BVS4GCR: BVS4GCR

In lb664/BVS4GCR: Bayesian Variable Selection for Gaussian Copula Regression models

BVS4GCR

Description

MCMC implementation of Bayesian Variable Selection for Gaussian Copula Regression models

Usage

BVS4GCR(
  Y,
  X,
  pFixed,
  responseType,
  EVgamma = c(2, 4),
  tau = 1,
  niter,
  burnin,
  thin,
  monitor,
  fullCov = FALSE,
  nonDecomp = FALSE,
  nTrial = NULL,
  negBinomParInit = NULL,
  link = "probit",
  alpha = 0.05,
  probAddDel = 0.9,
  probAdd = 0.5,
  probMix = 0.25,
  geomMean = EVgamma[1],
  graphParNiter = max(16, dim(Y)[2] * dim(Y)[2]),
  std = FALSE,
  seed = 31122021
)

Arguments

Y	(number of samples) times (number of responses) matrix of responses
X	(number of samples) times (number of predictors) matrix of predictors. If the model contains an intercept (for all responses), the first column is a vector of ones, followed by the covariates and the predictors on which Bayesian Variable Selection is performed
pFixed	Number of covariates (including the intercept)
responseType	Character vector specifying the response's type. Response's types currently supported are c("Gaussian", "binary", "ordinal", "binomial", "geometric", "negative binomial")
EVgamma	A priori number of expected predictors associated with each response and its variance. For details, see \insertCiteKohn2001;textualBVS4GCR and \insertCiteAlexopoulos2021;textualBVS4GCR. Default value is EVgamma = c(2, 2) which implies a priori range of important predictors between 0 and 6 for each response
tau	Prior variance of the beta (regression coefficient) prior with the default value set at 1. If std = TRUE, then tau = 1
niter	Number of MCMC iterations (including burn-in)
burnin	Number of iterations to be discarded as burn-in
thin	Store the outcome at every thin-th iteration
monitor	Display the monitored-th iteration
fullCov	Logical parameter to select full or sparse inverse correlation
nonDecomp	Logical, TRUE if a non-decomposable graphical model is selected
nTrial	Number of trials of binomial responses if "binomial" is included in responseType
negBinomParInit	Initial value for the over-dispersion parameter of the negative binomial distribution
link	Specification for the model link function with default value ("probit")
alpha	Level of significance (0.05 default) to select important predictors for each response by using one at-a-time response, univariable LM or GLM method. Used in the construction of the proposal distribution of gamma
probAddDel	Probability (0.9 default) of adding/deleting one predictor to be directly associated in a model with a response during Markov chain Monte Carlo exploration
probAdd	Probability (0.9 default) of adding one predictor to be directly associated in a model with a response during Markov chain Monte Monte Carlo exploration
probMix	Probability (0.25 default) of selecting a geometric proposal distribution that samples the index of the predictor being added/deleted. For details, see \insertCiteAlexopoulos2021;textualBVS4GCR
geomMean	Mean of the geometric proposal distribution that samples the index of the predictor being added/deleted. Default value is the number of expected predictors directly associated with each response and specified in EVgamma. For details, see \insertCiteAlexopoulos2021;textualBVS4GCR
graphParNiter	Number of Markov chain Monte Carlo internal iterations to sample the graphical model between responses. Default value is the max between 16 and the square of the number of responses
std	Logical parameter to standardise the predictors before the analysis. Default value is set at TRUE
seed	Seed used to initialise the Markov chain Monte Carlo algorithm

Details

For details regarding the model and the algorithm, see details in \insertCiteAlexopoulos2021;textualBVS4GCR. Types of responses currently supported: c("Gaussian", "binary", "ordinal", "binomial", "negative binomial"). For ordinal categorical responses it is assumed that the first category is labelled with zero. The maximum number of categories supported is currently 5.

Value

The value returned is a list object list(B, G, R, D, Cutoff, NBPar)

• Β Matrix of the (thinned) samples drawn from the posterior distribution of the regression coefficients

• G 3D array of the (thinned) samples drawn from the posterior distribution of the graphs

• R 3D array of the (thinned) samples drawn from the posterior distribution of the correlation matrix between the responses

• D Matrix of the (thinned) samples drawn from the posterior distribution of the standard deviations of the responses (non-identifiable parameters in the case of discrete responses)

• cutPoint List of the samples drawn from the posterior distribution of the cut-off points for ordinal categorical responses

• NBPar List of the samples drawn from the posterior distribution of the over-dispersion parameter for the negative responses

• postMean List of the posterior means of list(B, G, R, D)

• hyperPar List of the hyper-parameters list(tau) and the parameters of the beta-binomial distribution derived from EVgamma

• samplerPar List of parameters used in the Markov chain Monte Carlo algorithm list(alpha, probAddDel, probAdd, probMix, geomMean, graphParNiter)

• opt List of options used list(std, seed)

References

\insertAllCited

Examples

# Example 1: Toy example with only Gaussian responses and decomposable graphs specification

d <- Sim_GCRdata(m = 4, n = 500, p = 100, pFixed = 2, responseType = rep("Gaussian", 4), 
                 extras = list(rhoY = 0.8, rhoX = 0.7, fixedEffect = rep(1, 2), muEffect = 
                               rep(3, 4), varEffect = rep(0.5, 4), s1Lev = 0.05, s2Lev = 0.95, 
                               sdResponse = rep(1, 4)), 
                 seed = 28061971)

output <- BVS4GCR(d$Y, d$X, pFixed = 2, responseType = rep("Gaussian", 4), EVgamma = c(5, 3^2), 
                  niter = 250, burnin = 50, thin = 4, monitor = 50, seed = 280610971)

# Example 2: Toy example with only Gaussian responses and non-decomposable graphs 

d <- Sim_GCRdata(m = 8, n = 1000, p = 100, pFixed = 1, responseType = rep("Gaussian", 8), 
                 extras = list(rhoY = 0.8, rhoX = 0.7, fixedEffect = 1, muEffect = rep(3, 8), 
                               varEffect = rep(0.5, 8), s1Lev = 0.10, s2Lev = 0.90, 
                               sdResponse = rep(1, 8)), 
                 seed = 28061971)

output <- BVS4GCR(d$Y, d$X, pFixed = 1, responseType = rep("Gaussian", 8), EVgamma = c(5, 3^2), 
                  niter = 250, burnin = 50, thin = 4, monitor = 50, 
                  nonDecomp = TRUE, seed = 280610971)


# Example 3: Combination of Gaussian, binary and ordinal responses (similar to Scenario I in
# Alexopoulos and Bottolo (2021))

d <- Sim_GCRdata(m = 4, n = 50, p = 30, pFixed = 1, responseType = c("Gaussian", "binary", 
                 "ordinal", "ordinal"), 
                 extras = list(rhoY = 0.8, rhoX = 0.7, fixedEffect = 1, muEffect = 
                               c(1, rep(0.5, 3)), varEffect = c(1, rep(0.2, 3)), s1Lev = 0.15, 
                               s2Lev = 0.95, sdResponse = rep(1, 4), nCat = c(3, 4), cutPoint = 
                               cbind(c(0.5, NA, NA), c(0.5, 1.2, 2))), 
                 seed = 28061971)

output <- BVS4GCR(d$Y, d$X, pFixed = 1, responseType = c("Gaussian", "binary", "ordinal", 
                  "ordinal"), EVgamma =  c(5, 3^2), 
                  niter = 250, burnin = 50, thin = 4, monitor = 50, seed = 280610971)


# Example 4: Combination of Gaussian and count responses (Scenario IV in Alexopoulos and Bottolo 
# (2021))

d <- Sim_GCRdata(m = 4, n = 100, p = 100, pFixed = 2, responseType = c("Gaussian", "binomial", 
                 "negative binomial", "negative binomial"), 
                 extras = list(rhoY = 0.8, rhoX = 0.7, fixedEffect = c(-0.5, 0.5), muEffect = 
                               c(1, rep(0.5, 3)), varEffect = c(1, rep(0.2, 3)), s1Lev = 0.05, 
                               s2Lev = 0.95, sdResponse = rep(1, 4), nTrial = 10, negBinomPar = 
                               c(0.5, 0.75)), 
                 seed = 28061971)

output <- BVS4GCR(d$Y, d$X, pFixed = 2, responseType = c("Gaussian", "binomial", 
                  "negative binomial", "negative binomial"), EVgamma = c(5, 3 ^2), 
                  niter = 250, burnin = 50, thin = 4, monitor = 50, seed = 280610971, 
                  nTrial = 10, negBinomParInit = rep(1, 2))

lb664/BVS4GCR documentation built on Feb. 6, 2023, 11:21 p.m.