mcmcSAR | R Documentation |
mcmcSAR
implements the Bayesian estimator of the linear-in-mean SAR model when only the linking probabilities are available or can be estimated.
mcmcSAR(
formula,
contextual,
start,
G0.obs,
G0 = NULL,
mlinks = list(),
hyperparms = list(),
ctrl.mcmc = list(),
iteration = 2000L,
data
)
formula |
object of class formula: a symbolic description of the model. The |
contextual |
(optional) logical; if true, this means that all individual variables will be set as contextual variables. Set
|
start |
(optional) vector of starting value of the model parameter as |
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 |
list of sub-network matrices (or simply network matrix if there is only one sub-network). |
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 |
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.
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).
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).
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.
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
.
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.
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.
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 |
mlinks |
return value of |
accept.rate |
acceptance rates. |
prop.net |
proportion of observed network data. |
method.net |
network formation model specification. |
start |
starting values. |
formula |
input value of |
contextual |
input value of |
ctrl.mcmc |
return value of |
smmSAR
, sim.IV
# 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")
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.