mcmcSAR {PartialNetwork}

Bayesian Estimator of SAR model

Description

mcmcSAR implements the Bayesian estimator of the linear-in-mean SAR model when only the linking probabilities are available or can be estimated.

Usage

mcmcSAR(
  formula,
  contextual,
  start,
  G0.obs,
  G0 = NULL,
  mlinks = list(),
  hyperparms = list(),
  ctrl.mcmc = list(),
  iteration = 2000L,
  data
)

Arguments

formula

object of class formula: a symbolic description of the model. The formula should be as for example y ~ x1 + x2 | x1 + x2 where y is the endogenous vector, the listed variables before the pipe, x1, x2 are the individual exogenous variables and the listed variables after the pipe, x1, x2 are the contextual observable variables. Other formulas may be y ~ x1 + x2 for the model without contextual effects, y ~ -1 + x1 + x2 | x1 + x2 for the model without intercept, or y ~ x1 + x2 | x2 + x3 to allow the contextual variables to be different from the individual variables.

contextual

(optional) logical; if true, this means that all individual variables will be set as contextual variables. Set formula as y ~ x1 + x2 and contextual as TRUE is equivalent to set formula as y ~ x1 + x2 | x1 + x2.

start

(optional) vector of starting value of the model parameter as (\beta' ~ \gamma' ~ \alpha ~ \sigma^2)', where \beta is the individual variables parameter, \gamma is the contextual variables parameter, \alpha is the peer effect parameter and \sigma^2 the variance of the error term. If the start is missing, a Maximum Likelihood estimator will be used, where the network matrix is that given through the argument G0 (if provided) or generated from it distribution.

G0.obs

list of matrices (or simply matrix if the list contains only one matrix) indicating the part of the network data which is observed. If the (i,j)-th element of the m-th matrix is one, then the element at the same position in the network data will be considered as observed and will not be inferred in the MCMC. In contrast, if the (i,j)-th element of the m-th matrix is zero, the element at the same position in the network data will be considered as a starting value of the missing link which will be inferred. G0.obs can also take "none" when no part of the network data is observed (equivalent to the case where all the entries are zeros) and "all" when the network data is fully observed (equivalent to the case where all the entries are ones).

G0

list of sub-network matrices (or simply network matrix if there is only one sub-network). G0 is made up of starting values for the entries with missing network data and observed values for the entries with observed network data. G0 is optional when G0.obs = "none".

mlinks

list specifying the network formation model (see Section Network formation model in Details).

hyperparms

(optional) is a list of hyperparameters (see Section Hyperparameters in Details).

ctrl.mcmc

list of MCMC controls (see Section MCMC control in Details).

iteration

number of MCMC steps to be performed.

data

optional data frame, list or environment (or object coercible by as.data.frame to a data frame) containing the variables in the model. If missing, the variables are taken from environment(formula), typically the environment from which mcmcSAR is called.

Details

Outcome model

The model is given by

\mathbf{y} = \mathbf{X}\beta + \mathbf{G}\mathbf{X}\gamma + \alpha \mathbf{G}\mathbf{y} + \epsilon.

where

\epsilon \sim N(0, \sigma^2).

The parameters to estimate in this model are the matrix \mathbf{G}, the vectors \beta, \gamma and the scalar \alpha, \sigma. Prior distributions are assumed on \mathbf{A}, the adjacency matrix in which \mathbf{A}_{ij} = 1 if i is connected to j and \mathbf{A}_{ij} = 0 otherwise, and on \beta, \gamma, \alpha and \sigma^2.

\mathbf{A}_{ij} \sim Bernoulli(\mathbf{P}_{ij})

(\beta' ~ \gamma')'|\sigma^2 \sim \mathcal{N}(\mu_{\theta}, \sigma^2\Sigma_{\theta})

\zeta = \log\left(\frac{\alpha}{1 - \alpha}\right) \sim \mathcal{N}(\mu_{\zeta}, \sigma_{\zeta}^2)

\sigma^2 \sim IG(\frac{a}{2}, \frac{b}{2})

where \mathbf{P} is the linking probability. The linking probability is an hyperparameters that can be set fixed or updated using a network formation model.

Network formation model

The linking probability can be set fixed or updated using a network formation model. Information about how \mathbf{P} should be handled in in the MCMC can be set through the argument mlinks which should be a list with named elements. Divers specifications of network formation model are possible. The list assigned to mlist should include an element named model. The expected values of model are "none" (default value), "logit", "probit", and "latent space".

• "none" means that the network distribution \mathbf{P} is set fixed throughout the MCMC,

• "probit" or "logit" implies that the network distribution \mathbf{P} will be updated using a Probit or Logit model,

• "latent spate" means that \mathbf{P} will be updated following Breza et al. (2020).

Fixed network distribution

To set \mathbf{P} fixed, mlinks could contain,

• dnetwork, a list, where the m-th elements is the matrix of link probability in the m-th sub-network.

• model = "none" (optional as "none" is the default value).

Probit and Logit models

For the Probit and Logit specification as network formation model, the following elements could be declared in mlinks.

• model = "probit" or model = "logit".

• mlinks.formula object of class formula: a symbolic description of the Logit or Probit model. The formula should only specify the explanatory variables, as for example ~ x1 + x2, the variables x1 and x2 are the dyadic observable characteristics. Each variable should verify length(x) == sum(N^2 - N), where N is a vector of the number of individual in each sub-network. Indeed, x will be associated with the entries (1, 2); (1, 3); (1, 4); ...; (2, 1); (2, 3); (2, 4); ... of the linking probability and as so, in all the sub-networks. Functions mat.to.vec and vec.to.mat can be used to convert a list of dyadic variable as in matrix form to a format that suits mlinks.formula.

• weights (optional) is a vector of weights of observed entries. This is important to address the selection problem of observed entries. Default is a vector of ones.

• estimates (optional when a part of the network is observed) is a list containing rho, a vector of the estimates of the Probit or Logit parameters, and var.rho the covariance matrix of the estimator. These estimates can be automatically computed when a part of the network data is available. In this case, rho and the unobserved part of the network are updated without using the observed part of the network. The latter is assumed non-stochastic in the MCMC. In addition, if G0.obs = "none", estimates should also include N, a vector of the number of individuals in each sub-network.

• prior (optional) is a list containing rho, a vector of the prior beliefs on rho, and var.rho the prior covariance matrix of rho. This input is relevant only when the observed part of the network is used to update rho, i.e. only when estimates = NULL (so, either estimates or prior should be NULL).
  To understand the difference between estimates and prior, note that estimates includes initial estimates of rho and var.rho, meaning that the observed part of the network is not used in the MCMC to update rho. In contrast, prior contains the prior beliefs of the user, and therefore, rho is updated using this prior and information from the observed part of the network. In addition, if G0.obs = "none", prior should also include N, a vector of the number of individuals in each sub-network.

• mlinks.data optional data frame, list or environment (or object coercible by as.data.frame to a data frame) containing the dyadic observable characteristics If missing, the variables will be taken from environment(mlinks.formula), typically the environment from which mcmcARD is called.

Latent space models

The following element could be declared in mlinks.

• model = "latent space".

• estimates a list of objects of class mcmcARD, where the m-th element is Breza et al. (2020) estimator as returned by the function mcmcARD in the m-th sub-network.

• mlinks.data (required only when ARD are partially observed) is a list of matrices, where the m-th element is the variable matrix to use to compute distance between individuals (could be the list of traits) in the m-th sub-network. The distances will be used to compute gregariousness and coordinates for individuals without ARD by k-nearest neighbors approach.

• obsARD (required only when ARD are partially observed) is a list of logical vectors, where the i-th entry of the m-th vector indicates by TRUE or FALSE if the i-th individual in the m-th sub-network has ARD or not.

• mARD (optional, default value is ⁠rep(1, M⁠)) is a vector indicating the number of neighbors to use in each sub-network.

• burninARD (optional) set the burn-in to summarize the posterior distribution in estimates.

Hyperparameters

All the hyperparameters can be defined through the argument hyperparms (a list) and should be named as follow.

• mutheta, the prior mean of (\beta' ~ \gamma')'|\sigma^2. The default value assumes that the prior mean is zero.

• invstheta as \Sigma_{\theta}^{-1}. The default value is a diagonal matrix with 0.01 on the diagonal.

• muzeta, the prior mean of \zeta. The default value is zero.

• invszeta, the inverse of the prior variance of \zeta with default value equal to 2.

• a and b which default values equal to 4.2 and 2.2 respectively. This means for example that the prior mean of \sigma^2 is 1.

Inverses are used for the prior variance through the argument hyperparms in order to allow non informative prior. Set the inverse of the prior variance to 0 is equivalent to assume a non informative prior.

MCMC control

During the MCMC, the jumping scales of \alpha and \rho are updated following Atchade and Rosenthal (2005) in order to target the acceptance rate to the target value. This requires to set a minimal and a maximal jumping scales through the parameter ctrl.mcmc. The parameter ctrl.mcmc is a list which can contain the following named components.

• target: the default value is c("alpha" = 0.44, "rho" = 0.234).

• jumpmin: the default value is c("alpha" = 1e-5, "rho" = 1e-5).

• jumpmax: the default value is c("alpha" = 10, "rho" = 10).

• print.level: an integer in {0, 1, 2} that indicates if the MCMC progression should be printed in the console. If 0, the MCMC progression is not be printed. If 1 (default value), the progression is printed and if 2, the simulations from the posterior distribution are printed.

• block.max: The maximal number of entries that can be updated simultaneously in \mathbf{A}. It might be more efficient to update simultaneously 2 or 3 entries (see Boucher and Houndetoungan, 2022).

If block.max > 1, several entries are randomly chosen from the same row and updated simultaneously. The number of entries chosen is randomly chosen between 1 and block.max. In addition, the entries are not chosen in order. For example, on the row i, the entries (i, 5) and (i, 9) can be updated simultaneously, then the entries (i, 1), (i, 3), (i, 8), and so on.

Value

A list consisting of:

n.group

number of groups.

N

vector of each group size.

time

elapsed time to run the MCMC in second.

iteration

number of MCMC steps performed.

posterior

matrix (or list of matrices) containing the simulations.

hyperparms

return value of hyperparms.

mlinks

return value of mlinks.

accept.rate

acceptance rates.

prop.net

proportion of observed network data.

method.net

network formation model specification.

start

starting values.

formula

input value of formula and mlinks.formula.

contextual

input value of contextual.

ctrl.mcmc

return value of ctrl.mcmc.

See Also

smmSAR, sim.IV

Examples

# We assume that the network is fully observed 
# See our vignette for examples where the network is partially observed 
# Number of groups
M             <- 50
# size of each group
N             <- rep(30,M)
# individual effects
beta          <- c(2,1,1.5)
# contextual effects
gamma         <- c(5,-3)
# endogenous effects
alpha         <- 0.4
# std-dev errors
se            <- 1
# prior distribution
prior         <- runif(sum(N*(N-1)))
prior         <- vec.to.mat(prior, N, normalise = FALSE)
# covariates
X             <- cbind(rnorm(sum(N),0,5),rpois(sum(N),7))
# true network
G0            <- sim.network(prior)
# normalise
G0norm        <- norm.network(G0)
# simulate dependent variable use an external package
y             <- CDatanet::simsar(~ X, contextual = TRUE, Glist = G0norm,
                                  theta = c(alpha, beta, gamma, se))
y             <- y$y
# dataset
dataset       <- as.data.frame(cbind(y, X1 = X[,1], X2 = X[,2]))
out.none1     <- mcmcSAR(formula = y ~ X1 + X2, contextual = TRUE, G0.obs = "all",
                         G0 = G0, data = dataset, iteration = 1e4)
summary(out.none1)
plot(out.none1)
plot(out.none1, plot.type = "dens")