uniNet: A function that generates samples for a univariate network model.

Description

This function generates samples for a univariate network model, which is given by

$$Y_{i_s}|\mu_{i_s} \sim f(y_{i_s}| \mu_{i_s}, \sigma_{e}^{2}) ~~~ i=1,\ldots, N_{s},~s=1,\ldots,S ,$$ $$g(\mu_{i_s}) = \boldsymbol{x}^\top_{i_s} \boldsymbol{\beta} + \sum_{j\in \textrm{net}(i_s)}w_{i_sj}u_{j} + w^{*}_{i_s}u^{*},$$ $$\boldsymbol{\beta} \sim \textrm{N}(\boldsymbol{0}, \alpha\boldsymbol{I}),$$ $$u_{j} \sim \textrm{N}(0, \sigma_{u}^{2}),$$ $$u^{*} \sim \textrm{N}(0, \sigma_{u}^{2}),$$ $$\sigma_{u}^{2} \sim \textrm{Inverse-Gamma}(\alpha_{2}, \xi_{2}),$$ $$\sigma_{e}^{2} \sim \textrm{Inverse-Gamma}(\alpha_{3}, \xi_{3}).$$

The covariates for the $i$th individual in the $s$th spatial unit or other grouping are included in a $p \times 1$ vector $\boldsymbol{x}_{i_s}$. The corresponding $p \times 1$ vector of fixed effect parameters are denoted by $\boldsymbol{\beta}$, which has an assumed multivariate Gaussian prior with mean $\boldsymbol{0}$ and diagonal covariance matrix $\alpha\boldsymbol{I}$ that can be chosen by the user. A conjugate Inverse-Gamma prior is specified for $\sigma_{e}^{2}$, and the corresponding hyperparamaterers ($\alpha_{3}$, $\xi_{3}$) can be chosen by the user.

The $J \times 1$ vector of alter random effects are denoted by $\boldsymbol{u} = (u_{1}, \ldots, u_{J})_{J \times 1}$ and modelled as independently Gaussian with mean zero and a constant variance, and due to the row standardised nature of $\boldsymbol{W}$, $\sum_{j \in \textrm{net}(i_s)} w_{i_sj} u_{j}$ represents the average (mean) effect that the peers of individual $i$ in spatial unit or group $s$ have on that individual. $w^{*}_{i_s}u^{*}$ is an isolation effect, which is an effect for individuals who don't nominate any friends. This is achieved by setting $w^{*}_{i_s}=1$ if individual $i_s$ nominates no peers and $w^{*}_{i_s}=0$ otherwise, and if $w^{*}_{i_s}=1$ then clearly $\sum_{j\in \textrm{net}(i_{s})}w_{i_{s}j}u_{jr}=0$ as $\textrm{net}(i_{s})$ is the empty set. A conjugate Inverse-Gamma prior is specified for the random effects variance $\sigma_{u}^{2}$, and the corresponding hyperparamaterers ($\alpha_{2}$, $\xi_{2}$) can be chosen by the user.

The exact specification of each of the likelihoods (binomial, Gaussian, and Poisson) are given below:

$$\textrm{Binomial:} ~ Y_{i_s} \sim \textrm{Binomial}(n_{i_s}, \theta_{i_s}) ~ \textrm{and} ~ g(\mu_{i_s}) = \textrm{ln}(\theta_{i_s} / (1 - \theta_{i_s})),$$ $$\textrm{Gaussian:} ~ Y_{i_s} \sim \textrm{N}(\mu_{i_s}, \sigma_{e}^{2}) ~ \textrm{and} ~ g(\mu_{i_s}) = \mu_{i_s},$$ $$\textrm{Poisson:} ~ Y_{i_s} \sim \textrm{Poisson}(\mu_{i_s}) ~ \textrm{and} ~ g(\mu_{i_s}) = \textrm{ln}(\mu_{i_s}).$$

Usage

uniNet(formula, data, trials, family, W, numberOfSamples = 10, burnin = 0, thin = 1, 
seed = 1, trueBeta = NULL, trueURandomEffects = NULL, trueSigmaSquaredU = NULL,
trueSigmaSquaredE = NULL, covarianceBetaPrior = 10^5, a2 = 0.001, b2 = 0.001, 
a3 = 0.001, b3 = 0.001, centerURandomEffects = TRUE)

Value

call: The matched call.
y: The response used.
X: The design matrix used.
standardizedX: The standardized design matrix used.
W: The network matrix used.
samples: The matrix of simulated samples from the posterior distribution of each parameter in the model (excluding random effects).
betaSamples: The matrix of simulated samples from the posterior distribution of $\boldsymbol{\beta}$ parameters in the model.
sigmaSquaredUSamples: The vector of simulated samples from the posterior distribution of $\sigma_{u}^{2}$ in the model.
sigmaSquaredESamples: The vector of simulated samples from the posterior distribution of $\sigma_{e}^{2}$ in the model.
uRandomEffectsSamples: The matrix of simulated samples from the posterior distribution of network random effects $\boldsymbol{u}$ in the model.
acceptanceRates: The acceptance rates of parameters in the model (excluding random effects) from the MCMC sampling scheme .
uRandomEffectsAcceptanceRate: The acceptance rates of network random effects in the model from the MCMC sampling scheme.
timeTaken: The time taken for the model to run.
burnin: The number of MCMC samples to discard as the burn-in period.
thin: The value by which to thin $\texttt{numberOfSamples}$.
DBar: DBar for the model.
posteriorDeviance: The posterior deviance for the model.
posteriorLogLikelihood: The posterior log likelihood for the model.
pd: The number of effective parameters in the model.
DIC: The DIC for the model.

Arguments

formula: A formula for the covariate part of the model using a similar syntax to that used in the lm() function.
data: An optional data.frame containing the variables in the formula.
trials: A vector the same length as the response containing the total number of trials $n_{i_s}$. Only used if $\texttt{family}$=``binomial".
family: The data likelihood model that must be ``gaussian", ``poisson" or ``binomial".
W: A matrix $\boldsymbol{W}$ that encodes the social network structure and whose rows sum to 1.
numberOfSamples: The number of samples to generate pre-thin.
burnin: The number of MCMC samples to discard as the burn-in period.
thin: The value by which to thin $\texttt{numberOfSamples}$.
seed: A seed for the MCMC algorithm.
trueBeta: If available, the true value of $\boldsymbol{\beta}$.
trueURandomEffects: If available, the true value of $\boldsymbol{u}$.
trueSigmaSquaredU: If available, the true value $\sigma_{u}^{2}$.
trueSigmaSquaredE: If available, the true value $\sigma_{e}^{2}$.
covarianceBetaPrior: A scalar prior $\alpha$ for the covariance parameter of the beta prior, such that the covariance is $\alpha\boldsymbol{I}$.
a2: The shape parameter for the Inverse-Gamma distribution relating to the network random effects $\alpha_{2}$.
b2: The scale parameter for the Inverse-Gamma distribution relating to the network random effects $\xi_{2}$.
a3: The shape parameter for the Inverse-Gamma distribution relating to the error terms $\alpha_{3}$. Only used if $\texttt{family}$=``gaussian".
b3: The scale parameter for the Inverse-Gamma distribution relating to the error terms $\xi_{3}$. Only used if $\texttt{family}$=``gaussian".
centerURandomEffects: A choice to center the network random effects after each iteration of the MCMC sampler.

Author

George Gerogiannis

Examples

Run this code

  #################################################
  #### Run the model on simulated data
  #################################################
  #### Load other libraries required
  library(MCMCpack)
  
  #### Set up a network
  observations <- 200
  numberOfMultipleClassifications <- 50
  W <- matrix(rbinom(observations * numberOfMultipleClassifications, 1, 0.05), 
              ncol = numberOfMultipleClassifications)
  numberOfActorsWithNoPeers <- sum(apply(W, 1, function(x) { sum(x) == 0 }))
  peers <- sample(1:numberOfMultipleClassifications, numberOfActorsWithNoPeers,
                  TRUE)
  actorsWithNoPeers <- which(apply(W, 1, function(x) { sum(x) == 0 }))
  for(i in 1:numberOfActorsWithNoPeers) {
    W[actorsWithNoPeers[i], peers[i]] <- 1
  }
  W <- t(apply(W, 1, function(x) { x / sum(x) }))
  
  #### Generate the covariates and response data
  X <- matrix(rnorm(2 * observations), ncol = 2)
  colnames(X) <- c("x1", "x2")
  beta <- c(1, -0.5, 0.5)
  sigmaSquaredU <- 1
  uRandomEffects <- rnorm(numberOfMultipleClassifications, mean = 0, 
                          sd = sqrt(sigmaSquaredU))

  logTheta <- cbind(rep(1, observations), X) %*% beta + W %*% uRandomEffects
  Y <- rpois(n = observations, lambda = exp(logTheta))
  data <- data.frame(cbind(Y, X))
  
  #### Run the model
  formula <- Y ~ x1 + x2
  if (FALSE) model <- uniNet(formula = formula, data = data, family="poisson",
                           W = W, numberOfSamples = 10000, burnin = 10000, 
                           thin = 10, seed = 1)

Run the code above in your browser using DataLab