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
McCormick and Zheng (2015) write the likelihood of the model with respect to the spherical coordinate ,
the trait locations
, the degree
, the fraction of ties in the network that are
made with members of group k
, the trait intensity parameter
and
.
The following
prior distributions are defined.
Identification
For identification, some and
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 while
fixb
can be used to set the desired values for .
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,
,
,
and
is 0.44.
jumpmin
: The default value isc(0,1,1e-7,1e-7,1e-7)*1e-5
. The minimal jumping of everyis 0, every
is
, and every
,
and
is
.
jumpmax
: The default value isc(100,1,1,1,1)*20
. The maximal jumping scale is 20 except forwhich 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")
})
}