probe_one: Fitting PaRtitiOned empirical Bayes Ecm (PROBE) algorithm to sparse high-dimensional linear models.

Description

A wrapper function for the one-at-a-time variant of the PROBE algorithm.

Usage

probe_one(Y, X, ep = 0.001, maxit = 10000, Y_test = NULL, X_test = NULL, 
verbose = FALSE, signal = NULL, eta_i = NULL, alpha = 0.05, plot_ind = FALSE, 
order.method = "lasso", adj = 10, delta = 0.4, update_order= NULL, beta_start= NULL,
seed = NULL)

Value

A list including

beta_ast_hat MAP estimates of the regression coefficients (\(\beta^\ast\)),

beta_hat, beta_hat_var MAP estimates of the posterior expectation (beta_hat) and variance (beta_hat_var) of the prior mean (\(\beta\)) of the regression coefficients assuming \(\gamma=1\),

gamma_hat the posterior expectation of the latent \(\gamma\) variables,

sigma2_est MAP estimate of the residual variance,

E_step full results of the final E_step,

count the total number of iterations before convergence.

Arguments

Y: The outcome variable.
X: An n x M matrix of sparse predictors variables.
ep: Value against which to compare convergence criterion (default = 0.001).
maxit: Maximum number of iterations the algorithm will run for (default = 10000).
Y_test: (optional) Test Y data used plotting purposes only (doesn't impact results).
X_test: (optional) Test X data used plotting purposes only (doesn't impact results).
verbose: A logical (true/false) value whether to print algorithm iteration progress and summary quantities (default = FALSE).
signal: (optional) A vector of indicies of the true non-null coefficients. This is used to calculate the true and false discovery rates by iteration for simulated data. Used plotting purposes only (doesn't impact results).
eta_i: (optional) A vector of the true signal. This is used to calculate the MSE by iteration for simulated data. Used plotting purposes only (doesn't impact results).
alpha: (optional) significance level
plot_ind: A logical values (True/False) for whether to output plots on algorithm results and progress (default = FALSE)
order.method: Updating order and initial values of the algorithm. For lasso (default) or ridge, a lasso or a ridge regression model (fit with 10-fold CV) will be fitted and used. The update_order is defined by the absolute values of the coefficient and beta_start is the coefficient values. When using none, update_order and beta_start must be given. random will randomly select the updating order and use very small values for beta_start.
adj: Bandwidth parameter for empirical Bayes E-step. The bandwidth will be equal to adj times Silverman's 'rule of thumb' (default = 10).
delta: Learning rate for iteration t is (1 + t)^(-1 + delta) (default delta = 0.4).
update_order: Manual value for the updating order for when order.method = "none" is used.
beta_start: Manual value for the starting beta coefficients for when order.method = "none" is used.
seed: Seed value to ensure reproducibility when order.method = "lasso", order.method = "ridge", or order.method = "random".

References

McLain, A. C., Zgodic, A., & Bondell, H. (2022). Sparse high-dimensional linear regression with a partitioned empirical Bayes ECM algorithm. arXiv preprint arXiv:2209.08139..

Examples

Run this code

### Example
data(Sim_data)
data(Sim_data_test)
attach(Sim_data)
attach(Sim_data_test)
plot_ind <- TRUE
adj <- 10

# Run the analysis. Y_test and X_test are included for plotting purposes only
full_res <- probe_one( Y = Y, X = X, Y_test = Y_test, order.method = "random",
X_test = X_test, plot_ind = plot_ind, adj = adj)

# Predicting for test data
pred_res <- predict_probe_func(full_res, X = X_test)
sqrt(mean((Y_test - pred_res$Pred)^2))

# Estimate of the residual variance and true value
full_res$sigma2_est
sigma2_tr

# RMSE of estimated beta coefficients
beta_ast_est <- c(full_res$beta_ast_hat)
sqrt(mean((beta_ast_est - beta_tr)^2))

# Posterior expectation of gamma by true
gamma_est <- full_res$E_step$gamma
table(gamma_est > 0.5, beta_tr > 0)
sum(gamma_est)
sum(gamma_est[beta_tr>0])