mr.mash | R Documentation |
Performs multivariate multiple regression with mixture-of-normals prior.
mr.mash(
X,
Y,
S0,
w0 = rep(1/(length(S0)), length(S0)),
V = NULL,
mu1_init = matrix(0, nrow = ncol(X), ncol = ncol(Y)),
tol = 0.0001,
convergence_criterion = c("mu1", "ELBO"),
max_iter = 5000,
update_w0 = TRUE,
update_w0_method = "EM",
w0_penalty = rep(1, length(S0)),
update_w0_max_iter = Inf,
w0_threshold = 0,
compute_ELBO = TRUE,
standardize = TRUE,
verbose = TRUE,
update_V = FALSE,
update_V_method = c("full", "diagonal"),
version = c("Rcpp", "R"),
e = 1e-08,
ca_update_order = c("consecutive", "decreasing_logBF", "increasing_logBF", "random"),
nthreads = as.integer(NA)
)
X |
n x p matrix of covariates. |
Y |
n x r matrix of responses. |
S0 |
List of length K containing the desired r x r prior covariance matrices on the regression coefficients. |
w0 |
K-vector with prior mixture weights, each associated with
the respective covariance matrix in |
V |
r x r residual covariance matrix. |
mu1_init |
p x r matrix of initial estimates of the posterior
mean regression coefficients. These should be on the same scale as
the X provided. If |
tol |
Convergence tolerance. |
convergence_criterion |
Criterion to use for convergence check. |
max_iter |
Maximum number of iterations for the optimization algorithm. |
update_w0 |
If |
update_w0_method |
Method to update prior weights. Only EM is currently supported. |
w0_penalty |
K-vector of penalty terms (>=1) for each mixture component. Default is all components are unpenalized. |
update_w0_max_iter |
Maximum number of iterations for the update of w0. |
w0_threshold |
Drop mixture components with weight less than this value. Components are dropped at each iteration after 15 initial iterations. This is done to prevent from dropping some poetentially important components prematurely. |
compute_ELBO |
If |
standardize |
If |
verbose |
If |
update_V |
if |
update_V_method |
Method to update residual covariance. So far,
"full" and "diagonal" are supported. If |
version |
Whether to use R or C++ code to perform the coordinate ascent updates. |
e |
A small number to add to the diagonal elements of the prior matrices to improve numerical stability of the updates. |
ca_update_order |
The order with which coordinates are updated. So far, "consecutive", "decreasing_logBF", "increasing_logBF", "random" are supported. |
nthreads |
Number of RcppParallel threads to use for the
updates. When |
A mr.mash fit, stored as a list with some or all of the following elements:
mu1 |
p x r matrix of posterior means for the regression coeffcients. |
S1 |
r x r x p array of posterior covariances for the regression coeffcients. |
w1 |
p x K matrix of posterior assignment probabilities to the mixture components. |
V |
r x r residual covariance matrix |
w0 |
K-vector with (updated, if |
.
S0 |
r x r x K array of prior covariance matrices on the regression coefficients |
.
intercept |
r-vector containing posterior mean estimate of the intercept. |
fitted |
n x r matrix of fitted values. |
G |
r x r covariance matrix of fitted values. |
pve |
r-vector of proportion of variance explained by the covariates. |
ELBO |
Evidence Lower Bound (ELBO) at last iteration. |
progress |
A data frame including information regarding convergence criteria at each iteration. |
converged |
|
Y |
n x r matrix of responses at last iteration (only relevant when missing values are present in the input Y). |
###Set seed
set.seed(123)
###Simulate X and Y
##Set parameters
n <- 1000
p <- 100
p_causal <- 20
r <- 5
###Simulate data
out <- simulate_mr_mash_data(n, p, p_causal, r, pve=0.5, B_cor=1,
B_scale=1, X_cor=0, X_scale=1, V_cor=0)
###Split the data in training and test sets
Ytrain <- out$Y[-c(1:200), ]
Xtrain <- out$X[-c(1:200), ]
Ytest <- out$Y[c(1:200), ]
Xtest <- out$X[c(1:200), ]
###Specify the covariance matrices for the mixture-of-normals prior.
univ_sumstats <- compute_univariate_sumstats(Xtrain, Ytrain,
standardize=TRUE, standardize.response=FALSE)
grid <- autoselect.mixsd(univ_sumstats, mult=sqrt(2))^2
S0 <- compute_canonical_covs(ncol(Ytrain), singletons=TRUE,
hetgrid=c(0, 0.25, 0.5, 0.75, 1))
S0 <- expand_covs(S0, grid, zeromat=TRUE)
###Fit mr.mash
fit <- mr.mash(Xtrain, Ytrain, S0, update_V=TRUE)
# Compare the "fitted" values of Y against the true Y in the training set.
plot(fit$fitted,Ytrain,pch = 20,col = "darkblue",xlab = "true",
ylab = "fitted")
abline(a = 0,b = 1,col = "magenta",lty = "dotted")
# Predict the multivariate outcomes in the test set using the fitted model.
Ytest_est <- predict(fit,Xtest)
plot(Ytest_est,Ytest,pch = 20,col = "darkblue",xlab = "true",
ylab = "predicted")
abline(a = 0,b = 1,col = "magenta",lty = "dotted")
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.