The FDP for different \(g_1\) figure in cashr

Introduction

library(ggplot2)

boxplot.quantile <- function(x) {
  r <- quantile(x, probs = c(0.10, 0.25, 0.5, 0.75, 0.90))
  names(r) <- c("ymin", "lower", "middle", "upper", "ymax")
  return(r)
}

boxplot.quantile.sq <- function (x) {
  r <- sqrt(quantile(x^2, probs = c(0.10, 0.25, 0.5, 0.75, 0.90)))
  names(r) <- c("ymin", "lower", "middle", "upper", "ymax")
  return(r)
}

mean.sq <- function (x) {
  r <- sqrt(mean(x^2))
  return(r)
}

mysqrt_trans <- function() {
  scales::trans_new("mysqrt", 
                    transform = base::sqrt,
                    inverse = function(x) ifelse(x<0, 0, x^2),
                    domain = c(0, Inf))
}

q <- 0.1
method.name.FDR <- c("cashr", "BH", "qvalue", "ashr", "locfdr")
method.col.FDR <- scales::hue_pal()(length(method.name.FDR))[c(5, 1, 2, 4, 3)]

density.ggdata.normal <- readRDS("../output/paper/simulation/density.ggdata.normal.rds")
density.ggdata.nearnormal <- readRDS("../output/paper/simulation/density.ggdata.nearnormal.rds")
density.ggdata.spiky <- readRDS("../output/paper/simulation/density.ggdata.spiky.rds")
density.ggdata.flattop <- readRDS("../output/paper/simulation/density.ggdata.flattop.rds")
density.ggdata.skew <- readRDS("../output/paper/simulation/density.ggdata.skew.rds")
density.ggdata.bimodal <- readRDS("../output/paper/simulation/density.ggdata.bimodal.rds")

density.g.ggdata <- rbind.data.frame(
  density.ggdata.normal,
  density.ggdata.nearnormal,
  density.ggdata.spiky,
  density.ggdata.flattop,
  density.ggdata.skew,
  density.ggdata.bimodal
)

density.g.plot <- ggplot(data = density.g.ggdata, aes(x = plotx, y = ploty)) +
  geom_line() +
  facet_wrap(~g, nrow = 1) +
  labs(x = expression(theta), y = expression(g[1](theta))) +
  theme(# axis.title.x = element_text(size = 15),
        axis.title.x = element_blank(),
        axis.text.x = element_text(size = 10),
        axis.title.y = element_text(size = 15),
        axis.text.y = element_text(size = 10),
        strip.text = element_text(size = 15),
        legend.position = "none",
        legend.background = element_rect(color = "grey"),
        legend.text = element_text(size = 12))

FDP.q.all.ggdata.normal <- readRDS("../output/paper/simulation/FDP.q.all.ggdata.normal.rds")
FDP.q.all.ggdata.nearnormal <- readRDS("../output/paper/simulation/FDP.q.all.ggdata.nearnormal.rds")
FDP.q.all.ggdata.spiky <- readRDS("../output/paper/simulation/FDP.q.all.ggdata.spiky.rds")
FDP.q.all.ggdata.flattop <- readRDS("../output/paper/simulation/FDP.q.all.ggdata.flattop.rds")
FDP.q.all.ggdata.skew <- readRDS("../output/paper/simulation/FDP.q.all.ggdata.skew.rds")
FDP.q.all.ggdata.bimodal <- readRDS("../output/paper/simulation/FDP.q.all.ggdata.bimodal.rds")

FDP.q.g.ggdata <- rbind.data.frame(
  FDP.q.all.ggdata.normal,
  FDP.q.all.ggdata.nearnormal,
  FDP.q.all.ggdata.spiky,
  FDP.q.all.ggdata.flattop,
  FDP.q.all.ggdata.skew,
  FDP.q.all.ggdata.bimodal
)

FDP.q.g.ggdata$Method <- plyr::mapvalues(FDP.q.g.ggdata$Method, from = c("CASH", "BHq", "ASH"), to = c("cashr", "BH", "ashr"))

FDP.q.g.ggdata$Method <- factor(FDP.q.g.ggdata$Method, levels = c("cashr", "BH", "qvalue", "ashr", "locfdr"))

FDP.q.g.ggdata$pi0 <- as.numeric(levels(FDP.q.g.ggdata$pi0))[FDP.q.g.ggdata$pi0]

FDP.q.g.plot <- ggplot(data = FDP.q.g.ggdata, aes(x = Method, y = FDP, fill = Method, color = Method)) +
  stat_summary(fun.data = boxplot.quantile, geom = "boxplot", position = "dodge", aes(width = 0.75), show.legend = FALSE) +
  stat_summary(fun.y = mean, geom = "point", position = position_dodge(width = 0.9), show.legend = FALSE, shape = 13, size = 3) +
  scale_color_manual(labels = method.name.FDR, values = method.col.FDR) +
  scale_fill_manual(labels = method.name.FDR, values = alpha(method.col.FDR, 0.35)) +
  scale_y_continuous(trans = "mysqrt", breaks = c(0, 0.1, 0.2, 0.4, 0.6, 0.8)) +
  facet_grid(pi0 ~ g, labeller = label_bquote(rows = pi[0] == .(pi0))) +
  scale_x_discrete(limits = rev(levels(FDP.q.g.ggdata$Method))) +
  coord_flip() +
  geom_hline(yintercept = q, col = "black", linetype = "dashed", size = 1) +
  labs(y = "FDP", title = bquote(paste("Nominal FDR = ", .(q)))) +
  theme(plot.title = element_text(size = 12, hjust = 0),
        axis.title.y = element_blank(),
        axis.text.y = element_text(size = 15),
        axis.title.x = element_text(size = 15),
        axis.text.x = element_text(size = 10),
        strip.text = element_text(size = 15),
        legend.position = "bottom",
        legend.background = element_rect(color = "grey"),
        legend.text = element_text(size = 12))

Warning: Ignoring unknown aesthetics: width

FDP.q.g.plot.save <- gridExtra::arrangeGrob(
  density.g.plot + 
    theme(plot.margin = grid::unit(c(5.5, 30, 5.5, 10), "points")) +
    theme(plot.title = element_blank()),
  FDP.q.g.plot + theme(strip.text.x = element_blank(), plot.title = element_blank()),
  heights = c(1, 2),
  top = grid::textGrob(bquote(paste("Nominal FDR = ", .(q))), gp = grid::gpar(fontsize = 15), hjust = 2.35)
)

Warning: Removed 8 rows containing non-finite values (stat_summary).

Warning: Removed 8 rows containing non-finite values (stat_summary).

ggsave("../output/paper/FDP.q.g.pdf", FDP.q.g.plot.save, height = 6, width = 9)

pi0.plot <- 0.9

FDP.q.g.pi0.plot <- ggplot(data = FDP.q.g.ggdata[FDP.q.g.ggdata$pi0 == pi0.plot, ], aes(x = Method, y = FDP, fill = Method, color = Method)) +
  stat_summary(fun.data = boxplot.quantile, geom = "boxplot", position = "dodge", aes(width = 0.5), show.legend = FALSE) +
  stat_summary(fun.y = mean, geom = "point", position = position_dodge(width = 0.9), show.legend = FALSE) +
  scale_color_manual(labels = method.name.FDR, values = method.col.FDR) +
  scale_fill_manual(labels = method.name.FDR, values = alpha(method.col.FDR, 0.35)) +
  scale_y_continuous(trans = "mysqrt", breaks = c(0, 0.1, 0.2, 0.4, 0.6, 0.8)) +
  facet_wrap(~ g, nrow = 1) +
  scale_x_discrete(limits = rev(levels(FDP.q.g.ggdata$Method))) +
  coord_flip() +
  geom_hline(yintercept = q, col = "black", linetype = "dashed", size = 1) +
  labs(y = "FDP", title = bquote(paste("Nominal FDR = ", .(q), "(", pi[0] == .(pi0.plot), ")"))) +
  theme(plot.title = element_text(size = 12, hjust = 0),
        axis.title.y = element_blank(),
        axis.text.y = element_text(size = 15),
        axis.title.x = element_text(size = 15),
        axis.text.x = element_text(size = 10),
        strip.text = element_text(size = 15),
        legend.position = "bottom",
        legend.background = element_rect(color = "grey"),
        legend.text = element_text(size = 12))

Warning: Ignoring unknown aesthetics: width

FDP.q.g.pi0.plot.save <- gridExtra::arrangeGrob(
  density.g.plot,
  FDP.q.g.pi0.plot,
  heights = c(1, 1.1)
)

ggsave("../output/paper/FDP.q.g.pi0.pdf", FDP.q.g.pi0.plot.save, height = 4, width = 12)

Session information

sessionInfo()

R version 3.4.3 (2017-11-30)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: macOS High Sierra 10.13.6

Matrix products: default
BLAS: /Library/Frameworks/R.framework/Versions/3.4/Resources/lib/libRblas.0.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/3.4/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] ggplot2_2.2.1

loaded via a namespace (and not attached):
 [1] Rcpp_0.12.16      knitr_1.20        whisker_0.3-2    
 [4] magrittr_1.5      workflowr_1.1.1   munsell_0.4.3    
 [7] colorspace_1.3-2  rlang_0.2.0       stringr_1.3.1    
[10] plyr_1.8.4        tools_3.4.3       grid_3.4.3       
[13] gtable_0.2.0      R.oo_1.22.0       git2r_0.21.0     
[16] htmltools_0.3.6   yaml_2.1.19       lazyeval_0.2.1   
[19] rprojroot_1.3-2   digest_0.6.15     tibble_1.4.2     
[22] gridExtra_2.3     reshape2_1.4.3    R.utils_2.6.0    
[25] evaluate_0.10.1   rmarkdown_1.9     labeling_0.3     
[28] stringi_1.2.2     pillar_1.2.2      compiler_3.4.3   
[31] scales_0.5.0      backports_1.1.2   R.methodsS3_1.7.1