mvnBvs | R Documentation |
The function can be used to perform variable selection for multivariate normal responses incorporating not only information on the mean model, but also information on the variance-covariance structure of the outcomes. A multivariate prior is specified on the latent binary selection indicators to incorporate the dependence between outcomes into the variable selection procedure.
mvnBvs(Y, lin.pred, data, model = "unstructured", hyperParams, startValues, mcmcParams)
Y |
a data.frame containing q continuous multivariate outcomes from |
lin.pred |
a list containing two formula objects: the first formula specifies the p covariates for which variable selection is to be performed; the second formula specifies the confounders to be adjusted for (but on which variable selection is not to be performed) in the regression analysis. |
data |
a data.frame containing the variables named in the formulas in |
model |
a character that specifies the covariance structure of the model: either "unstructured" or "factor-analytic". |
hyperParams |
a list containing lists or vectors for hyperparameter values in hierarchical models. Components include,
|
startValues |
a numeric vector containing starting values for model parameters: c( |
mcmcParams |
a list containing variables required for MCMC sampling. Components include,
|
mvnBvs
returns an object of class mvnBvs
.
Kyu Ha Lee, Mahlet G. Tadesse, Brent A. Coull
Maintainer: Kyu Ha Lee <klee15239@gmail.com>
Lee, K. H., Tadesse, M. G., Baccarelli, A. A., Schwartz J., and Coull, B. A. (2017),
Multivariate Bayesian variable selection exploiting dependence structure among outcomes:
application to air pollution effects on DNA methylation, Biometrics, Volume 73, Issue 1, pages 232-241.
# loading a data set data(simData_cont) Y <- simData_cont$Y data <- simData_cont$X form1 <- as.formula( ~ cov.1+cov.2) form2 <- as.formula( ~ 1) lin.pred <- list(form1, form2) p <- dim(data)[2] p_adj <- 0 q <- dim(Y)[2] ##################### ## Hyperparameters ## ## Common hyperparameters ## eta = 0.1 v = rep(10, q) omega = rep(log(0.5/(1-0.5)), p-p_adj) common.beta0 <- c(rep(0, q), 10^6) ## Unstructured model ## rho0 <- q + 4 Psi0 <- diag(3, q) US.Sigma <- c(rho0, Psi0) ## Factor-analytic model ## FA.lam <- c(rep(0, q), 10^6) FA.sigSq <- c(2, 1) ## hyperParams <- list(eta=eta, v=v, omega=omega, beta0=common.beta0, US=list(US.Sigma=US.Sigma), FA=list(lambda=FA.lam, sigmaSq=FA.sigSq)) ################### ## MCMC SETTINGS ## ## Setting for the overall run ## numReps <- 50 thin <- 1 burninPerc <- 0.5 ## Tuning parameters for specific updates ## ## - those common to all models mhProp_beta_var <- matrix(0.5, p+p_adj, q) ## ## - those specific to the unstructured model mhrho_prop <- 1000 mhPsi_prop <- diag(1, q) ## ## - those specific to the factor-analytic model mhProp_lambda_var <- 0.5 ## mcmc.US <- list(run=list(numReps=numReps, thin=thin, burninPerc=burninPerc), tuning=list(mhProp_beta_var=mhProp_beta_var, mhrho_prop=mhrho_prop, mhPsi_prop=mhPsi_prop)) ## mcmc.FA <- list(run=list(numReps=numReps, thin=thin, burninPerc=burninPerc), tuning=list(mhProp_beta_var=mhProp_beta_var, mhProp_lambda_var=mhProp_lambda_var)) ##################### ## Starting Values ## ## - those common to all models beta0 <- rep(0, q) B <- matrix(sample(x=c(0.3, 0), size=q, replace = TRUE), p+p_adj, q) gamma <- B gamma[gamma != 0] <- 1 ## ## - those specific to the unstructured model Sigma <- diag(1, q) ## ## - those specific to the factor-analytic model lambda <- rep(0.5, q) sigmaSq <- 1 startValues <- as.vector(c(beta0, B, gamma, Sigma)) #################################### ## Fitting the unstructured model ## #################################### fit.us <- mvnBvs(Y, lin.pred, data, model="unstructured", hyperParams, startValues, mcmcParams=mcmc.US) fit.us summ.fit.us <- summary(fit.us); names(summ.fit.us) summ.fit.us ####################################### ## Fitting the factor-analytic model ## ####################################### startValues <- as.vector(c(beta0, B, gamma, sigmaSq, lambda)) fit.fa <- mvnBvs(Y, lin.pred, data, model="factor-analytic", hyperParams, startValues, mcmcParams=mcmc.FA) fit.fa summ.fit.fa <- summary(fit.fa); names(summ.fit.fa) summ.fit.fa
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.