An R Package for Bayesian Social Network Modelling

Description

Implements a class of univariate and multivariate spatio-network generalised linear mixed models, with inference in a Bayesian setting using Markov chain Monte Carlo (MCMC) simulation. The response variable can be binomial, Gaussian, and Poisson.

Details

Author(s)

George Gerogiannis <g.gerogiannis.1@research.gla.ac.uk>

Examples

  ## See the examples in the function specific help files.

A function that extracts valuable properties from a raw social network.

Description

This function transforms a network, which is a data.frame type in a specified format, in to a resultant n by n adjacency matrix, where a_{ij} = 0 if vertex i and j (i \neq j) are not adjacent i.e. vertex i and j are not the head/tail of an edge e and a_{ij} = 1 if vertex i and j (i \neq j) are adjacent i.e. vertex i and j are the head/tail of an edge e. a_{ij} = 0 when i = j.

Usage

getAdjacencyMatrix(rawNetwork)

Arguments

Value

Author(s)

George Gerogiannis

Examples

rawNetwork = matrix(NA, 4, 3)
rawNetwork = as.data.frame(rawNetwork)
colnames(rawNetwork)[1] = "labels"
rawNetwork[, 1] = c("A", "B", "C", "D")
rawNetwork[, 2] = c(0, "C", "D", 0)
rawNetwork[, 3] = c("B", 0, "A", "C")
getAdjacencyMatrix(rawNetwork)

rawNetwork = matrix(NA, 4, 3)
rawNetwork = as.data.frame(rawNetwork)
colnames(rawNetwork)[2] = "labels"
rawNetwork[, 1] = c(NA, "Charlie", "David", 0)
rawNetwork[, 2] = c("Alistar", "Bob", "Charlie", "David")
rawNetwork[, 3] = c("Bob", NA, "Alistar", "Charlie")
getAdjacencyMatrix(rawNetwork)

rawNetwork = matrix(NA, 4, 3)
rawNetwork = as.data.frame(rawNetwork)
colnames(rawNetwork)[1] = "labels"
rawNetwork[, 1] = c(245, 344, 234, 104)
rawNetwork[, 2] = c(NA, 234, 104, NA)
rawNetwork[, 3] = c(344, 0, 245, 234)
getAdjacencyMatrix(rawNetwork)

rawNetwork = matrix(NA, 4, 3)
rawNetwork = as.data.frame(rawNetwork)
colnames(rawNetwork)[1] = "labels"
rawNetwork[, 1] = c(245, 344, 234, 104)
rawNetwork[, 2] = c(32, 234, 104, 0)
rawNetwork[, 3] = c(344, 20, 245, 234)
getAdjacencyMatrix(rawNetwork)

rawNetwork = matrix(NA, 4, 3)
rawNetwork = as.data.frame(rawNetwork)
colnames(rawNetwork)[1] = "labels"
rawNetwork[, 1] = c("Alistar", "Bob", "Charlie", "David")
rawNetwork[, 2] = c(NA, "Charlie", "David", 0)
rawNetwork[, 3] = c("Bob", "Blaine", "Alistar", "Charlie")
getAdjacencyMatrix(rawNetwork)

rawNetwork = matrix(NA, 4, 3)
rawNetwork = as.data.frame(rawNetwork)
colnames(rawNetwork)[1] = "labels"
rawNetwork[, 1] = c("Alistar", "Bob", "Charlie", "David")
rawNetwork[, 2] = c(0, "Charlie", 0, 0)
rawNetwork[, 3] = c("Bob", "Blaine", "Alistar", 0)
getAdjacencyMatrix(rawNetwork)

rawNetwork = matrix(NA, 4, 3)
rawNetwork = as.data.frame(rawNetwork)
colnames(rawNetwork)[1] = "labels"
rawNetwork[, 1] = c(245, 344, 234, 104)
rawNetwork[, 2] = c(32, 0, 104, 0)
rawNetwork[, 3] = c(34, 0, 245, 234)
getAdjacencyMatrix(rawNetwork)

A function that generates a data.frame that is the membership matrix of the network.

Description

A function that generates a data.frame that is the membership matrix of the network given individual IDs and the alters that they have nominated.

Usage

  getMembershipMatrix(individualID, alters)

Arguments

Value

Author(s)

George Gerogiannis

Examples

  
  individualID = data.frame(c(1, 2, 3))
  alters = data.frame(c(5, 3, 2), c(5, 6, 1))
  getMembershipMatrix(individualID, alters)
  
  individualID = data.frame(c(1, 2, 3))
  alters = data.frame(c(NA, 3, 2), c(NA, NA, 1))
  getMembershipMatrix(individualID, alters)
  
  individualID = data.frame(c(1, 2, 3))
  alters = data.frame(c(NA, 3, NA), c(NA, NA, 1))
  getMembershipMatrix(individualID, alters)
  
  individualID = data.frame(c(1, 2, 3))
  alters = data.frame(c(NA, 3, NA), c(6, NA, 1))
  getMembershipMatrix(individualID, alters)
  

A function that generates a data.frame that stores the number of alters with a given level of a factor an individual has.

Description

This is a function that can be used to generates a data.frame that stores the number of alters with a given level of a factor an individual has.

Usage

getTotalAltersByStatus(individualID, status, alters)

Arguments

Value

Author(s)

George Gerogiannis

Examples

individualID = data.frame(c(1, 2, 3, 4))
status = data.frame(c(10, 20, 30, 20))
alters = data.frame(c(4, 3, 2, 1), c(3, 4, 1, 2), c(2, 1, 4, 3))
totalAltersByStatus = getTotalAltersByStatus(individualID, status, alters)

individualID = data.frame(c(1, 2, 3, 4))
status = data.frame(c("RegularSmoke", "Nonsmoker", "CasualSmoker", "Nonsmoker"))
alters = data.frame(c(4, 3, 2, 1), c(3, 4, 1, 2), c(5, 1, 5, 3))
totalAltersByStatus = getTotalAltersByStatus(individualID, status, alters)

individualID = data.frame(c(1, 2, 3, 4))
status = data.frame(c(NA, "Nonsmoker", "CasualSmoker", "Nonsmoker"))
alters = data.frame(c(4, 3, 2, 1), c(3, 4, 1, 2), c(5, 1, 5, 3))
totalAltersByStatus = getTotalAltersByStatus(individualID, status, alters)

individualID = data.frame(c(10, 20))
status = data.frame(c(NA, "Nonsmoker"))
alters = data.frame(c(NA, 10), c(20, NA))
totalAltersByStatus = getTotalAltersByStatus(individualID, status, alters)

individualID = data.frame(c(NA, 20))
status = data.frame(c("Smoker", "Nonsmoker"))
alters = data.frame(c(NA, 10), c(20, NA))
totalAltersByStatus = getTotalAltersByStatus(individualID, status, alters)

A function that generates samples for a multivariate fixed effects and network model.

Description

This function that generates samples for a multivariate fixed effects and network model, which is given by

Y_{i_sr}|\mu_{i_sr} \sim f(y_{i_sr}| \mu_{i_sr}, \sigma_{er}^{2}) ~~~ i=1,\ldots, N_{s},~s=1,\ldots,S ,~r=1,\ldots,R,

g(\mu_{i_sr}) = \boldsymbol{x}^\top_{i_s} \boldsymbol{\beta}_{r} + \sum_{j\in \textrm{net}(i_s)}w_{i_sj}u_{jr}+ w^{*}_{i_s}u^{*}_{r},

\boldsymbol{\beta}_{r} \sim \textrm{N}(\boldsymbol{0}, \alpha\boldsymbol{I})

\boldsymbol{u}_{j} = (u_{1j},\ldots, u_{Rj}) \sim \textrm{N}(\boldsymbol{0}, \boldsymbol{\Sigma}_{\boldsymbol{u}}),

\boldsymbol{u}^{*} = (u_{1}^*,\ldots, u_{R}^*) \sim \textrm{N}(\boldsymbol{0}, \boldsymbol{\Sigma}_{\boldsymbol{u}}),

\boldsymbol{\Sigma}_{\boldsymbol{u}} \sim \textrm{Inverse-Wishart}(\xi_{\boldsymbol{u}}, \boldsymbol{\Omega}_{\boldsymbol{u}}),

\sigma_{er}^{2} \sim \textrm{Inverse-Gamma}(\alpha_{3}, \xi_{3}).

The covariates for the ith individual in the sth 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 relating to the rth response are denoted by \boldsymbol{\beta}_{r}, 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_{er}^{2}, and the corresponding hyperparamaterers (\alpha_{3}, \xi_{3}) can be chosen by the user.

The R \times 1 vector of random effects for the jth alter is denoted by \boldsymbol{u}_{j} = (u_{j1}, \ldots, u_{jR})_{R \times 1}, while the R \times 1 vector of isolation effects for all R outcomes is denoted by \boldsymbol{u}^{*} = (u_{1}^*,\ldots, u_{R}^*), and both are assigned multivariate Gaussian prior distributions. The unstructured covariance matrix \boldsymbol{\Sigma}_{\boldsymbol{u}} captures the covariance between the R outcomes at the network level, and a conjugate Inverse-Wishart prior is specified for this covariance matrix \boldsymbol{\Sigma}_{\boldsymbol{u}}. The corresponding hyperparamaterers (\xi_{\boldsymbol{u}}, \boldsymbol{\Omega}_{\boldsymbol{u}}) 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_sr} \sim \textrm{Binomial}(n_{i_sr}, \theta_{i_sr}) ~ \textrm{and} ~ g(\mu_{i_sr}) = \textrm{ln}(\theta_{i_sr} / (1 - \theta_{i_sr})),

\textrm{Gaussian:} ~ Y_{i_sr} \sim \textrm{N}(\mu_{i_sr}, \sigma_{er}^{2}) ~ \textrm{and} ~ g(\mu_{i_sr}) = \mu_{i_sr},

\textrm{Poisson:} ~ Y_{i_sr} \sim \textrm{Poisson}(\mu_{i_sr}) ~ \textrm{and} ~ g(\mu_{i_sr}) = \textrm{ln}(\mu_{i_sr}).

Usage

multiNet(formula, data, trials, family, W, numberOfSamples = 10, burnin = 0, 
thin = 1, seed = 1, trueBeta = NULL, trueURandomEffects = NULL, 
trueVarianceCovarianceU = NULL, trueSigmaSquaredE = NULL, 
covarianceBetaPrior = 10^5, xi, omega, a3 = 0.001, b3 = 0.001, 
centerURandomEffects = TRUE)

Arguments

Value

Author(s)

George Gerogiannis

A function that generates samples for a multivariate fixed effects, spatial, and network model.

Description

This function that generates samples for a multivariate fixed effects, spatial, and network model, which is given by

Y_{i_sr}|\mu_{i_sr} \sim f(y_{i_sr}| \mu_{i_sr}, \sigma_{er}^{2}) ~~~ i=1,\ldots, N_{s},~s=1,\ldots,S ,~r=1,\ldots,R,

g(\mu_{i_sr}) = \boldsymbol{x}^\top_{i_s} \boldsymbol{\beta}_{r} + \phi_{sr} + \sum_{j\in \textrm{net}(i_s)}w_{i_sj}u_{jr}+ w^{*}_{i_s}u^{*}_{r},

\boldsymbol{\beta}_{r} \sim \textrm{N}(\boldsymbol{0}, \alpha\boldsymbol{I})

\boldsymbol{\phi}_{r} = (\phi_{1r},\ldots, \phi_{Sr}) \sim \textrm{N}(\boldsymbol{0}, \tau_{r}^{2}(\rho_{r}(\textrm{diag}(\boldsymbol{A1})-\boldsymbol{A})+(1-\rho_{r})\boldsymbol{I})^{-1}),

\boldsymbol{u}_{j} = (u_{1j},\ldots, u_{Rj}) \sim \textrm{N}(\boldsymbol{0}, \boldsymbol{\Sigma}_{\boldsymbol{u}}),

\boldsymbol{u}^{*} = (u_{1}^*,\ldots, u_{R}^*) \sim \textrm{N}(\boldsymbol{0}, \boldsymbol{\Sigma}_{\boldsymbol{u}}),

\tau_{r}^{2} \sim \textrm{Inverse-Gamma}(a_{1}, b_{1}),

\rho_{r} \sim \textrm{Uniform}(0, 1),

\boldsymbol{\Sigma}_{\boldsymbol{u}} \sim \textrm{Inverse-Wishart}(\xi_{\boldsymbol{u}}, \boldsymbol{\Omega}_{\boldsymbol{u}}),

\sigma_{er}^{2} \sim \textrm{Inverse-Gamma}(\alpha_{3}, \xi_{3}).

The covariates for the ith individual in the sth 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 relating to the rth response are denoted by \boldsymbol{\beta}_{r}, 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_{er}^{2}, and the corresponding hyperparamaterers (\alpha_{3}, \xi_{3}) can be chosen by the user.

Spatial correlation in these areal unit level random effects is most often modelled by a conditional autoregressive (CAR) prior distribution. Using this model spatial correlation is induced into the random effects via a non-negative spatial adjacency matrix \boldsymbol{A} = (a_{sl})_{S \times S}, which defines how spatially close the S areal units are to each other. The elements of \boldsymbol{A}_{S \times S} can be binary or non-binary, and the most common specification is that a_{sl} = 1 if a pair of areal units (\mathcal{G}_{s}, \mathcal{G}_{l}) share a common border or are considered neighbours by some other measure, and a_{sl} = 0 otherwise. Note, a_{ss} = 0 for all s. \tau^{2}_{r} measures the variance of these random effects for the rth response, where a conjugate Inverse-Gamma prior is specified for \tau^{2}_{r} and the corresponding hyperparamaterers (a_{1}, b_{1}) can be chosen by the user. \rho_{r} controls the level of spatial autocorrelation. A non-conjugate uniform prior is specified for \rho_{r}.

The R \times 1 vector of random effects for the jth alter is denoted by \boldsymbol{u}_{j} = (u_{j1}, \ldots, u_{jR})_{R \times 1}, while the R \times 1 vector of isolation effects for all R outcomes is denoted by \boldsymbol{u}^{*} = (u_{1}^*,\ldots, u_{R}^*), and both are assigned multivariate Gaussian prior distributions. The unstructured covariance matrix \boldsymbol{\Sigma}_{\boldsymbol{u}} captures the covariance between the R outcomes at the network level, and a conjugate Inverse-Wishart prior is specified for this covariance matrix \boldsymbol{\Sigma}_{\boldsymbol{u}}. The corresponding hyperparamaterers (\xi_{\boldsymbol{u}}, \boldsymbol{\Omega}_{\boldsymbol{u}}) 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_sr} \sim \textrm{Binomial}(n_{i_sr}, \theta_{i_sr}) ~ \textrm{and} ~ g(\mu_{i_sr}) = \textrm{ln}(\theta_{i_sr} / (1 - \theta_{i_sr})),

\textrm{Gaussian:} ~ Y_{i_sr} \sim \textrm{N}(\mu_{i_sr}, \sigma_{er}^{2}) ~ \textrm{and} ~ g(\mu_{i_sr}) = \mu_{i_sr},

\textrm{Poisson:} ~ Y_{i_sr} \sim \textrm{Poisson}(\mu_{i_sr}) ~ \textrm{and} ~ g(\mu_{i_sr}) = \textrm{ln}(\mu_{i_sr}).

Usage

multiNetLeroux(formula, data, trials, family, squareSpatialNeighbourhoodMatrix,
spatialAssignment, W, numberOfSamples = 10, burnin = 0, thin = 1, seed = 1,
trueBeta = NULL, trueSpatialRandomEffects = NULL, trueURandomEffects = NULL, 
trueSpatialTauSquared = NULL, trueSpatialRho = NULL, 
trueVarianceCovarianceU = NULL, trueSigmaSquaredE = NULL, 
covarianceBetaPrior = 10^5, a1 = 0.001, b1 = 0.001, xi, omega, a3 = 0.001, 
b3 = 0.001, centerSpatialRandomEffects = TRUE, centerURandomEffects = TRUE)

Arguments

Value

Author(s)

George Gerogiannis

A function that generates samples for a multivariate fixed effects, grouping, and network model.

Description

This function that generates samples for a multivariate fixed effects, grouping, and network model, which is given by

Y_{i_sr}|\mu_{i_sr} \sim f(y_{i_sr}| \mu_{i_sr}, \sigma_{er}^{2}) ~~~ i=1,\ldots, N_{s},~s=1,\ldots,S ,~r=1,\ldots,R,

g(\mu_{i_sr}) = \boldsymbol{x}^\top_{i_s} \boldsymbol{\beta}_{r} v_{sr} + \sum_{j\in \textrm{net}(i_s)}w_{i_sj}u_{jr}+ w^{*}_{i_s}u^{*}_{r},

\boldsymbol{\beta}_{r} \sim \textrm{N}(\boldsymbol{0}, \alpha\boldsymbol{I})

\boldsymbol{v}_{s} = (v_{s1},\ldots, v_{sR}) \sim \textrm{N}(\boldsymbol{0}, \boldsymbol{\Sigma}_{\boldsymbol{v}})\boldsymbol{v}_{s} = (v_{s1},\ldots, v_{sR}) \sim \textrm{N}(\boldsymbol{0}, \boldsymbol{\Sigma}_{\boldsymbol{v}}),

\boldsymbol{u}_{j} = (u_{1j},\ldots, u_{Rj}) \sim \textrm{N}(\boldsymbol{0}, \boldsymbol{\Sigma}_{\boldsymbol{u}}),

\boldsymbol{u}^{*} = (u_{1}^*,\ldots, u_{R}^*) \sim \textrm{N}(\boldsymbol{0}, \boldsymbol{\Sigma}_{\boldsymbol{u}}),

\boldsymbol{\Sigma}_{\boldsymbol{v}} \sim \textrm{Inverse-Wishart}(\xi_{\boldsymbol{v}}, \boldsymbol{\Omega}_{\boldsymbol{v}}),

\boldsymbol{\Sigma}_{\boldsymbol{u}} \sim \textrm{Inverse-Wishart}(\xi_{\boldsymbol{u}}, \boldsymbol{\Omega}_{\boldsymbol{u}}),

\sigma_{er}^{2} \sim \textrm{Inverse-Gamma}(\alpha_{3}, \xi_{3}).

The covariates for the ith individual in the sth 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 relating to the rth response are denoted by \boldsymbol{\beta}_{r}, 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_{er}^{2}, and the corresponding hyperparamaterers (\alpha_{3}, \xi_{3}) can be chosen by the user.

The R \times 1 vector of random effects for the $s$th group is denoted by \boldsymbol{v}_{s} = (v_{s1}, \ldots, v_{sR})_{R \times 1}, which is assigned a joint Gaussian prior distribution with an unstructured covariance matrix \boldsymbol{\Sigma}_{\boldsymbol{v}} that captures the covariance between the R outcomes. A conjugate Inverse-Wishart prior is specified for the random effects covariance matrix \boldsymbol{\Sigma}_{\boldsymbol{v}}. The corresponding hyperparamaterers (\xi_{\boldsymbol{v}}, \boldsymbol{\Omega}_{\boldsymbol{v}}) can be chosen by the user.

The R \times 1 vector of random effects for the jth alter is denoted by \boldsymbol{u}_{j} = (u_{j1}, \ldots, u_{jR})_{R \times 1}, while the R \times 1 vector of isolation effects for all R outcomes is denoted by \boldsymbol{u}^{*} = (u_{1}^*,\ldots, u_{R}^*), and both are assigned multivariate Gaussian prior distributions. The unstructured covariance matrix \boldsymbol{\Sigma}_{\boldsymbol{u}} captures the covariance between the R outcomes at the network level, and a conjugate Inverse-Wishart prior is specified for this covariance matrix \boldsymbol{\Sigma}_{\boldsymbol{u}}. The corresponding hyperparamaterers (\xi_{\boldsymbol{u}}, \boldsymbol{\Omega}_{\boldsymbol{u}}) 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_sr} \sim \textrm{Binomial}(n_{i_sr}, \theta_{i_sr}) ~ \textrm{and} ~ g(\mu_{i_sr}) = \textrm{ln}(\theta_{i_sr} / (1 - \theta_{i_sr})),

\textrm{Gaussian:} ~ Y_{i_sr} \sim \textrm{N}(\mu_{i_sr}, \sigma_{er}^{2}) ~ \textrm{and} ~ g(\mu_{i_sr}) = \mu_{i_sr},

\textrm{Poisson:} ~ Y_{i_sr} \sim \textrm{Poisson}(\mu_{i_sr}) ~ \textrm{and} ~ g(\mu_{i_sr}) = \textrm{ln}(\mu_{i_sr}).

Usage

multiNetRand(formula, data, trials, family, V, W, numberOfSamples = 10, burnin = 0, 
thin = 1, seed = 1, trueBeta = NULL, trueVRandomEffects = NULL, 
trueURandomEffects = NULL, trueVarianceCovarianceV = NULL, 
trueVarianceCovarianceU = NULL, trueSigmaSquaredE = NULL, 
covarianceBetaPrior = 10^5, xiV, omegaV, xi, omega, a3 = 0.001, 
b3 = 0.001, centerVRandomEffects = TRUE, centerURandomEffects = TRUE)

Arguments

Value

Author(s)

George Gerogiannis

A function that plots visual MCMC diagnostics of the fitted model.

Description

This function takes a netcmc object of samples from the posterior distribution of a parameter(s) and returns a visual convergence diaagnostics in the form of a density plot, trace plot, and ACF plot.

Usage

  ## S3 method for class 'netcmc'
plot(x, ...)

Arguments

Value

Returns a trace plot, density plot and ACF plot for the posterior distribution of a parameter(s) in a netcmc object.

Author(s)

George Gerogiannis

A function that gets a summary of the fitted model.

Description

This function takes a netcmc object and returns a summary of the fitted model. The summary includes, for selected parameters, posterior medians and 95 percent credible intervals, the effective number of independent samples and the Geweke convergence diagnostic in the form of a Z-score.

Usage

  ## S3 method for class 'netcmc'
print(x, ...)

Arguments

Value

Returns a model summary for a netcmc object.

Author(s)

George Gerogiannis

A function that gets a summary of the fitted model.

Description

This function takes a netcmc object and returns a summary of the fitted model. The summary includes, for selected parameters, posterior medians and 95 percent credible intervals, the effective number of independent samples and the Geweke convergence diagnostic in the form of a Z-score.

Usage

  ## S3 method for class 'netcmc'
summary(object, ...)

Arguments

Value

Returns a model summary for a netcmc object.

Author(s)

George Gerogiannis

A function that generates samples for a univariate fixed effects model.

Description

This function generates samples for a univariate fixed effects 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},

\boldsymbol{\beta} \sim \textrm{N}(\boldsymbol{0}, \alpha\boldsymbol{I}),

\sigma_{e}^{2} \sim \textrm{Inverse-Gamma}(\alpha_{3}, \xi_{3}).

The covariates for the ith individual in the sth 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 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

uni(formula, data, trials, family, numberOfSamples = 10, burnin = 0, thin = 1, seed = 1,
trueBeta = NULL, trueSigmaSquaredE = NULL, covarianceBetaPrior = 10^5, 
a3 = 0.001, b3 = 0.001)

Arguments

Value

Author(s)

George Gerogiannis

Examples

  #################################################
  #### Run the model on simulated data
  #################################################
  
  #### Generate the covariates and response data
  observations <- 100
  X <- matrix(rnorm(2 * observations), ncol = 2)
  colnames(X) <- c("x1", "x2")
  beta <- c(2, -2, 2)
  logit <- cbind(rep(1, observations), X) %*% beta
  prob <- exp(logit) / (1 + exp(logit))
  trials <- rep(50, observations)
  Y <- rbinom(n = observations, size = trials, prob = prob)
  data <- data.frame(cbind(Y, X))
  
  #### Run the model
  formula <- Y ~ x1 + x2
  ## Not run: model <- uni(formula = formula, data = data, family="binomial", 
                        trials = trials, numberOfSamples = 10000, 
                        burnin = 10000, thin = 10, seed = 1)
## End(Not run)
                        

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 ith individual in the sth 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)

Arguments

Value

Author(s)

George Gerogiannis

Examples

  #################################################
  #### 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
  ## Not run: model <- uniNet(formula = formula, data = data, family="poisson",
                           W = W, numberOfSamples = 10000, burnin = 10000, 
                           thin = 10, seed = 1)
## End(Not run)

A function that generates samples for a univariate network Leroux model.

Description

This function generates samples for a univariate network Leroux 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} + \phi_{s} + \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}),

\phi_{s} | \boldsymbol{\phi}_{-s} \sim \textrm{N}\bigg(\frac{ \rho \sum_{l = 1}^{S} a_{sl} \phi_{l} }{ \rho \sum_{l = 1}^{S} a_{sl} + 1 - \rho }, \frac{ \tau^{2} }{ \rho \sum_{l = 1}^{S} a_{sl} + 1 - \rho } \bigg),

u_{j} \sim \textrm{N}(0, \sigma_{u}^{2}),

u^{*} \sim \textrm{N}(0, \sigma_{u}^{2}),

\tau^{2} \sim \textrm{Inverse-Gamma}(\alpha_{1}, \xi_{1}),

\rho \sim \textrm{Uniform}(0, 1),

\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 ith individual in the sth 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.

Spatial correlation in these areal unit level random effects is most often modelled by a conditional autoregressive (CAR) prior distribution. Using this model spatial correlation is induced into the random effects via a non-negative spatial adjacency matrix \boldsymbol{A} = (a_{sl})_{S \times S}, which defines how spatially close the S areal units are to each other. The elements of \boldsymbol{A}_{S \times S} can be binary or non-binary, and the most common specification is that a_{sl} = 1 if a pair of areal units (\mathcal{G}_{s}, \mathcal{G}_{l}) share a common border or are considered neighbours by some other measure, and a_{sl} = 0 otherwise. Note, a_{ss} = 0 for all s. \boldsymbol{\phi}_{-s}=(\phi_1,\ldots,\phi_{s-1}, \phi_{s+1},\ldots,\phi_{S}). Here \tau^{2} is a measure of the variance relating to the spatial random effects \boldsymbol{\phi}, while \rho controls the level of spatial autocorrelation, with values close to one and zero representing strong autocorrelation and independence respectively. A non-conjugate uniform prior on the unit interval is specified for the single level of spatial autocorrelation \rho. In contrast, a conjugate Inverse-Gamma prior is specified for the random effects variance \tau^{2}, and corresponding hyperparamaterers (\alpha_{1}, \xi_{1}) 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

uniNetLeroux(formula, data, trials, family,
squareSpatialNeighbourhoodMatrix, spatialAssignment, W, numberOfSamples = 10, 
burnin = 0, thin = 1, seed = 1, trueBeta = NULL, 
trueSpatialRandomEffects = NULL, trueURandomEffects = NULL, 
trueSpatialTauSquared = NULL, trueSpatialRho = NULL, trueSigmaSquaredU = NULL,
trueSigmaSquaredE = NULL, covarianceBetaPrior = 10^5, a1 = 0.001, b1 = 0.001, 
a2 = 0.001, b2 = 0.001, a3 = 0.001, b3 = 0.001, 
centerSpatialRandomEffects = TRUE, centerURandomEffects = TRUE)

Arguments

Value

Author(s)

George Gerogiannis

Examples

  #################################################
  #### 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) }))
  
  #### Set up a spatial structure
  numberOfSpatialAreas <- 100
  factor = sample(1:numberOfSpatialAreas, observations, TRUE)
  spatialAssignment = matrix(NA, ncol = numberOfSpatialAreas, 
                             nrow = observations)
  for(i in 1:length(factor)){
    for(j in 1:numberOfSpatialAreas){
      if(factor[i] == j){
        spatialAssignment[i, j] = 1
      } else {
        spatialAssignment[i, j] = 0
      }
    }
  }
  
  gridAxis = sqrt(numberOfSpatialAreas)
  easting = 1:gridAxis
  northing = 1:gridAxis
  grid = expand.grid(easting, northing)
  numberOfRowsInGrid = nrow(grid)
  distance = as.matrix(dist(grid))
  squareSpatialNeighbourhoodMatrix = array(0, c(numberOfRowsInGrid, 
                                                numberOfRowsInGrid))
  squareSpatialNeighbourhoodMatrix[distance==1] = 1

  #### Generate the covariates and response data
  X <- matrix(rnorm(2 * observations), ncol = 2)
  colnames(X) <- c("x1", "x2")
  beta <- c(2, -2, 2)
  
  spatialRho <- 0.5
  spatialTauSquared <- 2
  spatialPrecisionMatrix = spatialRho * 
    (diag(apply(squareSpatialNeighbourhoodMatrix, 1, sum)) -
     squareSpatialNeighbourhoodMatrix) + (1 - spatialRho) * 
     diag(rep(1, numberOfSpatialAreas))
  spatialCovarianceMatrix = solve(spatialPrecisionMatrix)
  spatialPhi = mvrnorm(n = 1, mu = rep(0, numberOfSpatialAreas), 
                       Sigma = (spatialTauSquared * spatialCovarianceMatrix))
  
  sigmaSquaredU <- 2
  uRandomEffects <- rnorm(numberOfMultipleClassifications, mean = 0, 
                          sd = sqrt(sigmaSquaredU))
  
  logit <- cbind(rep(1, observations), X) %*% beta + 
    spatialAssignment %*% spatialPhi + W %*% uRandomEffects
  prob <- exp(logit) / (1 + exp(logit))
  trials <- rep(50, observations)
  Y <- rbinom(n = observations, size = trials, prob = prob)
  data <- data.frame(cbind(Y, X))
  
  #### Run the model
  formula <- Y ~ x1 + x2
  ## Not run: model <- uniNetLeroux(formula = formula, data = data, 
    family="binomial",  W = W,
    spatialAssignment = spatialAssignment, 
    squareSpatialNeighbourhoodMatrix = squareSpatialNeighbourhoodMatrix,
    trials = trials, numberOfSamples = 10000, 
    burnin = 10000, thin = 10, seed = 1)
## End(Not run)

A function that generates samples for a univariate network group model.

Description

This function generates samples for a univariate network group 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} + v_{s} + \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}),

v_{s} \sim \textrm{N}(0, \tau^{2}),

u_{j} \sim \textrm{N}(0, \sigma_{u}^{2}),

u^{*} \sim \textrm{N}(0, \sigma_{u}^{2}),

\tau^{2} \sim \textrm{Inverse-Gamma}(\alpha_{1}, \xi_{1}),

\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 ith individual in the sth 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 S \times 1 vector of random effects for the groups are collectively denoted by \boldsymbol{v} = (v_{1}, \ldots, v_{S})_{S \times 1}, and each element is assigned an independent zero-mean Gaussian prior distribution with a constant variance \tau^{2}. A conjugate Inverse-Gamma prior is specified for \tau^{2}. The corresponding hyperparamaterers (\alpha_{1}, \xi_{1}) 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

uniNetRand(formula, data, trials, family, groupAssignment, W, numberOfSamples = 10, 
burnin = 0, thin = 1, seed = 1, trueBeta = NULL, 
trueGroupRandomEffects = NULL, trueURandomEffects = NULL, 
trueTauSquared = NULL, trueSigmaSquaredU = NULL, 
trueSigmaSquaredE = NULL, covarianceBetaPrior = 10^5, a1 = 0.001, b1 = 0.001, 
a2 = 0.001, b2 = 0.001, a3 = 0.001, b3 = 0.001, 
centerGroupRandomEffects = TRUE, centerURandomEffects = TRUE)

Arguments

Value

Author(s)

George Gerogiannis

Examples

  #################################################
  #### 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) }))
  
  #### Set up a single level classification
  numberOfSingleClassifications <- 20
  factor = sample(1:numberOfSingleClassifications, observations, TRUE)
  V = matrix(NA, ncol = numberOfSingleClassifications, nrow = observations)
  for(i in 1:length(factor)){
    for(j in 1:numberOfSingleClassifications){
      if(factor[i] == j){
        V[i, j] = 1
      } else {
        V[i, j] = 0
      }
    }
  }
  
  #### 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)
  tauSquared <- 0.5
  vRandomEffects <- rnorm(numberOfSingleClassifications, mean = 0, 
                          sd = sqrt(tauSquared))
  sigmaSquaredU <- 1
  uRandomEffects <- rnorm(numberOfMultipleClassifications, mean = 0, 
                          sd = sqrt(sigmaSquaredU))

  logTheta <- cbind(rep(1, observations), X) %*% beta + V %*% vRandomEffects 
            + W %*% uRandomEffects
  Y <- rpois(n = observations, lambda = exp(logTheta))
  data <- data.frame(cbind(Y, X))
  
  #### Run the model
  formula <- Y ~ x1 + x2
  ## Not run: model <- uniNetRand(formula = formula, data = data, family="poisson", 
                               W = W, groupAssignment = V, 
                               numberOfSamples = 10000, burnin = 10000,
                               thin = 10, seed = 1)
## End(Not run)