mcmcARD {PartialNetwork} | R Documentation |
Estimate network model using ARD
Description
mcmcARD
estimates the network model proposed by Breza et al. (2020).
Usage
mcmcARD(
Y,
traitARD,
start,
fixv,
consb,
iteration = 2000L,
sim.d = TRUE,
sim.zeta = TRUE,
hyperparms = NULL,
ctrl.mcmc = list()
)
Arguments
Y |
is a matrix of ARD. The entry (i, k) is the number of i's friends having the trait k. |
traitARD |
is the matrix of traits for individuals with ARD. The entry (i, k) is equal to 1 if i has the trait k and 0 otherwise. |
start |
is a list containing starting values of |
fixv |
is a vector setting which location parameters are fixed for identifiability. These fixed positions are used to rotate the latent surface back to a common orientation at each iteration using a Procrustes transformation (see Section Identification in Details). |
consb |
is a vector of the subset of |
iteration |
is the number of MCMC steps to be performed. |
sim.d |
is logical indicating whether the degree |
sim.zeta |
is logical indicating whether the degree |
hyperparms |
is an 8-dimensional vector of hyperparameters (in this order) |
ctrl.mcmc |
is a list of MCMC controls (see Section MCMC control in Details). |
Details
The linking probability is given by
Model
P_{ij} \propto \exp(\nu_i + \nu_j + \zeta\mathbf{z}_i\mathbf{z}_j).
McCormick and Zheng (2015) write the likelihood of the model with respect to the spherical coordinate \mathbf{z}_i
,
the trait locations \mathbf{v}_k
, the degree d_i
, the fraction of ties in the network that are
made with members of group k b_k
, the trait intensity parameter \eta_k
and \zeta
.
The following
prior distributions are defined.
\mathbf{z}_i \sim Uniform ~ von ~ Mises-Fisher
\mathbf{v}_k \sim Uniform ~ von ~ Mises-Fisher
d_i \sim log-\mathcal{N}(\mu_d, \sigma_d)
b_k \sim log-\mathcal{N}(\mu_b, \sigma_b)
\eta_k \sim Gamma(\alpha_{\eta}, \beta_{\eta})
\zeta \sim Gamma(\alpha_{\zeta}, \beta_{\zeta})
Identification
For identification, some \mathbf{v}_k
and b_k
need to be exogenously fixed around their given starting value
(see McCormick and Zheng, 2015 for more details). The parameter fixv
can be used
to set the desired value for \mathbf{v}_k
while fixb
can be used to set the desired values for b_k
.
MCMC control
During the MCMC, the jumping scales are updated following Atchade and Rosenthal (2005) in order to target the acceptance rate of each parameter to the target
values. This
requires to set minimal and 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 isrep(0.44, 5)
. The target of every\mathbf{z}_i
,d_i
,b_k
,\eta_k
and\zeta
is 0.44.jumpmin
: The default value isc(0,1,1e-7,1e-7,1e-7)*1e-5
. The minimal jumping of every\mathbf{z}_i
is 0, everyd_i
is10^{-5}
, and everyb_k
,\eta_k
and\zeta
is10^{-12}
.jumpmax
: The default value isc(100,1,1,1,1)*20
. The maximal jumping scale is 20 except for\mathbf{z}_i
which is set to 2000.print
: A logical value which indicates if the MCMC progression should be printed in the console. The default value isTRUE
.
Value
A list consisting of:
n |
dimension of the sample with ARD. |
K |
number of traits. |
p |
hypersphere dimension. |
time |
elapsed time in second. |
iteration |
number of MCMC steps performed. |
simulations |
simulations from the posterior distribution. |
hyperparms |
return value of hyperparameters (updated and non updated). |
accept.rate |
list of acceptance rates. |
start |
starting values. |
ctrl.mcmc |
return value of |
Examples
# Sample size
N <- 500
# ARD parameters
genzeta <- 1
mu <- -1.35
sigma <- 0.37
K <- 12 # number of traits
P <- 3 # Sphere dimension
# Generate z (spherical coordinates)
genz <- rvMF(N,rep(0,P))
# Generate nu from a Normal distribution with parameters mu and sigma (The gregariousness)
gennu <- rnorm(N,mu,sigma)
# compute degrees
gend <- N*exp(gennu)*exp(mu+0.5*sigma^2)*exp(logCpvMF(P,0) - logCpvMF(P,genzeta))
# Link probabilities
Probabilities <- sim.dnetwork(gennu,gend,genzeta,genz)
# Adjacency matrix
G <- sim.network(Probabilities)
# Generate vk, the trait location
genv <- rvMF(K,rep(0,P))
# set fixed some vk distant
genv[1,] <- c(1,0,0)
genv[2,] <- c(0,1,0)
genv[3,] <- c(0,0,1)
# eta, the intensity parameter
geneta <-abs(rnorm(K,2,1))
# Build traits matrix
densityatz <- matrix(0,N,K)
for(k in 1:K){
densityatz[,k] <- dvMF(genz,genv[k,]*geneta[k])
}
trait <- matrix(0,N,K)
NK <- floor(runif(K, 0.8, 0.95)*colSums(densityatz)/apply(densityatz, 2, max))
for (k in 1:K) {
trait[,k] <- rbinom(N, 1, NK[k]*densityatz[,k]/sum(densityatz[,k]))
}
# print a percentage of people having a trait
colSums(trait)*100/N
# Build ARD
ARD <- G %*% trait
# generate b
genb <- numeric(K)
for(k in 1:K){
genb[k] <- sum(G[,trait[,k]==1])/sum(G)
}
############ ARD Posterior distribution ###################
# initialization
d0 <- exp(rnorm(N)); b0 <- exp(rnorm(K)); eta0 <- rep(1,K);
zeta0 <- 05; z0 <- matrix(rvMF(N,rep(0,P)),N); v0 <- matrix(rvMF(K,rep(0,P)),K)
# We need to fix some of the vk and bk for identification (see Breza et al. (2020) for details).
vfixcolumn <- 1:6
bfixcolumn <- c(3, 5)
b0[bfixcolumn] <- genb[bfixcolumn]
v0[vfixcolumn,] <- genv[vfixcolumn,]
start <- list("z" = z0, "v" = v0, "d" = d0, "b" = b0, "eta" = eta0, "zeta" = zeta0)
# MCMC
out <- mcmcARD(Y = ARD, traitARD = trait, start = start, fixv = vfixcolumn,
consb = bfixcolumn, iteration = 5000)
# plot simulations
# plot d
plot(out$simulations$d[,100], type = "l", col = "blue", ylab = "")
abline(h = gend[100], col = "red")
# plot coordinates of individuals
i <- 123 # individual 123
{
lapply(1:3, function(x) {
plot(out$simulations$z[i, x,] , type = "l", ylab = "", col = "blue", ylim = c(-1, 1))
abline(h = genz[i, x], col = "red")
})
}
# plot coordinates of traits
k <- 8
{
lapply(1:3, function(x) {
plot(out$simulations$v[k, x,] , type = "l", ylab = "", col = "blue", ylim = c(-1, 1))
abline(h = genv[k, x], col = "red")
})
}