Movatterモバイル変換

The CDatanet Package

Description

TheCDatanet package simulates and estimates peer effect models and network formation models. The peer effect models include linear-in-means models (Lee, 2004; Lee et al., 2010),Tobit models (Xu and Lee, 2015), and discrete numerical data models (Houndetoungan, 2024).The network formation models include pairwise regressions with degree heterogeneity (Graham, 2017; Yan et al., 2019) and exponential random graph models (Mele, 2017).To enhance computation speed,CDatanet usesC++ via theRcpp package (Eddelbuettel et al., 2011).

Author(s)

Maintainer: Aristide Houndetounganahoundetoungan@ecn.ulaval.ca

References

Eddelbuettel, D., & Francois, R. (2011).Rcpp: SeamlessR andC++ integration.Journal of Statistical Software, 40(8), 1-18,doi:10.18637/jss.v040.i08.

Houndetoungan, E. A. (2025). Count Data Models with Heterogeneous Peer Effects. Available at arXiv:2405.17290,doi:10.48550/arXiv.2405.17290.

Lee, L. F. (2004). Asymptotic distributions of quasi-maximum likelihood estimators for spatial autoregressive models.Econometrica, 72(6), 1899-1925,doi:10.1111/j.1468-0262.2004.00558.x.

Lee, L. F., Liu, X., & Lin, X. (2010). Specification and estimation of social interaction models with network structures. The Econometrics Journal, 13(2), 145-176,doi:10.1111/j.1368-423X.2010.00310.x

Xu, X., & Lee, L. F. (2015). Maximum likelihood estimation of a spatial autoregressive Tobit model.Journal of Econometrics, 188(1), 264-280,doi:10.1016/j.jeconom.2015.05.004.

Graham, B. S. (2017). An econometric model of network formation with degree heterogeneity.Econometrica, 85(4), 1033-1063,doi:10.3982/ECTA12679.

Mele, A. (2017). A structural model of dense network formation.Econometrica, 85(3), 825-850,doi:10.3982/ECTA10400.

Yan, T., Jiang, B., Fienberg, S. E., & Leng, C. (2019). Statistical inference in a directed network model with covariates.Journal of the American Statistical Association, 114(526), 857-868,doi:10.1080/01621459.2018.1448829.

See Also

Useful links:

• https://github.com/ahoundetoungan/CDatanet

• Report bugs athttps://github.com/ahoundetoungan/CDatanet/issues

Estimating Count Data Models with Social Interactions under Rational Expectations Using the NPL Method

Description

cdnet estimates count data models with social interactions under rational expectations using the NPL algorithm (see Houndetoungan, 2024).

Usage

cdnet(  formula,  Glist,  group,  Rmax,  Rbar,  starting = list(lambda = NULL, Gamma = NULL, delta = NULL),  Ey0 = NULL,  ubslambda = 1L,  optimizer = "fastlbfgs",  npl.ctr = list(),  opt.ctr = list(),  cov = TRUE,  data)

Arguments

Details

Model

The count variabley_i takes the valuer with probability.

P_{ir} = F(\sum_{s = 1}^S \lambda_s \bar{y}_i^{e,s} + \mathbf{z}_i'\Gamma - a_{h(i),r}) - F(\sum_{s = 1}^S \lambda_s \bar{y}_i^{e,s} + \mathbf{z}_i'\Gamma - a_{h(i),r + 1}).

In this equation,\mathbf{z}_i is a vector of control variables;F is the distribution function of the standard normal distribution;\bar{y}_i^{e,s} is the average ofE(y) among peers using thes-th network definition;a_{h(i),r} is ther-th cut-point in the cost grouph(i).

The following identification conditions have been introduced:\sum_{s = 1}^S \lambda_s > 0,a_{h(i),0} = -\infty,a_{h(i),1} = 0, anda_{h(i),r} = \infty for anyr \geq R_{\text{max}} + 1. The last condition implies thatP_{ir} = 0 for anyr \geq R_{\text{max}} + 1.For anyr \geq 1, the distance between two cut-points isa_{h(i),r+1} - a_{h(i),r} = \delta_{h(i),r} + \sum_{s = 1}^S \lambda_s.As the number of cut-points can be large, a quadratic cost function is considered forr \geq \bar{R}_{h(i)}, where\bar{R} = (\bar{R}_{1}, ..., \bar{R}_{L}).With the semi-parametric cost function,a_{h(i),r + 1} - a_{h(i),r} = \bar{\delta}_{h(i)} + \sum_{s = 1}^S \lambda_s.

The model parameters are:\lambda = (\lambda_1, ..., \lambda_S)',\Gamma, and\delta = (\delta_1', ..., \delta_L')',where\delta_l = (\delta_{l,2}, ..., \delta_{l,\bar{R}_l}, \bar{\delta}_l)' forl = 1, ..., L.The number of single parameters in\delta_l depends onR_{\text{max}} and\bar{R}_l. The components\delta_{l,2}, ..., \delta_{l,\bar{R}_l} or/and\bar{\delta}_l must be removed in certain cases.
IfR_{\text{max}} = \bar{R}_l \geq 2, then\delta_l = (\delta_{l,2}, ..., \delta_{l,\bar{R}_l})'.
IfR_{\text{max}} = \bar{R}_l = 1 (binary models), then\delta_l must be empty.
IfR_{\text{max}} > \bar{R}_l = 1, then\delta_l = \bar{\delta}_l.

npl.ctr

The model parameters are estimated using the Nested Partial Likelihood (NPL) method. This approachbegins with an initial guess for\theta andE(y) and iteratively refines them.The solution converges when the\ell_1-distance between two consecutive estimates of\theta andE(y) is smaller than a specified tolerance.

The argumentnpl.ctr must include the following parameters:

tol

the tolerance level for the NPL algorithm (default is 1e-4).

maxit

the maximum number of iterations allowed (default is 500).

print

a boolean value indicating whether the estimates should be printed at each step.

S

the number of simulations performed to compute the integral in the covariance using importance sampling.

Value

A list consisting of:

References

Houndetoungan, A. (2024). Count Data Models with Heterogeneous Peer Effects. Available at SSRN 3721250,doi:10.2139/ssrn.3721250.

See Also

sart,sar,simcdnet.

Examples

set.seed(123)M      <- 5 # Number of sub-groupsnvec   <- round(runif(M, 100, 200))n      <- sum(nvec)# Adjacency matrixA      <- list()for (m in 1:M) {  nm           <- nvec[m]  Am           <- matrix(0, nm, nm)  max_d        <- 30 #maximum number of friends  for (i in 1:nm) {    tmp        <- sample((1:nm)[-i], sample(0:max_d, 1))    Am[i, tmp] <- 1  }  A[[m]]       <- Am}Anorm  <- norm.network(A) #Row-normalization# XX      <- cbind(rnorm(n, 1, 3), rexp(n, 0.4))# Two group:group  <- 1*(X[,1] > 0.95)# Networks# length(group) = 2 and unique(sort(group)) = c(0, 1)# The networks must be defined as to capture:# peer effects of `0` on `0`, peer effects of `1` on `0`# peer effects of `0` on `1`, and peer effects of `1` on `1`G        <- list()cums     <- c(0, cumsum(nvec))for (m in 1:M) {  tp     <- group[(cums[m] + 1):(cums[m + 1])]  Am     <- A[[m]]  G[[m]] <- norm.network(list(Am * ((1 - tp) %*% t(1 - tp)),                              Am * ((1 - tp) %*% t(tp)),                              Am * (tp %*% t(1 - tp)),                              Am * (tp %*% t(tp))))}# Parameterslambda <- c(0.2, 0.3, -0.15, 0.25) Gamma  <- c(4.5, 2.2, -0.9, 1.5, -1.2)delta  <- rep(c(2.6, 1.47, 0.85, 0.7, 0.5), 2) # Datadata   <- data.frame(X, peer.avg(Anorm, cbind(x1 = X[,1], x2 =  X[,2])))colnames(data) = c("x1", "x2", "gx1", "gx2")ytmp   <- simcdnet(formula = ~ x1 + x2 + gx1 + gx2, Glist = G, Rbar = rep(5, 2),                   lambda = lambda, Gamma = Gamma, delta = delta, group = group,                   data = data)y      <- ytmp$yhist(y, breaks = max(y) + 1)table(y)# Estimationest    <- cdnet(formula = y ~ x1 + x2 + gx1 + gx2, Glist = G, Rbar = rep(5, 2), group = group,                optimizer = "fastlbfgs", data = data,                opt.ctr = list(maxit = 5e3, eps_f = 1e-11, eps_g = 1e-11))summary(est)

Converting Data between Directed Network Models and Symmetric Network Models.

Description

homophili.data converts the matrix of explanatory variables between directed network models and symmetric network models.

Usage

homophili.data(data, nvec, to = c("lower", "upper", "symmetric"))

Arguments

Value

The transformeddata.frame.

Estimating Network Formation Models with Degree Heterogeneity: the Fixed Effect Approach

Description

homophily.fe implements a Logit estimator for a network formation model with homophily. The model includes degree heterogeneity using fixed effects (see details).

Usage

homophily.fe(  network,  formula,  data,  symmetry = FALSE,  fe.way = 1,  init = NULL,  method = c("L-BFGS", "Block-NRaphson", "Mix"),  ctr = list(maxit.opt = 10000, maxit.nr = 50, eps_f = 1e-09, eps_g = 1e-09, tol = 1e-04),  print = TRUE)

Arguments

Details

Letp_{ij} be the probability for a link to go from individuali to individualj.This probability is specified for two-way effect models (fe.way = 2) as

p_{ij} = F(\mathbf{x}_{ij}'\beta + \mu_i + \nu_j),

whereF is the cumulative distribution function of the standard logistic distribution. Unobserved degree heterogeneity is captured by\mu_i and\nu_j. These are treated as fixed effects (seehomophily.re for random effect models).As shown by Yan et al. (2019), the estimator of the parameter\beta is biased. A bias correction is necessary but not implemented in this version. However,the estimators of\mu_i and\nu_j are consistent.

For one-way fixed effect models (fe.way = 1),\nu_j = \mu_j. For symmetric models, the network is not directed, and the fixed effects need to be one-way.

Value

A list consisting of:

References

Yan, T., Jiang, B., Fienberg, S. E., & Leng, C. (2019). Statistical inference in a directed network model with covariates.Journal of the American Statistical Association, 114(526), 857-868,doi:10.1080/01621459.2018.1448829.

See Also

homophily.re.

Examples

set.seed(1234)M            <- 2 # Number of sub-groupsnvec         <- round(runif(M, 20, 50))beta         <- c(.1, -.1)Glist        <- list()dX           <- matrix(0, 0, 2)mu           <- list()nu           <- list()Emunu        <- runif(M, -1.5, 0) # Expectation of mu + nusmu2         <- 0.2snu2         <- 0.2for (m in 1:M) {  n          <- nvec[m]  mum        <- rnorm(n, 0.7*Emunu[m], smu2)  num        <- rnorm(n, 0.3*Emunu[m], snu2)  X1         <- rnorm(n, 0, 1)  X2         <- rbinom(n, 1, 0.2)  Z1         <- matrix(0, n, n)    Z2         <- matrix(0, n, n)    for (i in 1:n) {    for (j in 1:n) {      Z1[i, j] <- abs(X1[i] - X1[j])      Z2[i, j] <- 1*(X2[i] == X2[j])    }  }    Gm           <- 1*((Z1*beta[1] + Z2*beta[2] +                       kronecker(mum, t(num), "+") + rlogis(n^2)) > 0)  diag(Gm)     <- 0  diag(Z1)     <- NA  diag(Z2)     <- NA  Z1           <- Z1[!is.na(Z1)]  Z2           <- Z2[!is.na(Z2)]    dX           <- rbind(dX, cbind(Z1, Z2))  Glist[[m]]   <- Gm  mu[[m]]      <- mum  nu[[m]]      <- num}mu  <- unlist(mu)nu  <- unlist(nu)out   <- homophily.fe(network =  Glist, formula = ~ -1 + dX, fe.way = 2)muhat <- out$estimate$munuhat <- out$estimate$nuplot(mu, muhat)plot(nu, nuhat)

Estimating Network Formation Models with Degree Heterogeneity: the Bayesian Random Effect Approach

Description

homophily.re implements a Bayesian Probit estimator for network formation model with homophily. The model includes degree heterogeneity using random effects (see details).

Usage

homophily.re(  network,  formula,  data,  symmetry = FALSE,  group.fe = FALSE,  re.way = 1,  init = list(),  iteration = 1000,  print = TRUE)

Arguments

Details

Letp_{ij} be a probability for a link to go from the individuali to the individualj.This probability is specified for two-way effect models (re.way = 2) as

p_{ij} = F(\mathbf{x}_{ij}'\beta + \mu_i + \nu_j),

whereF is the cumulative of the standard normal distribution. Unobserved degree heterogeneity is captured by\mu_i and\nu_j. The latter are treated as random effects (seehomophily.fe for fixed effect models).
For one-way random effect models (re.way = 1),\nu_j = \mu_j. For symmetric models, the network is not directed and therandom effects need to be one way.

Value

A list consisting of:

See Also

homophily.fe.

Examples

set.seed(1234)library(MASS)M            <- 4 # Number of sub-groupsnvec         <- round(runif(M, 100, 500))beta         <- c(.1, -.1)Glist        <- list()dX           <- matrix(0, 0, 2)mu           <- list()nu           <- list()cst          <- runif(M, -1.5, 0)smu2         <- 0.2snu2         <- 0.2rho          <- 0.8Smunu        <- matrix(c(smu2, rho*sqrt(smu2*snu2), rho*sqrt(smu2*snu2), snu2), 2)for (m in 1:M) {  n          <- nvec[m]  tmp        <- mvrnorm(n, c(0, 0), Smunu)  mum        <- tmp[,1] - mean(tmp[,1])  num        <- tmp[,2] - mean(tmp[,2])  X1         <- rnorm(n, 0, 1)  X2         <- rbinom(n, 1, 0.2)  Z1         <- matrix(0, n, n)    Z2         <- matrix(0, n, n)    for (i in 1:n) {    for (j in 1:n) {      Z1[i, j] <- abs(X1[i] - X1[j])      Z2[i, j] <- 1*(X2[i] == X2[j])    }  }    Gm           <- 1*((cst[m] + Z1*beta[1] + Z2*beta[2] +                       kronecker(mum, t(num), "+") + rnorm(n^2)) > 0)  diag(Gm)     <- 0  diag(Z1)     <- NA  diag(Z2)     <- NA  Z1           <- Z1[!is.na(Z1)]  Z2           <- Z2[!is.na(Z2)]    dX           <- rbind(dX, cbind(Z1, Z2))  Glist[[m]]   <- Gm  mu[[m]]      <- mum  nu[[m]]      <- num}mu  <- unlist(mu)nu  <- unlist(nu)out   <- homophily.re(network =  Glist, formula = ~ dX, group.fe = TRUE,                       re.way = 2, iteration = 1e3)# plot simulationsplot(out$posterior$beta[,1], type = "l")abline(h = cst[1], col = "red")plot(out$posterior$beta[,2], type = "l")abline(h = cst[2], col = "red")plot(out$posterior$beta[,3], type = "l")abline(h = cst[3], col = "red")plot(out$posterior$beta[,4], type = "l")abline(h = cst[4], col = "red")plot(out$posterior$beta[,5], type = "l")abline(h = beta[1], col = "red")plot(out$posterior$beta[,6], type = "l")abline(h = beta[2], col = "red")plot(out$posterior$sigma2_mu, type = "l")abline(h = smu2, col = "red")plot(out$posterior$sigma2_nu, type = "l")abline(h = snu2, col = "red")plot(out$posterior$rho, type = "l")abline(h = rho, col = "red")i <- 10plot(out$posterior$mu[,i], type = "l")abline(h = mu[i], col = "red")plot(out$posterior$nu[,i], type = "l")abline(h = nu[i], col = "red")

Marginal Effects for Count Data Models and Tobit Models with Social Interactions

Description

meffects computes marginal effects for count data and Tobit models with social interactions.It is a generic function which means that new printing methods can be easily added for new classes.

Usage

meffects(model, ...)## S3 method for class 'cdnet'meffects(  model,  Glist,  cont.var,  bin.var,  type.var,  Glist.contextual,  data,  tol = 1e-10,  maxit = 500,  boot = 1000,  progress = TRUE,  ncores = 1,  ...)## S3 method for class 'summary.cdnet'meffects(  model,  Glist,  cont.var,  bin.var,  type.var,  Glist.contextual,  data,  tol = 1e-10,  maxit = 500,  boot = 1000,  progress = TRUE,  ncores = 1,  ...)## S3 method for class 'sart'meffects(  model,  Glist,  cont.var,  bin.var,  type.var,  Glist.contextual,  data,  tol = 1e-10,  maxit = 500,  boot = 1000,  progress = TRUE,  ncores = 1,  ...)## S3 method for class 'summary.sart'meffects(  model,  Glist,  cont.var,  bin.var,  type.var,  Glist.contextual,  data,  tol = 1e-10,  maxit = 500,  boot = 1000,  progress = TRUE,  ncores = 1,  ...)

Arguments

Value

A list containing:

info

General information about the model.

estimate

The Maximum Likelihood (ML) estimates of the parameters.

Ey

E(y), the expected values of the endogenous variable.

GEy

The average ofE(y) among peers.

cov

A list containing covariance matrices (ifcov = TRUE).

details

Additional outputs returned by the optimizer.

meffects

A list containing the marginal effects.

Examples

#' set.seed(123)M      <- 5 # Number of sub-groupsnvec   <- round(runif(M, 100, 200))n      <- sum(nvec)# Adjacency matrixA      <- list()for (m in 1:M) {  nm           <- nvec[m]  Am           <- matrix(0, nm, nm)  max_d        <- 30 #maximum number of friends  for (i in 1:nm) {    tmp        <- sample((1:nm)[-i], sample(0:max_d, 1))    Am[i, tmp] <- 1  }  A[[m]]       <- Am}Anorm  <- norm.network(A) #Row-normalization# XX      <- cbind(rnorm(n, 1, 3), rexp(n, 0.4))# Two group:group  <- 1*(X[,1] > 0.95)# Networks# length(group) = 2 and unique(sort(group)) = c(0, 1)# The networks must be defined as to capture:# peer effects of `0` on `0`, peer effects of `1` on `0`# peer effects of `0` on `1`, and peer effects of `1` on `1`G        <- list()cums     <- c(0, cumsum(nvec))for (m in 1:M) {  tp     <- group[(cums[m] + 1):(cums[m + 1])]  Am     <- A[[m]]  G[[m]] <- norm.network(list(Am * ((1 - tp) %*% t(1 - tp)),                              Am * ((1 - tp) %*% t(tp)),                              Am * (tp %*% t(1 - tp)),                              Am * (tp %*% t(tp))))}# Parameterslambda <- c(0.2, 0.3, -0.15, 0.25) Gamma  <- c(4.5, 2.2, -0.9, 1.5, -1.2)delta  <- rep(c(2.6, 1.47, 0.85, 0.7, 0.5), 2) # Datadata   <- data.frame(X, peer.avg(Anorm, cbind(x1 = X[,1], x2 =  X[,2])))colnames(data) = c("x1", "x2", "gx1", "gx2")ytmp   <- simcdnet(formula = ~ x1 + x2 + gx1 + gx2, Glist = G, Rbar = rep(5, 2),                   lambda = lambda, Gamma = Gamma, delta = delta, group = group,                   data = data)y      <- ytmp$yhist(y, breaks = max(y) + 1)table(y)# Estimationest    <- cdnet(formula = y ~ x1 + x2 + gx1 + gx2, Glist = G, Rbar = rep(5, 2), group = group,                optimizer = "fastlbfgs", data = data,                opt.ctr = list(maxit = 5e3, eps_f = 1e-11, eps_g = 1e-11))meffects(est, Glist = G, data = data, cont.var = c("x1", "x2", "gx1", "gx2"),         type.var = list(c("x1", "gx1"), c("x2", "gx2")), Glist.contextual = Anorm,         boot = 100, ncores = 2)

Creating Objects for Network Models

Description

Thevec.to.mat function creates a list of square matrices from a given vector.Elements of the generated matrices are taken from the vector and placed column-wise or row-wise, progressing from the first matrix in the list to the last.The diagonals of the generated matrices are set to zeros.
Themat.to.vec function creates a vector from a given list of square matrices.Elements of the generated vector are taken column-wise or row-wise, starting from the first matrix in the list to the last, excluding diagonal entries.
Thenorm.network function row-normalizes matrices in a given list.

Usage

norm.network(W)vec.to.mat(u, N, normalise = FALSE, byrow = FALSE)mat.to.vec(W, ceiled = FALSE, byrow = FALSE)

Arguments

Value

A vector of sizesum(N * (N - 1)) or a list oflength(N) square matrices, with matrix sizes determined by⁠N[1], N[2], ...⁠.

See Also

simnetwork,peer.avg.

Examples

# Generate a list of adjacency matrices## Sub-network sizesN <- c(250, 370, 120)  ## Rate of friendshipp <- c(0.2, 0.15, 0.18)   ## Network datau <- unlist(lapply(1:3, function(x) rbinom(N[x] * (N[x] - 1), 1, p[x])))W <- vec.to.mat(u, N)# Convert W into a list of row-normalized matricesG <- norm.network(W)# Recover uv <- mat.to.vec(G, ceiled = TRUE)all.equal(u, v)

Computing Peer Averages

Description

Thepeer.avg function computes peer average values using network data (provided as a list of adjacency matrices) and observable characteristics.

Usage

peer.avg(Glist, V, export.as.list = FALSE)

Arguments

Value

The matrix productdiag(Glist[[1]], Glist[[2]], ...) %*% V, wherediag() represents the block diagonal operator.

See Also

simnetwork,vec.to.mat

Examples

# Generate a list of adjacency matrices## Sub-network sizesN <- c(250, 370, 120)  ## Rate of friendshipp <- c(0.2, 0.15, 0.18)   ## Network datau <- unlist(lapply(1:3, function(x) rbinom(N[x] * (N[x] - 1), 1, p[x])))G <- vec.to.mat(u, N, normalise = TRUE)# Generate a vector yy <- rnorm(sum(N))# Compute G %*% yGy <- peer.avg(Glist = G, V = y)

Printing the Average Expected Outcomes for Count Data Models with Social Interactions

Description

Summary and print methods for the classsimcdEy as returned by the functionsimcdEy.

Usage

## S3 method for class 'simcdEy'print(x, ...)## S3 method for class 'simcdEy'summary(object, ...)## S3 method for class 'summary.simcdEy'print(x, ...)

Arguments

Value

A list of the same objects inobject.

Removing Identifiers with NA from Adjacency Matrices Optimally

Description

Theremove.ids function removes identifiers with missing values (NA) from adjacency matrices in an optimal way.Multiple combinations of rows and columns can be deleted to eliminate NAs, but this function ensures that the smallestnumber of rows and columns are removed to retain as much data as possible.

Usage

remove.ids(network, ncores = 1L)

Arguments

Value

A list containing:

network

A list of adjacency matrices without missing values.

id

A list of vectors indicating the indices of retained rows and columns for each matrix.

Examples

# Example 1: Small adjacency matrixA <- matrix(1:25, 5)A[1, 1] <- NAA[4, 2] <- NAremove.ids(A)# Example 2: Larger adjacency matrix with multiple NAsB <- matrix(1:100, 10)B[1, 1] <- NAB[4, 2] <- NAB[2, 4] <- NAB[, 8] <- NAremove.ids(B)

Estimating Linear-in-mean Models with Social Interactions

Description

sar computes quasi-maximum likelihood estimators for linear-in-mean models with social interactions (see Lee, 2004 and Lee et al., 2010).

Usage

sar(  formula,  Glist,  lambda0 = NULL,  fixed.effects = FALSE,  optimizer = "optim",  opt.ctr = list(),  print = TRUE,  cov = TRUE,  cinfo = TRUE,  data)

Arguments

Details

In the complete information model, the outcomey_i for individuali is defined as:

y_i = \lambda \bar{y}_i + \mathbf{z}_i'\Gamma + \epsilon_i,

where\bar{y}_i represents the average outcomey among individuali's peers,\mathbf{z}_i is a vector of control variables, and\epsilon_i \sim N(0, \sigma^2) is the error term.In the case of incomplete information models with rational expectations, the outcomey_i is defined as:

y_i = \lambda E(\bar{y}_i) + \mathbf{z}_i'\Gamma + \epsilon_i,

whereE(\bar{y}_i) is the expected average outcome ofi's peers, as perceived by individuali.

Value

A list consisting of:

References

Lee, L. F. (2004). Asymptotic distributions of quasi-maximum likelihood estimators for spatial autoregressive models.Econometrica, 72(6), 1899-1925,doi:10.1111/j.1468-0262.2004.00558.x.

Lee, L. F., Liu, X., & Lin, X. (2010). Specification and estimation of social interaction models with network structures. The Econometrics Journal, 13(2), 145-176,doi:10.1111/j.1368-423X.2010.00310.x

See Also

sart,cdnet,simsar.

Examples

# Groups' sizeset.seed(123)M      <- 5 # Number of sub-groupsnvec   <- round(runif(M, 100, 1000))n      <- sum(nvec)# Parameterslambda <- 0.4Gamma  <- c(2, -1.9, 0.8, 1.5, -1.2)sigma  <- 1.5theta  <- c(lambda, Gamma, sigma)# XX      <- cbind(rnorm(n, 1, 1), rexp(n, 0.4))# NetworkG      <- list()for (m in 1:M) {  nm           <- nvec[m]  Gm           <- matrix(0, nm, nm)  max_d        <- 30  for (i in 1:nm) {    tmp        <- sample((1:nm)[-i], sample(0:max_d, 1))    Gm[i, tmp] <- 1  }  rs           <- rowSums(Gm); rs[rs == 0] <- 1  Gm           <- Gm/rs  G[[m]]       <- Gm}# datadata   <- data.frame(X, peer.avg(G, cbind(x1 = X[,1], x2 =  X[,2])))colnames(data) <- c("x1", "x2", "gx1", "gx2")ytmp    <- simsar(formula = ~ x1 + x2 + gx1 + gx2, Glist = G,                   theta = theta, data = data) data$y  <- ytmp$yout     <- sar(formula = y ~ x1 + x2 + gx1 + gx2, Glist = G,                optimizer = "optim", data = data)summary(out)

Estimating Tobit Models with Social Interactions

Description

sart estimates Tobit models with social interactions based on the framework of Xu and Lee (2015).The method allows for modeling both complete and incomplete information scenarios in networks, incorporating rational expectations in the latter case.

Usage

sart(  formula,  Glist,  starting = NULL,  Ey0 = NULL,  optimizer = "fastlbfgs",  npl.ctr = list(),  opt.ctr = list(),  cov = TRUE,  cinfo = TRUE,  data)

Arguments

Details

For a complete information model, the outcomey_i is defined as:

\begin{cases}y_i^{\ast} = \lambda \bar{y}_i + \mathbf{z}_i'\Gamma + \epsilon_i, \\ y_i = \max(0, y_i^{\ast}),\end{cases}

where\bar{y}_i is the average ofy among peers,\mathbf{z}_i is a vector of control variables,and\epsilon_i \sim N(0, \sigma^2).

In the case of incomplete information models with rational expectations,y_i is defined as:

\begin{cases}y_i^{\ast} = \lambda E(\bar{y}_i) + \mathbf{z}_i'\Gamma + \epsilon_i, \\ y_i = \max(0, y_i^{\ast}).\end{cases}

Value

A list containing:

info

General information about the model.

estimate

The Maximum Likelihood (ML) estimates of the parameters.

Ey

E(y), the expected values of the endogenous variable.

GEy

The average ofE(y) among peers.

cov

A list including covariance matrices (ifcov = TRUE).

details

Additional outputs returned by the optimizer.

References

Xu, X., & Lee, L. F. (2015). Maximum likelihood estimation of a spatial autoregressive Tobit model.Journal of Econometrics, 188(1), 264-280,doi:10.1016/j.jeconom.2015.05.004.

See Also

sar,cdnet,simsart.

Examples

# Group sizesset.seed(123)M      <- 5 # Number of sub-groupsnvec   <- round(runif(M, 100, 200))n      <- sum(nvec)# Parameterslambda <- 0.4Gamma  <- c(2, -1.9, 0.8, 1.5, -1.2)sigma  <- 1.5theta  <- c(lambda, Gamma, sigma)# Covariates (X)X      <- cbind(rnorm(n, 1, 1), rexp(n, 0.4))# Network creationG      <- list()for (m in 1:M) {  nm           <- nvec[m]  Gm           <- matrix(0, nm, nm)  max_d        <- 30  for (i in 1:nm) {    tmp        <- sample((1:nm)[-i], sample(0:max_d, 1))    Gm[i, tmp] <- 1  }  rs           <- rowSums(Gm); rs[rs == 0] <- 1  Gm           <- Gm / rs  G[[m]]       <- Gm}# Data creationdata   <- data.frame(X, peer.avg(G, cbind(x1 = X[, 1], x2 = X[, 2])))colnames(data) <- c("x1", "x2", "gx1", "gx2")## Complete information gameytmp    <- simsart(formula = ~ x1 + x2 + gx1 + gx2, Glist = G, theta = theta,                    data = data, cinfo = TRUE)data$yc <- ytmp$y## Incomplete information gameytmp    <- simsart(formula = ~ x1 + x2 + gx1 + gx2, Glist = G, theta = theta,                    data = data, cinfo = FALSE)data$yi <- ytmp$y# Complete information estimation for ycoutc1   <- sart(formula = yc ~ x1 + x2 + gx1 + gx2, optimizer = "nlm",                Glist = G, data = data, cinfo = TRUE)summary(outc1)# Complete information estimation for yioutc1   <- sart(formula = yi ~ x1 + x2 + gx1 + gx2, optimizer = "nlm",                Glist = G, data = data, cinfo = TRUE)summary(outc1)# Incomplete information estimation for ycouti1   <- sart(formula = yc ~ x1 + x2 + gx1 + gx2, optimizer = "nlm",                Glist = G, data = data, cinfo = FALSE)summary(outi1)# Incomplete information estimation for yiouti1   <- sart(formula = yi ~ x1 + x2 + gx1 + gx2, optimizer = "nlm",                Glist = G, data = data, cinfo = FALSE)summary(outi1)

Counterfactual Analyses with Count Data Models and Social Interactions

Description

simcdpar computes the average expected outcomes for count data models with social interactions and standard errors using the Delta method.This function can be used to examine the effects of changes in the network or in the control variables.

Usage

simcdEy(object, Glist, data, group, tol = 1e-10, maxit = 500, S = 1000)

Arguments

Value

A list consisting of:

See Also

simcdnet

Simulating Count Data Models with Social Interactions Under Rational Expectations

Description

simcdnet simulates the count data model with social interactions under rational expectations developed by Houndetoungan (2024).

Usage

simcdnet(  formula,  group,  Glist,  parms,  lambda,  Gamma,  delta,  Rmax,  Rbar,  cont.var,  bin.var,  tol = 1e-10,  maxit = 500,  data)

Arguments

Details

The count variabley_i takes the valuer with probability.

P_{ir} = F(\sum_{s = 1}^S \lambda_s \bar{y}_i^{e,s} + \mathbf{z}_i'\Gamma - a_{h(i),r}) - F(\sum_{s = 1}^S \lambda_s \bar{y}_i^{e,s} + \mathbf{z}_i'\Gamma - a_{h(i),r + 1}).

In this equation,\mathbf{z}_i is a vector of control variables;F is the distribution function of the standard normal distribution;\bar{y}_i^{e,s} is the average ofE(y) among peers using thes-th network definition;a_{h(i),r} is ther-th cut-point in the cost grouph(i).

The following identification conditions have been introduced:\sum_{s = 1}^S \lambda_s > 0,a_{h(i),0} = -\infty,a_{h(i),1} = 0, anda_{h(i),r} = \infty for anyr \geq R_{\text{max}} + 1. The last condition implies thatP_{ir} = 0 for anyr \geq R_{\text{max}} + 1.For anyr \geq 1, the distance between two cut-points isa_{h(i),r+1} - a_{h(i),r} = \delta_{h(i),r} + \sum_{s = 1}^S \lambda_s.As the number of cut-points can be large, a quadratic cost function is considered forr \geq \bar{R}_{h(i)}, where\bar{R} = (\bar{R}_{1}, ..., \bar{R}_{L}).With the semi-parametric cost function,a_{h(i),r + 1} - a_{h(i),r} = \bar{\delta}_{h(i)} + \sum_{s = 1}^S \lambda_s.

The model parameters are:\lambda = (\lambda_1, ..., \lambda_S)',\Gamma, and\delta = (\delta_1', ..., \delta_L')',where\delta_l = (\delta_{l,2}, ..., \delta_{l,\bar{R}_l}, \bar{\delta}_l)' forl = 1, ..., L.The number of single parameters in\delta_l depends onR_{\text{max}} and\bar{R}_l. The components\delta_{l,2}, ..., \delta_{l,\bar{R}_l} or/and\bar{\delta}_l must be removed in certain cases.
IfR_{\text{max}} = \bar{R}_l \geq 2, then\delta_l = (\delta_{l,2}, ..., \delta_{l,\bar{R}_l})'.
IfR_{\text{max}} = \bar{R}_l = 1 (binary models), then\delta_l must be empty.
IfR_{\text{max}} > \bar{R}_l = 1, then\delta_l = \bar{\delta}_l.

Value

A list consisting of:

References

Houndetoungan, A. (2024). Count Data Models with Heterogeneous Peer Effects. Available at SSRN 3721250,doi:10.2139/ssrn.3721250.

See Also

cdnet,simsart,simsar.

Examples

set.seed(123)M      <- 5 # Number of sub-groupsnvec   <- round(runif(M, 100, 200)) # Random group sizesn      <- sum(nvec) # Total number of individuals# Adjacency matrix for each groupA      <- list()for (m in 1:M) {  nm           <- nvec[m] # Size of group m  Am           <- matrix(0, nm, nm) # Empty adjacency matrix  max_d        <- 30 # Maximum number of friends  for (i in 1:nm) {    tmp        <- sample((1:nm)[-i], sample(0:max_d, 1)) # Sample friends    Am[i, tmp] <- 1 # Set friendship links  }  A[[m]]       <- Am # Add to the list}Anorm  <- norm.network(A) # Row-normalization of the adjacency matrices# Covariates (X)X      <- cbind(rnorm(n, 1, 3), rexp(n, 0.4)) # Random covariates# Two groups based on first covariategroup  <- 1 * (X[,1] > 0.95) # Assign to groups based on x1# Networks: Define peer effects based on group membership# The networks should capture:# - Peer effects of `0` on `0`# - Peer effects of `1` on `0`# - Peer effects of `0` on `1`# - Peer effects of `1` on `1`G        <- list()cums     <- c(0, cumsum(nvec)) # Cumulative indices for groupsfor (m in 1:M) {  tp     <- group[(cums[m] + 1):(cums[m + 1])] # Group membership for group m  Am     <- A[[m]] # Adjacency matrix for group m  # Define networks based on peer effects  G[[m]] <- norm.network(list(Am * ((1 - tp) %*% t(1 - tp)),                              Am * ((1 - tp) %*% t(tp)),                              Am * (tp %*% t(1 - tp)),                              Am * (tp %*% t(tp))))}# Parameters for the modellambda <- c(0.2, 0.3, -0.15, 0.25) Gamma  <- c(4.5, 2.2, -0.9, 1.5, -1.2)delta  <- rep(c(2.6, 1.47, 0.85, 0.7, 0.5), 2) # Repeated values for delta# Prepare data for the modeldata   <- data.frame(X, peer.avg(Anorm, cbind(x1 = X[,1], x2 = X[,2]))) colnames(data) = c("x1", "x2", "gx1", "gx2") # Set column names# Simulate outcomes using the `simcdnet` functionytmp   <- simcdnet(formula = ~ x1 + x2 + gx1 + gx2, Glist = G, Rbar = rep(5, 2),                   lambda = lambda, Gamma = Gamma, delta = delta, group = group,                   data = data)y      <- ytmp$y# Plot histogram of the simulated outcomeshist(y, breaks = max(y) + 1)# Display frequency table of the simulated outcomestable(y)

Simulating Network Data

Description

simnetwork generates adjacency matrices based on specified probabilities.

Usage

simnetwork(dnetwork, normalise = FALSE)

Arguments

Value

A list of (row-normalized) adjacency matrices.

Examples

# Generate a list of adjacency matrices## Sub-network sizesN         <- c(250, 370, 120)  ## Probability distributionsdnetwork  <- lapply(N, function(x) matrix(runif(x^2), x))## Generate networksG         <- simnetwork(dnetwork)

Simulating Data from Linear-in-Mean Models with Social Interactions

Description

simsar simulates continuous variables under linear-in-mean models with social interactions, following the specifications describedin Lee (2004) and Lee et al. (2010). The model incorporates peer interactions, where the value of an individual’s outcome dependsnot only on their own characteristics but also on the average characteristics of their peers in the network.

Usage

simsar(formula, Glist, theta, cinfo = TRUE, data)

Arguments

Details

In the complete information model, the outcomey_i for individuali is defined as:

y_i = \lambda \bar{y}_i + \mathbf{z}_i'\Gamma + \epsilon_i,

where\bar{y}_i represents the average outcomey among individuali's peers,\mathbf{z}_i is a vector of control variables, and\epsilon_i \sim N(0, \sigma^2) is the error term.In the case of incomplete information models with rational expectations, the outcomey_i is defined as:

y_i = \lambda E(\bar{y}_i) + \mathbf{z}_i'\Gamma + \epsilon_i,

whereE(\bar{y}_i) is the expected average outcome ofi's peers, as perceived by individuali.

Value

A list containing the following elements:

References

Lee, L. F. (2004). Asymptotic distributions of quasi-maximum likelihood estimators for spatial autoregressive models.Econometrica, 72(6), 1899-1925,doi:10.1111/j.1468-0262.2004.00558.x.

Lee, L. F., Liu, X., & Lin, X. (2010). Specification and estimation of social interaction models with network structures. The Econometrics Journal, 13(2), 145-176,doi:10.1111/j.1368-423X.2010.00310.x

See Also

sar,simsart,simcdnet.

Examples

# Groups' sizeset.seed(123)M      <- 5 # Number of sub-groupsnvec   <- round(runif(M, 100, 1000))n      <- sum(nvec)# Parameterslambda <- 0.4Gamma  <- c(2, -1.9, 0.8, 1.5, -1.2)sigma  <- 1.5theta  <- c(lambda, Gamma, sigma)# XX      <- cbind(rnorm(n, 1, 1), rexp(n, 0.4))# NetworkG      <- list()for (m in 1:M) {  nm           <- nvec[m]  Gm           <- matrix(0, nm, nm)  max_d        <- 30  for (i in 1:nm) {    tmp        <- sample((1:nm)[-i], sample(0:max_d, 1))    Gm[i, tmp] <- 1  }  rs           <- rowSums(Gm); rs[rs == 0] <- 1  Gm           <- Gm/rs  G[[m]]       <- Gm}# datadata   <- data.frame(X, peer.avg(G, cbind(x1 = X[,1], x2 =  X[,2])))colnames(data) <- c("x1", "x2", "gx1", "gx2")ytmp    <- simsar(formula = ~ x1 + x2 + gx1 + gx2, Glist = G,                   theta = theta, data = data) y       <- ytmp$y

Simulating Data from Tobit Models with Social Interactions

Description

simsart simulates censored data with social interactions (see Xu and Lee, 2015).

Usage

simsart(  formula,  Glist,  theta,  cont.var,  bin.var,  tol = 1e-15,  maxit = 500,  cinfo = TRUE,  data)

Arguments

Details

For a complete information model, the outcomey_i is defined as:

\begin{cases}y_i^{\ast} = \lambda \bar{y}_i + \mathbf{z}_i'\Gamma + \epsilon_i, \\ y_i = \max(0, y_i^{\ast}),\end{cases}

where\bar{y}_i is the average ofy among peers,\mathbf{z}_i is a vector of control variables,and\epsilon_i \sim N(0, \sigma^2).

In the case of incomplete information models with rational expectations,y_i is defined as:

\begin{cases}y_i^{\ast} = \lambda E(\bar{y}_i) + \mathbf{z}_i'\Gamma + \epsilon_i, \\ y_i = \max(0, y_i^{\ast}).\end{cases}

Value

A list consisting of:

yst

y^{\ast}, the latent variable.

y

The observed censored variable.

Ey

E(y), the expected value ofy.

Gy

The average ofy among peers.

GEy

The average ofE(y) among peers.

meff

A list including average and individual marginal effects.

iteration

The number of iterations performed per sub-network in the fixed-point iteration method.

References

Xu, X., & Lee, L. F. (2015). Maximum likelihood estimation of a spatial autoregressive Tobit model.Journal of Econometrics, 188(1), 264-280,doi:10.1016/j.jeconom.2015.05.004.

See Also

sart,simsar,simcdnet.

Examples

# Define group sizesset.seed(123)M      <- 5 # Number of sub-groupsnvec   <- round(runif(M, 100, 200)) # Number of nodes per sub-groupn      <- sum(nvec) # Total number of nodes# Define parameterslambda <- 0.4Gamma  <- c(2, -1.9, 0.8, 1.5, -1.2)sigma  <- 1.5theta  <- c(lambda, Gamma, sigma)# Generate covariates (X)X      <- cbind(rnorm(n, 1, 1), rexp(n, 0.4))# Construct network adjacency matricesG      <- list()for (m in 1:M) {  nm           <- nvec[m] # Nodes in sub-group m  Gm           <- matrix(0, nm, nm) # Initialize adjacency matrix  max_d        <- 30 # Maximum degree  for (i in 1:nm) {    tmp        <- sample((1:nm)[-i], sample(0:max_d, 1)) # Random connections    Gm[i, tmp] <- 1  }  rs           <- rowSums(Gm) # Normalize rows  rs[rs == 0]  <- 1  Gm           <- Gm / rs  G[[m]]       <- Gm}# Prepare datadata   <- data.frame(X, peer.avg(G, cbind(x1 = X[, 1], x2 = X[, 2])))colnames(data) <- c("x1", "x2", "gx1", "gx2") # Add column names# Complete information game simulationytmp    <- simsart(formula = ~ x1 + x2 + gx1 + gx2,                    Glist = G, theta = theta,                    data = data, cinfo = TRUE)data$yc <- ytmp$y # Add simulated outcome to the dataset# Incomplete information game simulationytmp    <- simsart(formula = ~ x1 + x2 + gx1 + gx2,                    Glist = G, theta = theta,                    data = data, cinfo = FALSE)data$yi <- ytmp$y # Add simulated outcome to the dataset

Summary for the Estimation of Count Data Models with Social Interactions under Rational Expectations

Description

Summary and print methods for the classcdnet as returned by the functioncdnet.

Usage

## S3 method for class 'cdnet'summary(object, Glist, data, S = 1000L, ...)## S3 method for class 'summary.cdnet'print(x, ...)## S3 method for class 'cdnet'print(x, ...)

Arguments

Value

A list of the same objects inobject.

Summary for the Estimation of Linear-in-mean Models with Social Interactions

Description

Summary and print methods for the classsar as returned by the functionsar.

Usage

## S3 method for class 'sar'summary(object, ...)## S3 method for class 'summary.sar'print(x, ...)## S3 method for class 'sar'print(x, ...)

Arguments

Value

A list of the same objects inobject.

Summary for the Estimation of Tobit Models with Social Interactions

Description

Summary and print methods for the classsart as returned by the functionsart.

Usage

## S3 method for class 'sart'summary(object, Glist, data, ...)## S3 method for class 'summary.sart'print(x, ...)## S3 method for class 'sart'print(x, ...)

Arguments

Value

A list of the same objects inobject.