Figures_for_the_paper

Introduction

This .Rmd documentation is to reproduce figures in the paper/manuscript.

library(Matrix); library(ggplot2); library(cowplot); library(susieR); library(BGLR);
library(glmnet); library(mr.ash.alpha); library(ncvreg); library(L0Learn); library(varbvs);
standardize = FALSE
source('code/method_wrapper.R')
source('code/sim_wrapper.R')

Popular shrinkage operators

1. Soft thresholding

\[S_{{\rm soft}, \lambda}(b) = {\rm sign}(b) \max \{|b| - \lambda, 0\} = \left\{\begin{array}{ll} 0, & |b| \leq \lambda \\ b - \lambda, & b > \lambda \\ b + \lambda, & b < -\lambda \end{array} \right.\]

2. Hard thresholding

\[S_{{\rm hard}, \lambda}(b) = {\rm sign}(b) (\max \{|b|, \lambda\} - \lambda) = \left\{\begin{array}{ll} 0, & |b| \leq \lambda \\ b, & b > \lambda \\ b, & b < -\lambda \end{array} \right.\]

3. SCAD (Smoothly Clipped Absolute Deviation)

\[S_{{\rm scad}, \lambda, \gamma}(b) = \left\{\begin{array}{ll} S_{{\rm soft}, \lambda}(b), & |b| \leq 2 \lambda \\ \frac{S_{{\rm soft}, \gamma\lambda/(\gamma - 1)}(b)}{1 - (1 / (\gamma - 1))}, & 2\lambda < |b| \leq \gamma \lambda \\ b, & \gamma \lambda < |b| \end{array}\right.\]

4. MCP (Minimax Concave Penalty)

\[S_{{\rm mcp}, \lambda, \gamma}(b) = \left\{\begin{array}{ll} \frac{S_{{\rm soft}, \lambda}(b)}{1 - (1 / (\gamma - 1))}, & |b| \leq \gamma \lambda \\ b, & \gamma \lambda < |b| \end{array}\right.\]

scad = function(b, lambda = 1, gamma = 3) {
  b1 = pmax(abs(b) - lambda, 0) * sign(b) * (abs(b) <= 2 * lambda)
  b2 = ((gamma-1) * b - sign(b) * gamma * lambda) / (gamma-2) * (abs(b) <= gamma * lambda) * (abs(b) > 2 * lambda)
  b3 = b * (abs(b) > gamma * lambda)
  return (b1+b2+b3)
}

mcp = function(b, lambda = 1, gamma = 2) {
  b1 = pmax(abs(b) - lambda, 0) * sign(b) * (b <= gamma * lambda) * gamma / (gamma-1)
  b2 = b * (abs(b) > gamma * lambda)
  return (b1+b2)
}

softt = function(b, lambda = 1) {
  return (pmax(abs(b) - lambda, 0) * sign(b))
}

hardt = function(b, lambda = 1) {
  return (b * (abs(b) >= lambda))
}

MR.ASH shrinkage operator may resemble the popular shrinkage operators

postmean = function(b, pi, sa2 = 0:100, sigma2 = 1){
  phi = outer(b^2, 1 / 2 / (1 + 1 / sa2) / sigma2);
  phi = exp(phi - apply(phi, 1, max))
  
  phi = t(pi * t(phi) / sqrt(1 + sa2));
  phi = phi / rowSums(phi);
  out = c(colSums(t(phi) / (1 + 1 / sa2))) * b
  return (out)
}

sa2 = 0:100; b = seq(0,10,0.001)

pi = exp(-sa2)
pi = pi / sum(pi)
b_softt = postmean(b, pi)

pi = double(101)
pi[1] = 0.8; pi[101] = 0.2
b_mcp = postmean(b, pi)

pi = double(101)
pi[1] = 0.6; pi[101] = 0.4
b_hardt = postmean(b, pi, sa2 = sa2^5)

pi = double(101)
pi[1] = 0.5; pi[2:51] = 0.2/50; pi[101] = 0.3
b_scad = postmean(b, pi, sa2 = sa2)
df = data.frame(b = c(b,b,b,b), sb = c(b_softt, b_hardt, b_scad, b_mcp),
                sb2 = c(softt(b, lambda = 1.43),
                        hardt(b, lambda = 5),
                        scad(b, lambda = 1, gamma = 4),
                        mcp(b, lambda = 2)),
                operator = rep(c("soft","hard","scad","mcp"), each = length(b)))
df$operator = factor(df$operator, levels = c("soft","hard","scad","mcp"))

p1 = ggplot(df) + geom_line(aes(x = b, y = sb, color = operator)) +
  theme_cowplot(font_size = 14) + theme(axis.line = element_blank()) +
  labs(y = "shrinakge of b (S(b))", title = "MR.ASH shrinkage operators") +
  theme(legend.position = "none")
p2 = ggplot(df) + geom_line(aes(x = b, y = sb2, color = operator)) +
  theme_cowplot(font_size = 14) + theme(axis.line = element_blank()) +
  labs(title = "123")  +
  labs(y = "shrinakge of b (S(b))", title = "Well-known shrinkage/thresholding operators")

fig = plot_grid(p1,p2, nrow = 1, rel_widths = c(0.463,0.5))
title = ggdraw() + draw_label("MR.ASH Shrinkage operator may resemble the well-known shrinkage operators", fontface = 'bold', size = 20) 
fig = plot_grid(title, fig, nrow = 2, rel_heights = c(0.06,0.95))
ggsave("figures/figure1_paper.pdf", fig, width = 18, height = 9)

R.hardt = function(b, lambda = 1) {
  out = b^2/2
  out[b > lambda] = lambda^2/2
  exp(-out)
}

R.softt = function(b, lambda = 1) {
  out = lambda^2/2 + (b - lambda) * lambda
  out[b <= lambda] = b[b <= lambda]^2/2
  exp(-out)
}

R.scad = function(b, lambda = 1, gamma = 3) {
  out = lambda^2/2 + (b - lambda) * lambda
  out[b <= lambda] = b[b <= lambda]^2/2
  out[b > 2 * lambda] = lambda^2 * 3/2 + (gamma - 2) * lambda^2/2 * 
    (1 - (gamma - b[b > 2 * lambda] / lambda)^2 / (gamma - 2)^2)
  out[b > gamma * lambda] = lambda^2 * 3/2 + (gamma - 2) * lambda^2/2
  exp(-out)
}

R.mcp = function(b, lambda = 1, gamma = 2) {
  out = lambda^2/2 + (gamma - 1) * lambda^2 / 2 * 
  (1 - (gamma - b / lambda)^2 / (gamma - 1)^2)
  out[b > gamma * lambda] = lambda^2 / 2 + (gamma - 1) * lambda^2 / 2 
  out[b <= lambda] = b[b <= lambda]^2/2
  exp(-out)
}

marginal_shape = function(b, lambda = 1, gamma = 2, operator = "ridge") {
  if (operator == "ridge") {
    return (exp(-lambda * b^2/2))
  } else if (operator == "soft") {
    out = b
    out[b >= 0] = R.softt(b[b >= 0], lambda = lambda)
    out[b < 0] = R.softt(-b[b < 0], lambda = lambda)
    return (out)
  } else if (operator == "hard") {
    out = b
    out[b >= 0] = R.hardt(b[b >= 0], lambda = lambda)
    out[b < 0] = R.hardt(-b[b < 0], lambda = lambda)
    return (out)
  } else if (operator == "scad") {
    out = b
    out[b >= 0] = R.scad(b[b >= 0], lambda = lambda, gamma = gamma)
    out[b < 0] = R.scad(-b[b < 0], lambda = lambda, gamma = gamma)
    return (out)
  } else if (operator == "mcp") {
    out = b
    out[b >= 0] = R.mcp(b[b >= 0], lambda = lambda, gamma = gamma)
    out[b < 0] = R.mcp(-b[b < 0], lambda = lambda, gamma = gamma)
    return (out)
  }
}

b = seq(-5,5,0.001)
df = data.frame(b = c(b,b,b,b,b), sb = c(marginal_shape(b, operator = "ridge"),
                                         marginal_shape(b, operator = "soft"),
                                         marginal_shape(b, operator = "hard"),
                                         marginal_shape(b, operator = "scad"),
                                         marginal_shape(b, operator = "mcp")),
                operator = rep(c("normal","soft","hard","scad","mcp"), each = length(b)))
df$operator = factor(df$operator, levels = c("normal","soft","hard","scad","mcp"))
ggplot(df) + geom_line(aes(x = b, y = sb, color = operator)) +
  theme_cowplot(font_size = 14) + theme(axis.line = element_blank())

Past versions of fig2-1.png

res_df       = readRDS("results/ridge_pve0.5.RDS")
p_list       = c(50,100,200,500,1000,2000)
method_list  = c("Mr.ASH","VarBVS","BayesB","Blasso","SuSiE","E-NET","Lasso","Ridge","SCAD","MCP","L0Learn","Ridge.opt")
col          = gg_color_hue(13)[1:11]
sdat = data.frame()
for (i in 1:6) {
  sdat       = rbind(sdat, data.frame(pred = colMeans(matrix(res_df[[i]]$pred,20,12)),
                                      time = colMeans(matrix(res_df[[i]]$time,20,12)),
                                      p = p_list[i],
                                      fit = method_list))
}
shape       = c(19,17,24,25,9,3,11,4,5,7,8)
sdat$fit    = factor(sdat$fit, levels =  c("Mr.ASH","E-NET","Lasso","Ridge",
                                           "SCAD","MCP","L0Learn",
                                           "VarBVS","BayesB","Blasso","SuSiE",
                                           "Ridge.opt"))
sdat1       = sdat[sdat$fit %in% c("Mr.ASH","E-NET","Lasso","Ridge","SCAD","MCP","L0Learn","Ridge.opt"),]

p1 = ggplot(sdat1) + geom_line(aes(x = p, y = pred, color = fit)) +
  geom_point(aes(x = p, y = pred, color = fit, shape = fit), size = 2.5) +
  theme_cowplot(font_size = 14) +
  scale_x_continuous(trans = "log10", breaks = p_list) +
  labs(y = "predictior error (rmse / sigma)", x = "number of coefficients (p)") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5)) +
  scale_color_manual(values = c(col[1:7],"gray50")) +
  scale_shape_manual(values = c(shape[1:7],15)) +
  scale_y_continuous(trans = "log10", limits = c(1.04,1.46), breaks = c(1.1,1.2,1.3,1.4))
fig_main  = p1
title     = ggdraw() + draw_label("Prediction Error (log-scale)", fontface = 'bold', size = 20) 
subtitle  = ggdraw() + draw_label("Scenario: IndepGauss + Normal, n = 500, p = 50,100,200,500,1000,2000, s = p, pve = 0.5", fontface  = 'bold', size = 18) 
fig       = plot_grid(title,subtitle,fig_main, ncol = 1, rel_heights = c(0.06,0.06,0.95))
fig

Past versions of fig3-1.png

res_df       = readRDS("results/signalshape_pve0.99.RDS")
method_list  = c("Mr.ASH","VarBVS","BayesB","Blasso","SuSiE","E-NET","Lasso","Ridge","SCAD","MCP","L0Learn")
method_level = c("Mr.ASH","E-NET","Lasso","Ridge",
                 "SCAD","MCP","L0Learn",
                 "VarBVS","BayesB","Blasso","SuSiE")
col          = gg_color_hue(13)[1:11]
for (i in 1:6) {
res_df[[i]]$fit   = rep(method_list, each = 20)
res_df[[i]]$fit   = factor(res_df[[i]]$fit, levels =  c("Mr.ASH","E-NET","Lasso","Ridge",
                                            "SCAD","MCP","L0Learn",
                                            "VarBVS","BayesB","Blasso","SuSiE"))
some              = c(1,2,3,5,6,7)
res_df[[i]] = res_df[[i]][res_df[[i]]$fit %in% method_level[some],]
}
pp = list()
signal_name = c("SparseLaplace","SparseT2","SparseT5","SparseNormal","SparseUnif","SparseConst")
for (i in 1:6) {
  d       = res_df[[i]]
  pp[[i]] = my.box2(d, "fit", "pred", cols = col[some], shapes = 1:6) +
  theme(axis.line    = element_blank(),
        axis.text.x  = element_text(angle = 45,hjust = 1),
        legend.position = "none") +
  geom_hline(yintercept = mean(d$pred[d$fit == "Mr.ASH"]), col = col[1],
             linetype = "dotted", size = 1.5) +
  scale_y_continuous(trans = "log10", breaks = c(1,1.2,1.4,1.6,1.8,2.0)) +
  coord_cartesian(ylim = c(0.95,1.8))
  subtitle  = ggdraw() + draw_label(paste(paste("Signal: ",signal_name[i], sep = ""),"", sep = ""),
                                    fontface  = 'bold', size = 18)
  pp[[i]]   = plot_grid(subtitle, pp[[i]], ncol = 1, rel_heights = c(0.06,0.95))
}
fig_main  = plot_grid(pp[[2]],pp[[3]],pp[[1]],pp[[4]],pp[[5]],pp[[6]], nrow = 2, rel_widths = c(0.3,0.3,0.3,0.3))
title     = ggdraw() + draw_label("Prediction Error (log-scale)", fontface = 'bold', size = 20) 
subtitle  = ggdraw() + draw_label("Scenario: LowdimIndepGauss, n = 500, p = 2000, s = 20, pve = 0.99", fontface = 'bold', size = 18)
fig       = plot_grid(title,subtitle,fig_main, ncol = 1, rel_heights = c(0.03,0.04,0.95))
fig

Past versions of fig4-1.png

res_df       = readRDS("results/sparsesignal.RDS")
sdat         = data.frame()
s_list       = c(1,5,20,100,500,2000)
col          = gg_color_hue(13)[1:11]
shape        = c(19,17,24,25,9,3,11,4,5,7,8)
for (i in 1:6) {
  sdat = rbind(sdat, data.frame(pred = colMeans(matrix(res_df[[i]]$pred, 20, 11)),
                                time = colMeans(matrix(res_df[[i]]$time, 20, 11)),
                                fit  = method_list,
                                s    = s_list[i]))
}
sdat$fit = factor(sdat$fit, levels = method_level)

p1 = ggplot(sdat) + geom_line(aes(x = s, y = pred, color = fit)) +
  geom_point(aes(x = s, y = pred, color = fit, shape = fit), size = 2.5) +
  scale_x_continuous(breaks = s_list, trans = "log10") +
  theme_cowplot(font_size = 14) +
  labs(y = "predictior error (rmse / sigma)", x = "number of nonzero coefficients (s)") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5)) +
  scale_color_manual(values = col) +
  scale_shape_manual(values = shape) +
  scale_y_continuous(trans = "log10", breaks = c(1,1.1,1.2,1.3,1.4)) +
  coord_cartesian(ylim = c(1,1.45))
title     = ggdraw() + draw_label("Prediction Error (log-scale)", fontface = 'bold', size = 20) 
subtitle  = ggdraw() + draw_label("Scenario: IndepGauss + SparseNormal, n = 500, p = 2000, s = 1,5,20,100,500,2000, pve = 0.5", fontface  = 'bold', size = 18) 
p0        = ggplot() + geom_blank() + theme_cowplot() + theme(axis.line = element_blank())
fig_main  = plot_grid(p0,p1,p0, nrow = 1, rel_widths = c(0.3,0.8,0.3))
fig       = plot_grid(title,subtitle,fig_main, ncol = 1, rel_heights = c(0.03,0.04,0.95))
fig

Past versions of fig5-1.png

res_df       = readRDS("results/diffpve.RDS")
method_list  = c("Mr.ASH","VarBVS","BayesB","Blasso","SuSiE","E-NET","Lasso","Ridge","SCAD","MCP","L0Learn")
method_level = c("Mr.ASH","E-NET","Lasso","Ridge",
                 "SCAD","MCP","L0Learn",
                 "VarBVS","BayesB","Blasso","SuSiE")
col          = gg_color_hue(13)[1:11]
shape        = c(19,17,24,25,9,3,11,4,5,7,8)
pve_list     = seq(0,0.9,0.1)
sdat         = data.frame()
for (i in 1:10) {
  sdat = rbind(sdat, data.frame(pred = colMeans(matrix(res_df[[i]]$pred, 20, 11)),
                                time = colMeans(matrix(res_df[[i]]$time, 20, 11)),
                                fit  = method_list,
                                pve  = pve_list[i]))
}
sdat$fit = factor(sdat$fit, levels = method_level)
sdat     = sdat[sdat$fit %in% c("Mr.ASH","E-NET","Lasso","Ridge",
                 "SCAD","MCP","L0Learn"),]
sdat$size = 0.5
sdat$size[1:10] = 1.2

p1 = ggplot(sdat) + geom_line(aes(x = pve, y = pred, color = fit), size = sdat$size) +
  geom_point(aes(x = pve, y = pred, color = fit, shape = fit), size = 2.5) +
  scale_x_continuous(breaks = pve_list) +
  theme_cowplot(font_size = 14) +
  labs(y = "predictior error (rmse / sigma)", x = "proportion of variance explained (pve)") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5),
        legend.position = "none") +
  scale_color_manual(values = col) +
  scale_shape_manual(values = shape) +
  scale_y_continuous(trans = "log10", breaks = c(1,1.1,1.2,1.3,1.4)) +
  coord_cartesian(ylim = c(1,1.3))

title     = ggdraw() + draw_label("Prediction Error (log-scale)", fontface = 'bold', size = 20) 
subtitle  = ggdraw() + draw_label("Scenario: IndepGauss + PointNormal, n = 500, p = 2000, s = 20, pve = 0 - 0.9", fontface  = 'bold', size = 14) 
p0        = ggplot() + geom_blank() + theme_cowplot() + theme(axis.line = element_blank())
fig_main  = plot_grid(p0,p1,p0, nrow = 1, rel_widths = c(0.3,0.8,0.3))
fig       = plot_grid(title,subtitle,fig_main, ncol = 1, rel_heights = c(0.03,0.04,0.95))
fig

Past versions of fig6-1.png

res_df       = readRDS("results/highdim_pve0.5.RDS")
method_list  = c("Mr.ASH","VarBVS","BayesB","Blasso","SuSiE","E-NET","Lasso","SCAD","MCP","L0Learn")
method_level = c("Mr.ASH","E-NET","Lasso",
                 "SCAD","MCP","L0Learn",
                 "VarBVS","BayesB","Blasso","SuSiE")
col          = gg_color_hue(13)[1:11][-4]
#shape        = c(19,4,5,7,8,17,24,25,9,3,11)[-8]
shape        = c(19,17,24,25,9,3,11,4,5,7,8)[-4]
p_list       = c(50,500,5000,50000)
sdat         = data.frame()
for (i in 1:4) {
  sdat = rbind(sdat, data.frame(pred = colMeans(matrix(res_df[[i]]$pred, 20, 10)),
                                time = colMeans(matrix(res_df[[i]]$time, 20, 10)),
                                fit  = method_list,
                                p = p_list[i]))
}
sdat$fit = factor(sdat$fit, levels = method_level)
sdat     = sdat[sdat$fit %in% c("Mr.ASH","E-NET","Lasso","Ridge",
                 "SCAD","MCP","L0Learn"),]
sdat$size = 0.5
sdat$size[1:4] = 1.2

p2 = ggplot(sdat) + geom_line(aes(x = p, y = pred, color = fit), size = sdat$size) +
  geom_point(aes(x = p, y = pred, color = fit, shape = fit), size = 2.5) +
  theme_cowplot(font_size = 14) +
  scale_x_continuous(trans = "log10", breaks = p_list) +
  labs(y = "predictior error (rmse / sigma)", x = "number of coefficients (p)") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5)) +
  scale_color_manual(values = col) +
  scale_shape_manual(values = shape) +
  scale_y_continuous(trans = "log10", breaks = c(1,1.1,1.2,1.3,1.4)) +
  coord_cartesian(ylim = c(1,1.3))

title     = ggdraw() + draw_label("Prediction Error (log-scale)", fontface = 'bold', size = 20) 
subtitle2 = ggdraw() + draw_label("Scenario: IndepGauss + SparseNormal, n = 500, p = 50, 500, 5000, 50000, pve = 0.5", fontface  = 'bold', size = 14) 
p0        = ggplot() + geom_blank() + theme_cowplot() + theme(axis.line = element_blank())
fig_main  = plot_grid(p0,p2,p0, nrow = 1, rel_widths = c(0.3,0.6,0.3))
fig       = plot_grid(title,subtitle2,fig_main, ncol = 1, rel_heights = c(0.06,0.06,0.95))
fig

Past versions of fig7-1.png

fig       = plot_grid(plot_grid(subtitle,p1, ncol = 1, rel_heights = c(0.05,0.95)),
                      plot_grid(subtitle2,p2, ncol = 1, rel_heights = c(0.05,0.95)), nrow = 1, rel_widths = c(0.5,0.55))
fig       = plot_grid(title, fig, ncol = 1, rel_heights = c(0.05,0.95))
ggsave("figures/figure6_paper.pdf", fig, width = 18, height = 9)

res_df       = readRDS("results/highdimdiffp2.RDS")
method_list  = c("Mr.ASH","E-NET","Lasso","Ridge","SCAD","MCP","L0Learn")
col          = gg_color_hue(13)[1:11]
shape        = c(19,17,24,25,9,3,11,4,5,7,8)
p_list       = c(50,500,5000,50000)
sdat         = data.frame()
for (i in 1:4) {
  res_df[[i]] = res_df[[i]][1:140,]
  sdat = rbind(sdat, data.frame(pred = colMeans(matrix(res_df[[i]]$pred, 20, 7)),
                                time = colMeans(matrix(res_df[[i]]$time, 20, 7)),
                                fit  = method_list,
                                p    = p_list[i]))
}
sdat$fit = factor(sdat$fit, levels = method_list)

p1 = ggplot(sdat) + geom_line(aes(x = p, y = time, color = fit)) +
  geom_point(aes(x = p, y = time, color = fit, shape = fit), size = 2.5) +
  theme_cowplot(font_size = 14) +
  scale_x_continuous(trans = "log10", breaks = p_list) +
  labs(y = "computation time (seconds)", x = "number of coefficients (p)") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5)) +
  scale_color_manual(values = col) +
  scale_shape_manual(values = shape) +
  scale_y_continuous(trans = "log10")

title     = ggdraw() + draw_label("Computation time (log-scale)", fontface = 'bold', size = 20) 
subtitle  = ggdraw() + draw_label("Scenario: IndepGauss + SparseNormal, n = 500, p = 50, 500, 5000, 50000, s = 20, pve = 0.5", fontface  = 'bold', size = 18) 
p0        = ggplot() + geom_blank() + theme_cowplot() + theme(axis.line = element_blank())
fig_main  = plot_grid(p0,p1,p0, nrow = 1, rel_widths = c(0.3,0.8,0.3))
fig       = plot_grid(title,subtitle,fig_main, ncol = 1, rel_heights = c(0.03,0.04,0.95))
fig

Past versions of fig71-1.png

res_df = readRDS("results/Ridge_for_paper.RDS")
out = matrix(0,9,12)
lower = matrix(0,9,12)
upper = matrix(0,9,12)
for (i in 1:9) {
  out[i,] = colMeans(matrix(res_df[[i]]$pred, 20, 12))
  lower[i,] = apply(matrix(res_df[[i]]$pred, 20, 12), 2, function(x) quantile(x, probs = 0.1))
  upper[i,] = apply(matrix(res_df[[i]]$pred, 20, 12), 2, function(x) quantile(x, probs = 0.9))
}
colnames(out) = c("Mr.ASH","VarBVS","BayesB","Blasso","SuSiE","E-NET","Lasso","Ridge","SCAD","MCP","L0Learn","Mr.ASH.opt")
ind = 1:12
out = out[,ind]
lower = lower[,ind]
upper = upper[,ind]
s_range = c(1,2,5,10,20,50,100,200,500)
col   = gg_color_hue(13)[1:12]
shape = c(19,17,24,25,9,3,11,4,5,7,8,1)
df = data.frame(s = rep(s_range, length(ind)), pred = c(out), fit = rep(colnames(out), each = 9),
                lower = c(lower), upper = c(upper))
df$fit = factor(df$fit, levels =  c("Mr.ASH","E-NET","Lasso","Ridge",
                                    "SCAD","MCP","L0Learn",
                                    "VarBVS","BayesB","Blasso","SuSiE","Mr.ASH.opt"))
df$size = 0.5
df$size[1:9] = 1.2
df = df[df$fit %in% c("Mr.ASH","E-NET","Lasso","Ridge",
                                    "SCAD","MCP","L0Learn"),]
p1 = ggplot(df) + geom_line(aes(x = s, y = pred, color = fit), size = df$size) +
  geom_point(aes(x = s, y = pred, color = fit, shape = fit), size = 3) +
  theme_cowplot(font_size = 14) +
  scale_x_continuous(trans = "log10", breaks = c(1,2,5,10,20,50,100,200,500)) +
  labs(y = "predictior error (rmse / sigma)", x = "number of nonzero coefficients (s)") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5),
        legend.position = "none") +
  scale_color_manual(values = col) +
  scale_shape_manual(values = shape) +
  scale_y_continuous(trans = "log10", breaks = c(1,1.1,1.2,1.3,1.4)) +
  coord_cartesian(ylim = c(1,1.46), xlim = c(1,600))
fig_main = p1
title     = ggdraw() + draw_label("Prediction Error (log-scale)", fontface = 'bold', size = 20) 
subtitle  = ggdraw() + draw_label("Scenario: IndepGauss + PointNormal, n = 1000, p = 500, s = 1-500, pve = 0.5", fontface  = 'bold', size = 14) 
fig1      = plot_grid(subtitle,fig_main, ncol = 1, rel_heights = c(0.05,0.95))

res_df       = readRDS("results/signalshape25.RDS")
method_list  = c("Mr.ASH","E-NET","Lasso","Ridge","SCAD","MCP","L0Learn")
method_level = c("Mr.ASH","E-NET","Lasso","Ridge","SCAD","MCP","L0Learn")
col          = gg_color_hue(13)[1:11]
x_range      = 2^(c(1,2,4,5,3,6,7,8) - 1)
df           = data.frame()
for (i in 1:8) {
res_df[[i]]       = res_df[[i]][1:140,]
res_df[[i]]$fit   = rep(method_list, each = 20)
res_df[[i]]$fit   = factor(res_df[[i]]$fit, levels = method_level)
df                = rbind(df, data.frame(pred = c(colMeans(matrix(res_df[[i]]$pred, 20, 7))),
                                         df   = x_range[i],
                                         fit  = method_list))
}
df$size = 0.5
df$size[1:8] = 1.2
df$fit = factor(df$fit, levels = method_level)
col   = gg_color_hue(13)[1:7]
shape = c(19,17,24,25,9,3,11,4,5,7,8,1)[1:7]
p1    = ggplot(df) + geom_line(aes(x = df, y = pred, color = fit), size = df$size) +
  geom_point(aes(x = df, y = pred, color = fit, shape = fit), size = 2.5) +
  theme_cowplot(font_size = 14) +
  scale_x_continuous(trans = "log10", breaks = c(1,2,4,8,16,32,64,128),
                     labels = c("t (df=1)","t (df=2)","Laplace","t (df=4)","t (df=8)",
                                "Normal","Uniform","Constant")) +
  labs(y = "predictior error (rmse / sigma)", x = "distribution of nonzero coefficients") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5),
        axis.text.x  = element_text(angle = 45,hjust = 1)) +
  scale_color_manual(values = col) +
  scale_shape_manual(values = shape) +
  scale_y_continuous(trans = "log10", breaks = c(1,1.1,1.2,1.3,1.4)) +
  coord_cartesian(ylim = c(1,sqrt(2)))
fig_main = p1
subtitle  = ggdraw() + draw_label("Scenario: IndepGauss + Spike-and-Slab, n = 1000, p = 500, s = 200, pve = 0.5", fontface  = 'bold', size = 14) 
fig2      = plot_grid(subtitle,fig_main, ncol = 1, rel_heights = c(0.05,0.95))
fig       = plot_grid(title, plot_grid(fig1, fig2, nrow = 1, rel_widths = c(0.5,0.55)),
                      ncol = 1, rel_heights = c(0.05,0.95))
ggsave("figures/figure4_paper.pdf", fig, width = 18, height = 9)
fig

Past versions of fig9-1.png

res_df       = readRDS("results/newhighdim.RDS")
method_list  = c("Mr.ASH","E-NET","Lasso","Ridge","SCAD","MCP","SCAD2","MCP2","L0Learn")
df           = data.frame()
p_range      = c(10,20,50,100,200,500,1000,2000,5000,10000)
for (i in 1:8) {
res_df[[i]]$fit   = rep(method_list, each = 20)
res_df[[i]]$fit   = factor(res_df[[i]]$fit, levels = method_list)
df                = rbind(df, data.frame(pred = c(colMeans(matrix(res_df[[i]]$pred, 20, 9))),
                                         time2 = apply(matrix(res_df[[i]]$time, 20, 9), 2, median),
                                         time = c(colMeans(matrix(res_df[[i]]$time, 20, 9))),
                                         p    = p_range[i],
                                         fit  = method_list))
}
df$fit = factor(df$fit, levels = method_list)
df = df[df$fit %in% c("Mr.ASH","E-NET","Lasso","Ridge","SCAD","MCP","L0Learn"),]
df = df[df$p %in% c(10,20,50,200,1000,2000),]
col   = gg_color_hue(13)[1:7]
shape = c(19,17,24,25,9,3,11,4,5,7,8,1,2,6)[1:7]
p1    = ggplot(df) + geom_line(aes(x = p, y = pred, color = fit)) +
  geom_point(aes(x = p, y = pred, color = fit, shape = fit), size = 2.5) +
  theme_cowplot(font_size = 14) +
  scale_x_continuous(trans = "log10", breaks = c(10,20,50,200,1000,2000)) +
  labs(y = "predictior error (rmse / sigma)", x = "number of nonzero coefficients (s)") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5),
        legend.position = "none") +
  scale_color_manual(values = col) +
  scale_shape_manual(values = shape) +
  scale_y_continuous(trans = "log10", breaks = c(1,1.1,1.2,1.3,1.4)) +
  coord_cartesian(ylim = c(1,1.4))
fig_main = p1
subtitle  = ggdraw() + draw_label("Prediction Error (log-scale)", fontface = 'bold', size = 14) 
title     = ggdraw() + draw_label("Scenario: IndepGauss + PointNormal, n = 500, p = 10-2000, s = 10, pve = 0.5", fontface  = 'bold', size = 20) 
fig1      = plot_grid(subtitle,fig_main, ncol = 1, rel_heights = c(0.05,0.95))

res_df       = readRDS("results/newhighdim.RDS")
method_list  = c("Mr.ASH","E-NET","Lasso","Ridge","SCAD","MCP","SCAD2","MCP2","L0Learn")
df           = data.frame()
p_range      = c(10,20,50,100,200,500,1000,2000,5000,10000)
for (i in 1:8) {
res_df[[i]]$fit   = rep(method_list, each = 20)
res_df[[i]]$fit   = factor(res_df[[i]]$fit, levels = method_list)
df                = rbind(df, data.frame(pred = c(colMeans(matrix(res_df[[i]]$pred, 20, 9))),
                                         time2 = apply(matrix(res_df[[i]]$time, 20, 9), 2, median),
                                         time = c(colMeans(matrix(res_df[[i]]$time, 20, 9))),
                                         p    = p_range[i],
                                         fit  = method_list))
}
df$fit = factor(df$fit, levels = method_list)
df = df[df$fit %in% c("Mr.ASH","E-NET","Lasso","Ridge","SCAD","MCP","L0Learn"),]
df = df[df$p %in% c(10,20,50,100,200,1000,2000),]
col   = gg_color_hue(13)[1:7]
shape = c(19,17,24,25,9,3,11,4,5,7,8,1,2,6)[1:7]
p1    = ggplot(df) + geom_line(aes(x = p, y = time, color = fit)) +
  geom_point(aes(x = p, y = time, color = fit, shape = fit), size = 2.5) +
  theme_cowplot(font_size = 14) +
  scale_x_continuous(trans = "log10", breaks = p_range) +
  labs(y = "computation time (seconds)", x = "number of coefficients (p)") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5)) +
  scale_color_manual(values = col) +
  scale_shape_manual(values = shape) +
  scale_y_continuous(trans = "log10", breaks = c(0.03,0.1,0.3,1,3,10,30,100))
fig_main = p1
subtitle  = ggdraw() + draw_label("Computation time (log-scale)", fontface = 'bold', size = 14) 
fig2      = plot_grid(subtitle,fig_main, ncol = 1, rel_heights = c(0.05,0.95))
fig       = plot_grid(title, plot_grid(fig1, fig2, nrow = 1, rel_widths = c(0.5,0.59)),
                      ncol = 1, rel_heights = c(0.05,0.95))
#ggsave("figures/figure5_paper.pdf", fig, width = 18, height = 8)
fig

Past versions of fig11-1.png

res_df       = readRDS("results/highdimdiffnpnew.RDS")
method_list  = c("Mr.ASH","E-NET","Lasso","Ridge","SCAD","MCP","SCAD2","MCP2","L0Learn")
df           = data.frame()
p_range      = c(20,50,100,200,500,1000,2000,5000,10000,20000)
for (i in 1:10) {
res_df[[i]]$fit   = rep(method_list, each = 20)
res_df[[i]]$fit   = factor(res_df[[i]]$fit, levels = method_list)
df                = rbind(df, data.frame(pred = c(colMeans(matrix(res_df[[i]]$pred, 20, 9))),
                                         time2 = apply(matrix(res_df[[i]]$time, 20, 9), 2, median),
                                         time = c(colMeans(matrix(res_df[[i]]$time, 20, 9))),
                                         p    = p_range[i],
                                         fit  = method_list))
}
df$fit = factor(df$fit, levels = method_list)
#df = df[df$fit %in% c("Mr.ASH","E-NET","Lasso","Ridge","SCAD","MCP","L0Learn"),]
df$size = 0.5
df$size[1:5] = 1.2; df$size[46:50] = 1.2;
col   = gg_color_hue(13)[1:9]
shape = c(19,17,24,25,9,3,11,4,5,7,8,1,2,6)[1:9]
p1    = ggplot(df[df$p <= 500, ]) + geom_line(aes(x = p, y = time, color = fit), size = df$size[df$p <= 500]) +
  geom_point(aes(x = p, y = time, color = fit, shape = fit), size = 2.5) +
  theme_cowplot(font_size = 14) +
  scale_x_continuous(trans = "log10", breaks = p_range) +
  labs(y = "compmutation time (seconds)", x = "number of coefficients (p)") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5),
        legend.position = "none") +
  scale_color_manual(values = col) +
  scale_shape_manual(values = shape) +
  scale_y_continuous(trans = "log10")
fig_main = p1
subtitle  = ggdraw() + draw_label("Computation time (log-scale) for p <= n", fontface = 'bold', size = 14) 
title     = ggdraw() + draw_label("Scenario: IndepGauss + PointNormal, n = 500, p = 20-20000, s = 20, pve = 0.5", fontface  = 'bold', size = 20) 
fig1      = plot_grid(subtitle,fig_main, ncol = 1, rel_heights = c(0.05,0.95))

p2    = ggplot(df[df$p > 500, ]) + geom_line(aes(x = p, y = time, color = fit), size = df$size[df$p > 500]) +
  geom_point(aes(x = p, y = time, color = fit, shape = fit), size = 2.5) +
  theme_cowplot(font_size = 14) +
  scale_x_continuous(trans = "log10", breaks = p_range) +
  labs(y = "compmutation time (seconds)", x = "number of coefficients (p)") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5)) +
  scale_color_manual(values = col) +
  scale_shape_manual(values = shape) +
  scale_y_continuous(trans = "log10")
fig_main = p2
subtitle  = ggdraw() + draw_label("Computation time (log-scale) for p > n", fontface = 'bold', size = 14) 
fig2      = plot_grid(subtitle,fig_main, ncol = 1, rel_heights = c(0.05,0.95))
fig13     = plot_grid(title, plot_grid(fig1, fig2, nrow = 1, rel_widths = c(0.5,0.55)),
                      ncol = 1, rel_heights = c(0.05,0.95))
ggsave("figures/figure5_paper.pdf", fig13, width = 18, height = 9)
fig13

res_df       = readRDS("results/highdimdiffnpnew.RDS")
method_list  = c("Mr.ASH","E-NET","Lasso","Ridge","SCAD","MCP","SCAD2","MCP2","L0Learn")
df           = data.frame()
p_range      = c(20,50,100,200,500,1000,2000,5000,10000,20000)
for (i in 1:10) {
res_df[[i]]$fit   = rep(method_list, each = 20)
res_df[[i]]$fit   = factor(res_df[[i]]$fit, levels = method_list)
df                = rbind(df, data.frame(pred = c(colMeans(matrix(res_df[[i]]$pred, 20, 9))),
                                         time2 = apply(matrix(res_df[[i]]$time, 20, 9), 2, median),
                                         time = c(colMeans(matrix(res_df[[i]]$time, 20, 9))),
                                         p    = p_range[i],
                                         fit  = method_list))
}
df$fit = factor(df$fit, levels = method_list)
#df = df[df$fit %in% c("Mr.ASH","E-NET","Lasso","Ridge","SCAD","MCP","L0Learn"),]
df$size = 0.5
df$size[1:7] = 1.2
df = df[df$p %in% c(20,50,200,500,1000,2000,10000,20000),]
col   = gg_color_hue(13)[1:9]
shape = c(19,17,24,25,9,3,11,4,5,7,8,1,2,6)[1:9]
p1    = ggplot(df) + geom_line(aes(x = p, y = pred, color = fit), size = df$size) +
  geom_point(aes(x = p, y = pred, color = fit, shape = fit), size = 2.5) +
  theme_cowplot(font_size = 14) +
  scale_x_continuous(trans = "log10", breaks = c(20,50,200,500,1000,2000,10000,20000)) +
  labs(y = "predictior error (rmse / sigma)", x = "number of coefficients (p)") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5)) +
  scale_color_manual(values = col) +
  scale_shape_manual(values = shape) +
  scale_y_continuous(trans = "log10")
fig_main = p1
#title     = ggdraw() + draw_label("Prediction Error (log-scale)", fontface = 'bold', size = 20) 
subtitle  = ggdraw() + draw_label("Scenario: IndepGauss + PointNormal, n = 500, p = 20-20000, s = 20, pve = 0.5", fontface  = 'bold', size = 20) 
fig14     = plot_grid(subtitle,fig_main, ncol = 1, rel_heights = c(0.05,0.95))

p2    = ggplot(df[df$p > 500, ]) + geom_line(aes(x = p, y = pred, color = fit), size = df$size[df$p > 500]) +
  geom_point(aes(x = p, y = pred, color = fit, shape = fit), size = 2.5) +
  theme_cowplot(font_size = 14) +
  scale_x_continuous(trans = "log10", breaks = p_range) +
  labs(y = "predictior error (rmse / sigma)", x = "number of coefficients (p)") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5)) +
  scale_color_manual(values = col) +
  scale_shape_manual(values = shape) +
  scale_y_continuous(trans = "log10")
fig_main = p2
subtitle  = ggdraw() + draw_label("Computation time (log-scale) for p > n", fontface = 'bold', size = 14) 
fig2      = plot_grid(subtitle,fig_main, ncol = 1, rel_heights = c(0.05,0.95))
fig       = plot_grid(title, plot_grid(fig1, fig2, nrow = 1, rel_widths = c(0.5,0.55)),
                      ncol = 1, rel_heights = c(0.05,0.95))
ggsave("figures/figure7_paper.pdf", fig, width = 18, height = 9)
fig

res_df = readRDS("results/newridge10.RDS")
out = matrix(0,9,15)
lower = matrix(0,9,15)
upper = matrix(0,9,15)
for (i in 1:9) {
  out[i,] = colMeans(matrix(res_df[[i]]$pred, 20, 15))
  lower[i,] = apply(matrix(res_df[[i]]$pred, 20, 15), 2, function(x) quantile(x, probs = 0.1))
  upper[i,] = apply(matrix(res_df[[i]]$pred, 20, 15), 2, function(x) quantile(x, probs = 0.9))
}
out = out[,c(1:11,15)]
colnames(out) = c("Mr.ASH","VarBVS","BayesB","Blasso","SuSiE","E-NET","Lasso","Ridge","SCAD","MCP","L0Learn","Mr.ASH.opt")
ind = c(6,7,8,9,10,11)
out = out[,ind]
lower = lower[,ind]
upper = upper[,ind]

col   = gg_color_hue(7)[1:7]
shape = c(19,17,24,25,9,3,11,4,5,7,8,1)[1:7]
df = data.frame(s = rep(s_range, length(ind)), pred = c(out), fit = rep(colnames(out), each = 9),
                lower = c(lower), upper = c(upper))
df$fit = factor(df$fit, levels =  c("Mr.ASH","E-NET","Lasso","Ridge",
                                    "SCAD","MCP","L0Learn",
                                    "VarBVS","BayesB","Blasso","SuSiE","Mr.ASH.opt"))
df = df[df$s %in% c(1,5,50,500), ]
p1 = ggplot(df) + geom_line(aes(x = s, y = pred, color = fit), size = 1) +
  geom_point(aes(x = s, y = pred, color = fit, shape = fit), size = 3) +
  theme_cowplot(font_size = 14) +
  scale_x_continuous(trans = "log10", breaks = c(1,5,50,500)) +
  labs(y = "predictior error (rmse / sigma)", x = "number of nonzero coefficients (s)") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5)) +
  scale_color_manual(values = col) +
  scale_shape_manual(values = shape) +
  scale_y_continuous(trans = "log10", breaks = c(1,1.1,1.2,1.3,1.4)) +
  coord_cartesian(ylim = c(1,1.46), xlim = c(1,600)) +
  geom_point(aes(x = 500, y = 1.1942835), size = 3, shape = 4) +
  geom_text(aes(x = 500, y = 1.1942835, label="Oracle "), vjust = 0.3, hjust = -0.2) +
  #geom_text(aes(x = 500, y = df$pred[4], label="Mr.ASH", col = "Mr.ASH"), vjust = 0.2, hjust = -0.2) +
  geom_text(aes(x = 500, y = df$pred[4], label="E-NET ", col = "E-NET"), vjust = -0.2, hjust = -0.2) +
  geom_text(aes(x = 500, y = df$pred[12], label="Ridge  ", col = "Ridge"), vjust = 1, hjust = -0.22) +
  geom_text(aes(x = 500, y = df$pred[8], label="Lasso  ", col = "Lasso"), vjust = 0.2, hjust = -0.22) +
  geom_text(aes(x = 500, y = df$pred[16], label="SCAD  ", col = "SCAD"), vjust = 0.2, hjust = -0.22) +
  geom_text(aes(x = 500, y = df$pred[20], label="MCP   ", col = "MCP"), vjust = 0.2, hjust = -0.22) +
  geom_text(aes(x = 500, y = df$pred[24], label="L0Learn", col = "L0Learn"), vjust = 0.2, hjust = -0.22)
fig_main = p1
title     = ggdraw() + draw_label("Prediction Error (log-scale)", fontface = 'bold', size = 20) 
subtitle  = ggdraw() + draw_label("Scenario: IndepGauss + PointNormal, n = 1000, p = 500, s = 1-500, pve = 0.5", fontface  = 'bold', size = 18) 
fig       = plot_grid(title,subtitle,fig_main, ncol = 1, rel_heights = c(0.05,0.05,0.95))
ggsave("PointNormal_Sparsity2_without_mrash.pdf", fig, width = 14, height = 8)
fig

Past versions of fig100-1.png

fig_main  = fig_main + 
    geom_ribbon(aes(x = s, ymin = lower, ymax = upper, color = fit, fill = fit),
                linetype = 2, size = 0, alpha = 0.1, colour = NA)
fig       = plot_grid(title,subtitle,fig_main, ncol = 1, rel_heights = c(0.05,0.05,0.95))
ggsave("PointNormal_Sparsity3.pdf", fig, width = 14, height = 8)
fig

Past versions of fig103-1.png

res_df = readRDS("results/newridge10.RDS")
out = matrix(0,9,15)
lower = matrix(0,9,15)
upper = matrix(0,9,15)
for (i in 1:9) {
  out[i,] = colMeans(matrix(res_df[[i]]$pred, 20, 15))
  lower[i,] = apply(matrix(res_df[[i]]$pred, 20, 15), 2, function(x) quantile(x, probs = 0.1))
  upper[i,] = apply(matrix(res_df[[i]]$pred, 20, 15), 2, function(x) quantile(x, probs = 0.9))
}
out = out[,c(1:11,15)]
colnames(out) = c("Mr.ASH","VarBVS","BayesB","Blasso","SuSiE","E-NET","Lasso","Ridge","SCAD","MCP","L0Learn","Mr.ASH.opt")
ind = c(6,7,8,9,10,11)
out = out[,ind]
lower = lower[,ind]
upper = upper[,ind]

col   = gg_color_hue(7)[2:7]
shape = c(19,17,24,25,9,3,11,4,5,7,8,1)[2:7]
df = data.frame(s = rep(s_range, length(ind)), pred = c(out), fit = rep(colnames(out), each = 9),
                lower = c(lower), upper = c(upper))
df$fit = factor(df$fit, levels =  c("Mr.ASH","E-NET","Lasso","Ridge",
                                    "SCAD","MCP","L0Learn",
                                    "VarBVS","BayesB","Blasso","SuSiE","Mr.ASH.opt"))
df = df[df$s %in% c(1,5,50,500), ]
df = df[df$fit %in% c("Lasso","Ridge","L0Learn"),]
p1 = ggplot(df) + geom_line(aes(x = s, y = pred, color = fit), size = 1) +
  geom_point(aes(x = s, y = pred, color = fit, shape = fit), size = 3) +
  theme_cowplot(font_size = 14) +
  scale_x_continuous(trans = "log10", breaks = c(1,5,50,500)) +
  labs(y = "predictior error (rmse / sigma)", x = "number of nonzero coefficients (s)") +
  theme(axis.line = element_blank(),
        plot.title = element_text(hjust = 0.5)) +
  scale_color_manual(values = col) +
  scale_shape_manual(values = shape) +
  scale_y_continuous(trans = "log10", breaks = c(1,1.1,1.2,1.3,1.4)) +
  coord_cartesian(ylim = c(1,1.46), xlim = c(1,600)) +
  geom_point(aes(x = 500, y = 1.1942835), size = 3, shape = 4) +
  geom_text(aes(x = 500, y = 1.1942835, label="Oracle "), vjust = 0.3, hjust = -0.2) +
  #geom_text(aes(x = 500, y = df$pred[4], label="Mr.ASH", col = "Mr.ASH"), vjust = 0.2, hjust = -0.2) +
  geom_text(aes(x = 500, y = df$pred[8], label="Ridge  ", col = "Ridge"), vjust = 1, hjust = -0.22) +
  geom_text(aes(x = 500, y = df$pred[4], label="Lasso  ", col = "Lasso"), vjust = 0.2, hjust = -0.22) +
  geom_text(aes(x = 500, y = df$pred[12], label="L0Learn", col = "L0Learn"), vjust = 0.2, hjust = -0.22)
fig_main = p1
title     = ggdraw() + draw_label("Prediction Error (log-scale)", fontface = 'bold', size = 20) 
subtitle  = ggdraw() + draw_label("Scenario: IndepGauss + PointNormal, n = 1000, p = 500, s = 1-500, pve = 0.5", fontface  = 'bold', size = 18) 
fig       = plot_grid(title,subtitle,fig_main, ncol = 1, rel_heights = c(0.05,0.05,0.95))
ggsave("PointNormal_Sparsity0_without_mrash.pdf", fig, width = 14, height = 8)
fig

Past versions of fig101-1.png

fig_main  = fig_main + 
    geom_ribbon(aes(x = s, ymin = lower, ymax = upper, color = fit, fill = fit),
                linetype = 2, size = 0, alpha = 0.1, colour = NA)
fig       = plot_grid(title,subtitle,fig_main, ncol = 1, rel_heights = c(0.05,0.05,0.95))
ggsave("PointNormal_Sparsity1.pdf", fig, width = 14, height = 8)
fig

Past versions of fig102-1.png

Session information

sessionInfo()

R version 3.5.3 (2019-03-11)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: macOS High Sierra 10.13.3

Matrix products: default
BLAS: /Library/Frameworks/R.framework/Versions/3.5/Resources/lib/libRblas.0.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/3.5/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] varbvs_2.6-5       L0Learn_1.1.0      ncvreg_3.11-1     
 [4] mr.ash.alpha_0.1-3 glmnet_2.0-16      foreach_1.4.4     
 [7] BGLR_1.0.8         susieR_0.7.1       cowplot_0.9.4     
[10] ggplot2_3.2.1      Matrix_1.2-17     

loaded via a namespace (and not attached):
 [1] Rcpp_1.0.3          RColorBrewer_1.1-2  plyr_1.8.4         
 [4] compiler_3.5.3      pillar_1.3.1        git2r_0.25.2       
 [7] workflowr_1.3.0     iterators_1.0.10    tools_3.5.3        
[10] digest_0.6.18       evaluate_0.13       tibble_2.1.1       
[13] gtable_0.3.0        lattice_0.20-38     pkgconfig_2.0.2    
[16] rlang_0.4.0         yaml_2.2.0          xfun_0.6           
[19] withr_2.1.2         stringr_1.4.0       dplyr_0.8.3        
[22] knitr_1.22          fs_1.3.0            rprojroot_1.3-2    
[25] grid_3.5.3          tidyselect_0.2.5    glue_1.3.1         
[28] R6_2.4.0            rmarkdown_1.12      latticeExtra_0.6-28
[31] reshape2_1.4.3      purrr_0.3.2         magrittr_1.5       
[34] whisker_0.3-2       codetools_0.2-16    backports_1.1.4    
[37] scales_1.0.0        htmltools_0.4.0     assertthat_0.2.1   
[40] colorspace_1.4-1    labeling_0.3        nor1mix_1.2-3      
[43] stringi_1.4.3       lazyeval_0.2.2      munsell_0.5.0      
[46] truncnorm_1.0-8     crayon_1.3.4