In alex-kalinka-cruk/crispRutils: Tools for working with CRISPR data.

Setup

num_reps <- 2

set.seed(1)

options(warn=-1)
suppressMessages(library(knitr))
suppressMessages(library(tidyr))
suppressMessages(library(magrittr))
suppressMessages(library(ggplot2))
suppressMessages(library(dplyr))
suppressMessages(library(EnvStats)) # For 'rpareto'.
suppressMessages(library(seqgendiff))

# Example to illustrate sensitivity issue for Mageck under strong selection + high prevalence.
load("../data/nz_lfc.rda")
load("../data/mg_lfc.rda")

# Data all from one screen.
c903 <- read.table("../data/HT29_c903.tsv", header=T, stringsAsFactors = F)[,c(1:7)] %>%
  dplyr::select(-ERS717283.plasmid) %>%
  dplyr::filter(rowSums(.[,3:6]) > 30)

# Data from 4 different screens.
cFour <- cbind(read.table("../data/HT29_c905.tsv", header=T, stringsAsFactors = F) %>%
                 select(sgRNA, gene, HT29_c905R4),
               read.table("../data/HT29_c906.tsv", header=T, stringsAsFactors = F) %>%
                 select(HT29_c906R8),
               read.table("../data/HT29_c907.tsv", header=T, stringsAsFactors = F) %>%
                 select(HT29_c907R7),
               read.table("../data/HT29_c908.tsv", header=T, stringsAsFactors = F) %>%
                 select(HT29_c908R6)) %>%
  dplyr::filter(rowSums(.[,3:6]) > 30)

## Common functions.
write_data <- function(data){
  # Write data for algorithms to use.
  write.table(data, file = "bthin_input.txt", col.names = T, row.names = F, sep="\t", quote=F)
}

run_mageck <- function(){
  cmd <- "source ~/miniconda3/bin/activate; mageck test -k bthin_input.txt -t treat.1,treat.2 -c control.1,control.2 -n mageck_out"
  system(cmd)
}

run_drugz <- function(){
  cmd <- "source ~/miniconda3/bin/activate; export PATH=~/git-projects/drugz:$PATH; drugz.py -i bthin_input.txt -c control.1,control.2 -x treat.1,treat.2 -o drugz-out.txt -unpaired"
  system(cmd)
}

read_mageck <- function(){
  return(read.table("mageck_out.gene_summary.txt", header=T, stringsAsFactors = F))
}

read_drugz <- function(){
  return(read.table("drugz-out.txt", header=T, stringsAsFactors = F))
}

# Process output to add expected logFC.
process_output <- function(lfc, type){
  if(type == "mageck"){
    data <- read_mageck()
    gene_col <- "id"
    col.1 <- "neg.fdr"
    col.2 <- "pos.fdr"
  }else{
    data <- read_drugz()
    gene_col <- "GENE"
    col.1 <- "fdr_synth"
    col.2 <- "fdr_supp"
  }
  ret <- data %>%
    mutate(lfc_bthin = ifelse(!!sym(gene_col) %in% names(lfc), 
                              lfc[match(!!sym(gene_col), names(lfc))], 0)) %>%
    rowwise() %>%
    mutate(sign.1 = (!!sym(col.1) < 0.1 || !!sym(col.2) < 0.1),
           # LOD cut-offs for genes to ignore when calculating performance estimates.
           ignore.1.4 = (abs(lfc_bthin) > 0 && abs(lfc_bthin) < 1.4),
           ignore.1.2 = (abs(lfc_bthin) > 0 && abs(lfc_bthin) < 1.2),
           ignore.1 = (abs(lfc_bthin) > 0 && abs(lfc_bthin) < 1),
           ignore.0.8 = (abs(lfc_bthin) > 0 && abs(lfc_bthin) < 0.8),
           ignore.0.6 = (abs(lfc_bthin) > 0 && abs(lfc_bthin) < 0.6),
           ignore.0.4 = (abs(lfc_bthin) > 0 && abs(lfc_bthin) < 0.4),
           ignore.0.2 = (abs(lfc_bthin) > 0 && abs(lfc_bthin) < 0.2),
           ignore.0 = F,
           true_neg = lfc_bthin==0, true_pos = abs(lfc_bthin)>0,
           false_pos = (sign.1 && lfc_bthin==0),
           false_neg = (!sign.1 && abs(lfc_bthin)>0)) %>%
    ungroup()
  return(ret)
}

# Estimate performance.
estimate_performance <- function(data){
  pf <- NULL
  for(i in seq(0,1.4,by=0.2)){
    pf <- rbind(pf,
                data %>%
                  filter(! (!!sym(paste("ignore.",i,sep="")))) %>%
                  summarise(logFC.LOD = i,
                            Specificity = 1-sum(false_pos)/sum(true_neg),
                            Sensitivity = 1-sum(false_neg)/sum(true_pos))
    )
  }
  return(pf)
}

# Randomly sample logFC.
rand_logfc <- function(fun, params, num){
  # Params expected to be in same order as function definition.
  return(fun(num, params[[1]], params[[2]]))
}

# Replicate signal additions.
add_signal_repl <- function(data, num_reps, num_genes, signal_fun, signal_params, 
                            thin_library = F, add_signal = T){
  # 'data' expected to have plasmid column(s) removed in advance.
  ret_perf <- NULL
  ret_sig_mageck <- ret_sig_drugz <- list()
  for(i in 1:num_reps){
    # 1. Randomly sample genes to receive a signal sampled from a particular distribution.
    if(add_signal){
      lfc <- rand_logfc(signal_fun, signal_params, num_genes)
      names(lfc) <- sample(unique(data$gene), num_genes, replace=F)
    }else{
      lfc <- rep(0,length(unique(data$gene)))
      names(lfc) <- unique(data$gene)
    }

    # 2. Build coef matrix.
    coef_mat <- as.matrix(data %>%
                        dplyr::select(gene) %>%
                        dplyr::mutate(gene_indicator = ifelse(gene %in% names(lfc),1,0)) %>%
                        dplyr::group_by(gene) %>%
                        dplyr::mutate(gene_indicator = ifelse(gene[1] %in% names(lfc),
                                                              rep(lfc[names(lfc)==gene[1]],n()),
                                                              rep(0,n()))) %>%
                        dplyr::ungroup() %>%
                        dplyr::select(gene_indicator))

    # 3. Build design matrix.
    design_cols <- rep(0,4)
    treat_cols <- sample(1:4,2,replace=F)
    design_cols[treat_cols] <- 1
    design_mat <- matrix(design_cols)
    colnames(design_mat) <- "treatment"

    # 4A. Thin two of the libraries randomly by 1/4.
    if(thin_library){
      thin_cols <- sample(1:4,2,replace=F)
      scaling_factor <- rep(0,4)
      scaling_factor[thin_cols] <- 0.5
      thout_lib <- thin_lib(mat = as.matrix(data[,3:ncol(data)]),
                            thinlog2 = scaling_factor)
      data <- data.frame(data[,1:2], thout_lib$mat, stringsAsFactors = F)
    }

    if(add_signal){
      # 4B. Add signal to randomly sampled genes and their gRNAs.
      thout <- thin_diff(mat = as.matrix(data[,3:ncol(data)]), 
                        design_fixed = design_mat, 
                        coef_fixed = coef_mat)
    }else{
      thout <- list()
      thout$mat <- as.matrix(data[,3:ncol(data)])
    }

    # 5. Assemble input data for algorithms.
    data_input <- data.frame(data[,1:2], thout$mat, stringsAsFactors = F)
    colnames(data_input)[setdiff(1:4,treat_cols)+2] <- paste("control.",1:2,sep="")
    colnames(data_input)[treat_cols+2] <- paste("treat.",1:2,sep="")
    write_data(data_input)

    # 6. Run algorithms.
    run_mageck()
    run_drugz()

    # 7. Add expected logFC.
    mageck_out <- process_output(lfc, type = "mageck")
    drugz_out <- process_output(lfc, type = "drugz")

    # 8. Estimate performance.
    pf <- rbind(data.frame(Algorithm="Mageck",estimate_performance(mageck_out)),
                data.frame(Algorithm="drugZ",estimate_performance(drugz_out)))

    ret_perf <- rbind(ret_perf, pf)
    ret_sig_mageck[[i]] <- mageck_out
    ret_sig_drugz[[i]] <- drugz_out
  }
  ret <- list()
  ret$signal_mageck <- ret_sig_mageck
  ret$signal_drugz <- ret_sig_drugz
  ret$performance <- ret_perf
  return(ret)
}

# Plot sensitivity across LODs.
plot_sens_lod <- function(data){
  pl <- ggplot(data$performance, aes(logFC.LOD,Sensitivity,color=Algorithm)) +
    geom_point() +
    geom_smooth(se=F) +
    facet_wrap(~dataset)
  print(pl)
}

setwd("../tmp")

Example

From a single high prevalence (10%), strong selection (SD 1.5) example.

# Compare normZ for logFC groups (expected).
nz_lfc <- rbind(p10.s$signal_drugz[[1]] %>%
                  mutate(selection = "Strong"),
                p10.m$signal_drugz[[1]] %>%
                  mutate(selection = "Moderate"),
                p10.w$signal_drugz[[1]] %>%
                  mutate(selection = "Weak")) %>%
  filter(abs(lfc_bthin) > 0) %>%
  mutate(selection = factor(selection, levels=c("Weak","Moderate","Strong")),
         normZ.abs = abs(normZ),
         lfc_group.abs = case_when(between(abs(lfc_bthin),0.01,0.2499999) ~ "[0.01,0.25)",
                                   between(abs(lfc_bthin),0.25,0.4999999) ~ "[0.25,0.5)",
                               between(abs(lfc_bthin),0.5,0.9999999) ~ "[0.5,1)",
                               between(abs(lfc_bthin),1,1.4999999) ~ "[1,1.5)",
                               between(abs(lfc_bthin),1.5,1.9999999) ~ "[1.5,2)",
                               abs(lfc_bthin) >= 2 ~ ">=2")) %>%
  filter(!is.na(lfc_group.abs))

# Mageck comparison using both RRA p-value score and observed logFC.
mg_lfc <- rbind(p10.s$signal_mageck[[1]] %>%
                  mutate(selection = "Strong"),
                p10.m$signal_mageck[[1]] %>%
                  mutate(selection = "Moderate"),
                p10.w$signal_mageck[[1]] %>%
                  mutate(selection = "Weak")) %>%
  rowwise() %>%
  mutate(selection = factor(selection, levels=c("Weak","Moderate","Strong")),
         logfc.obs = max(c(abs(neg.lfc),abs(pos.lfc))),
         log.RAA.pval = -log(min(c(neg.score,pos.score))),
         lfc_group.abs = case_when(between(abs(lfc_bthin),0.01,0.2499999) ~ "[0.01,0.25)",
                                   between(abs(lfc_bthin),0.25,0.4999999) ~ "[0.25,0.5)",
                               between(abs(lfc_bthin),0.5,0.9999999) ~ "[0.5,1)",
                               between(abs(lfc_bthin),1,1.4999999) ~ "[1,1.5)",
                               between(abs(lfc_bthin),1.5,1.9999999) ~ "[1.5,2)",
                               abs(lfc_bthin) >= 2 ~ ">=2")) %>%
  ungroup() %>%
  filter(!is.na(lfc_group.abs))

ggplot(nz_lfc,aes(lfc_group.abs,abs(sumZ),color=selection)) +
  geom_point(position=position_dodge(width=0.75)) +
  geom_boxplot() +
  ylim(0,45) +
  ggtitle("drugZ sumZ by added signal groups: 10% prevalence")

ggplot(nz_lfc,aes(lfc_group.abs,normZ.abs,color=selection)) +
  geom_point(position=position_dodge(width=0.75)) +
  geom_boxplot() +
  ggtitle("drugZ normZ by added signal groups: 10% prevalence")

ggplot(mg_lfc,aes(lfc_group.abs,logfc.obs,color=selection)) +
  geom_point(position=position_dodge(width=0.75)) +
  geom_boxplot() +
  ggtitle("Mageck logfc_obs by added signal groups: 10% prevalence")

ggplot(mg_lfc,aes(lfc_group.abs,log.RAA.pval,color=selection)) +
  geom_point(position=position_dodge(width=0.75)) +
  geom_boxplot() +
  ggtitle("Mageck log(RAA.pval) by added signal groups: 10% prevalence")

Symmetric fold changes

Prevalence: 10%

No library thinning

# Strong selection.
p10.s <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p10.m <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1)), stringsAsFactors = F))
# Weak selection.
p10.w <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(data.frame(Selection = "Strong", Prevalence = "10%", thin_lib = F, p10.s),
             data.frame(Selection = "Moderate", Prevalence = "10%", thin_lib = F, p10.m),
             data.frame(Selection = "Weak", Prevalence = "10%", thin_lib = F, p10.w))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p10.s)
plot_sens_lod(p10.m)
plot_sens_lod(p10.w)

Library thinning

# Strong selection.
p10.s.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p10.m.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p10.w.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "10%", thin_lib = T, p10.s.tl),
             data.frame(Selection = "Moderate", Prevalence = "10%", thin_lib = T, p10.m.tl),
             data.frame(Selection = "Weak", Prevalence = "10%", thin_lib = T, p10.w.tl))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p10.s.tl)
plot_sens_lod(p10.m.tl)
plot_sens_lod(p10.w.tl)

Prevalence: 5%

No library thinning

# Strong selection.
p5.s <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p5.m <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1)), stringsAsFactors = F))
# Weak selection.
p5.w <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "5%", thin_lib = F, p5.s),
             data.frame(Selection = "Moderate", Prevalence = "5%", thin_lib = F, p5.m),
             data.frame(Selection = "Weak", Prevalence = "5%", thin_lib = F, p5.w))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p5.s)
plot_sens_lod(p5.m)
plot_sens_lod(p5.w)

Library thinning

# Strong selection.
p5.s.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p5.m.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p5.w.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "5%", thin_lib = T, p5.s.tl),
             data.frame(Selection = "Moderate", Prevalence = "5%", thin_lib = T, p5.m.tl),
             data.frame(Selection = "Weak", Prevalence = "5%", thin_lib = T, p5.w.tl))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p5.s.tl)
plot_sens_lod(p5.m.tl)
plot_sens_lod(p5.w.tl)

Prevalence: 2.5%

No library thinning

# Strong selection.
p2.5.s <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p2.5.m <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1)), stringsAsFactors = F))
# Weak selection.
p2.5.w <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "2.5%", thin_lib = F, p2.5.s),
             data.frame(Selection = "Moderate", Prevalence = "2.5%", thin_lib = F, p2.5.m),
             data.frame(Selection = "Weak", Prevalence = "2.5%", thin_lib = F, p2.5.w))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p2.5.s)
plot_sens_lod(p2.5.m)
plot_sens_lod(p2.5.w)

Library thinning

# Strong selection.
p2.5.s.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p2.5.m.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p2.5.w.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "2.5%", thin_lib = T, p2.5.s.tl),
             data.frame(Selection = "Moderate", Prevalence = "2.5%", thin_lib = T, p2.5.m.tl),
             data.frame(Selection = "Weak", Prevalence = "2.5%", thin_lib = T, p2.5.w.tl))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p2.5.s.tl)
plot_sens_lod(p2.5.m.tl)
plot_sens_lod(p2.5.w.tl)

Prevalence: 1%

No library thinning

# Strong selection.
p1.s <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p1.m <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1)), stringsAsFactors = F))
# Weak selection.
p1.w <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "1%", thin_lib = F, p1.s),
             data.frame(Selection = "Moderate", Prevalence = "1%", thin_lib = F, p1.m),
             data.frame(Selection = "Weak", Prevalence = "1%", thin_lib = F, p1.w))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p1.s)
plot_sens_lod(p1.m)
plot_sens_lod(p1.w)

Library thinning

# Strong selection.
p1.s.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p1.m.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p1.w.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "1%", thin_lib = T, p1.s.tl),
             data.frame(Selection = "Moderate", Prevalence = "1%", thin_lib = T, p1.m.tl),
             data.frame(Selection = "Weak", Prevalence = "1%", thin_lib = T, p1.w.tl))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p1.s.tl)
plot_sens_lod(p1.m.tl)
plot_sens_lod(p1.w.tl)

Prevalence: 0.1%

No library thinning

# Strong selection.
p01.s <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p01.m <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1)), stringsAsFactors = F))
# Weak selection.
p01.w <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "0.1%", thin_lib = F, p01.s),
             data.frame(Selection = "Moderate", Prevalence = "0.1%", thin_lib = F, p01.m),
             data.frame(Selection = "Weak", Prevalence = "0.1%", thin_lib = F, p01.w))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p01.s)
plot_sens_lod(p01.m)
plot_sens_lod(p01.w)

Library thinning

# Strong selection.
p01.s.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p01.m.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p01.w.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = 0, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "0.1%", thin_lib = T, p01.s.tl),
             data.frame(Selection = "Moderate", Prevalence = "0.1%", thin_lib = T, p01.m.tl),
             data.frame(Selection = "Weak", Prevalence = "0.1%", thin_lib = T, p01.w.tl))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p01.s.tl)
plot_sens_lod(p01.m.tl)
plot_sens_lod(p01.w.tl)

Asymmetric fold changes

Gaussian

gauss_mean <- 0.1

Prevalence: 10%

No library thinning

# Strong selection.
p10.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p10.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F))
# Weak selection.
p10.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "10%", thin_lib = F, p10.s.as),
             data.frame(Selection = "Moderate", Prevalence = "10%", thin_lib = F, p10.m.as),
             data.frame(Selection = "Weak", Prevalence = "10%", thin_lib = F, p10.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p10.s.as)
plot_sens_lod(p10.m.as)
plot_sens_lod(p10.w.as)

Library thinning

# Strong selection.
p10.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p10.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p10.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "10%", thin_lib = T, p10.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "10%", thin_lib = T, p10.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "10%", thin_lib = T, p10.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p10.s.tl.as)
plot_sens_lod(p10.m.tl.as)
plot_sens_lod(p10.w.tl.as)

Prevalence: 5%

No library thinning

# Strong selection.
p5.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p5.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F))
# Weak selection.
p5.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "5%", thin_lib = F, p5.s.as),
             data.frame(Selection = "Moderate", Prevalence = "5%", thin_lib = F, p5.m.as),
             data.frame(Selection = "Weak", Prevalence = "5%", thin_lib = F, p5.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p5.s.as)
plot_sens_lod(p5.m.as)
plot_sens_lod(p5.w.as)

Library thinning

# Strong selection.
p5.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p5.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p5.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "5%", thin_lib = T, p5.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "5%", thin_lib = T, p5.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "5%", thin_lib = T, p5.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p5.s.tl.as)
plot_sens_lod(p5.m.tl.as)
plot_sens_lod(p5.w.tl.as)

Prevalence: 2.5%

No library thinning

# Strong selection.
p2.5.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p2.5.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F))
# Weak selection.
p2.5.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "2.5%", thin_lib = F, p2.5.s.as),
             data.frame(Selection = "Moderate", Prevalence = "2.5%", thin_lib = F, p2.5.m.as),
             data.frame(Selection = "Weak", Prevalence = "2.5%", thin_lib = F, p2.5.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p2.5.s.as)
plot_sens_lod(p2.5.m.as)
plot_sens_lod(p2.5.w.as)

Library thinning

# Strong selection.
p2.5.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p2.5.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p2.5.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "2.5%", thin_lib = T, p2.5.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "2.5%", thin_lib = T, p2.5.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "2.5%", thin_lib = T, p2.5.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p2.5.s.tl.as)
plot_sens_lod(p2.5.m.tl.as)
plot_sens_lod(p2.5.w.tl.as)

Prevalence: 1%

No library thinning

# Strong selection.
p1.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p1.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F))
# Weak selection.
p1.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "1%", thin_lib = F, p1.s.as),
             data.frame(Selection = "Moderate", Prevalence = "1%", thin_lib = F, p1.m.as),
             data.frame(Selection = "Weak", Prevalence = "1%", thin_lib = F, p1.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p1.s.as)
plot_sens_lod(p1.m.as)
plot_sens_lod(p1.w.as)

Library thinning

# Strong selection.
p1.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p1.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p1.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "1%", thin_lib = T, p1.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "1%", thin_lib = T, p1.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "1%", thin_lib = T, p1.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p1.s.tl.as)
plot_sens_lod(p1.m.tl.as)
plot_sens_lod(p1.w.tl.as)

Prevalence: 0.1%

No library thinning

# Strong selection.
p01.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p01.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F))
# Weak selection.
p01.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "0.1%", thin_lib = F, p01.s.as),
             data.frame(Selection = "Moderate", Prevalence = "0.1%", thin_lib = F, p01.m.as),
             data.frame(Selection = "Weak", Prevalence = "0.1%", thin_lib = F, p01.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p01.s.as)
plot_sens_lod(p01.m.as)
plot_sens_lod(p01.w.as)

Library thinning

# Strong selection.
p01.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p01.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p01.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "0.1%", thin_lib = T, p01.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "0.1%", thin_lib = T, p01.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "0.1%", thin_lib = T, p01.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p01.s.tl.as)
plot_sens_lod(p01.m.tl.as)
plot_sens_lod(p01.w.tl.as)

Gamma

gamma_params.strong <- list(shape = 1, scale = 1)
gamma_params.moderate <- list(shape = 1, scale = 2)
gamma_params.weak <- list(shape = 1, scale = 3)

Prevalence: 10%

No library thinning

# Strong selection.
p10.s <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong), stringsAsFactors = F))

# Moderate selection.
p10.m <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rgamma, 
                         signal_params =gamma_params.moderate), stringsAsFactors = F))
# Weak selection.
p10.w <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak), stringsAsFactors = F))

pall <- rbind(pall,
              data.frame(Selection = "Strong", Prevalence = "10%", thin_lib = F, p10.s),
             data.frame(Selection = "Moderate", Prevalence = "10%", thin_lib = F, p10.m),
             data.frame(Selection = "Weak", Prevalence = "10%", thin_lib = F, p10.w))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p10.s)
plot_sens_lod(p10.m)
plot_sens_lod(p10.w)

Library thinning

# Strong selection.
p10.s.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong, thin_library = T), stringsAsFactors = F))

# Moderate selection.
p10.m.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate, thin_library = T), stringsAsFactors = F))
# Weak selection.
p10.w.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak, thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "10%", thin_lib = T, p10.s.tl),
             data.frame(Selection = "Moderate", Prevalence = "10%", thin_lib = T, p10.m.tl),
             data.frame(Selection = "Weak", Prevalence = "10%", thin_lib = T, p10.w.tl))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p10.s.tl)
plot_sens_lod(p10.m.tl)
plot_sens_lod(p10.w.tl)

Prevalence: 5%

No library thinning

# Strong selection.
p5.s <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong), stringsAsFactors = F))

# Moderate selection.
p5.m <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate), stringsAsFactors = F))
# Weak selection.
p5.w <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "5%", thin_lib = F, p5.s),
             data.frame(Selection = "Moderate", Prevalence = "5%", thin_lib = F, p5.m),
             data.frame(Selection = "Weak", Prevalence = "5%", thin_lib = F, p5.w))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p5.s)
plot_sens_lod(p5.m)
plot_sens_lod(p5.w)

Library thinning

# Strong selection.
p5.s.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong, thin_library = T), stringsAsFactors = F))

# Moderate selection.
p5.m.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate, thin_library = T), stringsAsFactors = F))
# Weak selection.
p5.w.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak, thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "5%", thin_lib = T, p5.s.tl),
             data.frame(Selection = "Moderate", Prevalence = "5%", thin_lib = T, p5.m.tl),
             data.frame(Selection = "Weak", Prevalence = "5%", thin_lib = T, p5.w.tl))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p5.s.tl)
plot_sens_lod(p5.m.tl)
plot_sens_lod(p5.w.tl)

Prevalence: 2.5%

No library thinning

# Strong selection.
p2.5.s <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong), stringsAsFactors = F))

# Moderate selection.
p2.5.m <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate), stringsAsFactors = F))
# Weak selection.
p2.5.w <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "2.5%", thin_lib = F, p2.5.s),
             data.frame(Selection = "Moderate", Prevalence = "2.5%", thin_lib = F, p2.5.m),
             data.frame(Selection = "Weak", Prevalence = "2.5%", thin_lib = F, p2.5.w))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p2.5.s)
plot_sens_lod(p2.5.m)
plot_sens_lod(p2.5.w)

Library thinning

# Strong selection.
p2.5.s.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong, thin_library = T), stringsAsFactors = F))

# Moderate selection.
p2.5.m.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate, thin_library = T), stringsAsFactors = F))
# Weak selection.
p2.5.w.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak, thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "2.5%", thin_lib = T, p2.5.s.tl),
             data.frame(Selection = "Moderate", Prevalence = "2.5%", thin_lib = T, p2.5.m.tl),
             data.frame(Selection = "Weak", Prevalence = "2.5%", thin_lib = T, p2.5.w.tl))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p2.5.s.tl)
plot_sens_lod(p2.5.m.tl)
plot_sens_lod(p2.5.w.tl)

Prevalence: 1%

No library thinning

# Strong selection.
p1.s <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong), stringsAsFactors = F))

# Moderate selection.
p1.m <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate), stringsAsFactors = F))
# Weak selection.
p1.w <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "1%", thin_lib = F, p1.s),
             data.frame(Selection = "Moderate", Prevalence = "1%", thin_lib = F, p1.m),
             data.frame(Selection = "Weak", Prevalence = "1%", thin_lib = F, p1.w))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p1.s)
plot_sens_lod(p1.m)
plot_sens_lod(p1.w)

Library thinning

# Strong selection.
p1.s.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong, thin_library = T), stringsAsFactors = F))

# Moderate selection.
p1.m.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate, thin_library = T), stringsAsFactors = F))
# Weak selection.
p1.w.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak, thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "1%", thin_lib = T, p1.s.tl),
             data.frame(Selection = "Moderate", Prevalence = "1%", thin_lib = T, p1.m.tl),
             data.frame(Selection = "Weak", Prevalence = "1%", thin_lib = T, p1.w.tl))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p1.s.tl)
plot_sens_lod(p1.m.tl)
plot_sens_lod(p1.w.tl)

Prevalence: 0.1%

No library thinning

# Strong selection.
p01.s <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong), stringsAsFactors = F))

# Moderate selection.
p01.m <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate), stringsAsFactors = F))
# Weak selection.
p01.w <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "0.1%", thin_lib = F, p01.s),
             data.frame(Selection = "Moderate", Prevalence = "0.1%", thin_lib = F, p01.m),
             data.frame(Selection = "Weak", Prevalence = "0.1%", thin_lib = F, p01.w))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p01.s)
plot_sens_lod(p01.m)
plot_sens_lod(p01.w)

Library thinning

# Strong selection.
p01.s.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.strong, thin_library = T), stringsAsFactors = F))

# Moderate selection.
p01.m.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.moderate, thin_library = T), stringsAsFactors = F))
# Weak selection.
p01.w.tl <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rgamma, 
                         signal_params = gamma_params.weak, thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "0.1%", thin_lib = T, p01.s.tl),
             data.frame(Selection = "Moderate", Prevalence = "0.1%", thin_lib = T, p01.m.tl),
             data.frame(Selection = "Weak", Prevalence = "0.1%", thin_lib = T, p01.w.tl))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p01.s.tl)
plot_sens_lod(p01.m.tl)
plot_sens_lod(p01.w.tl)

Gaussian - bimodal

gauss_mean <- 1.5

Prevalence: 10%

No library thinning

# Strong selection.
p10.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p10.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F))
# Weak selection.
p10.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "10%", thin_lib = F, p10.s.as),
             data.frame(Selection = "Moderate", Prevalence = "10%", thin_lib = F, p10.m.as),
             data.frame(Selection = "Weak", Prevalence = "10%", thin_lib = F, p10.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p10.s.as)
plot_sens_lod(p10.m.as)
plot_sens_lod(p10.w.as)

Library thinning

# Strong selection.
p10.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p10.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p10.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "10%", thin_lib = T, p10.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "10%", thin_lib = T, p10.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "10%", thin_lib = T, p10.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p10.s.tl.as)
plot_sens_lod(p10.m.tl.as)
plot_sens_lod(p10.w.tl.as)

Prevalence: 5%

No library thinning

# Strong selection.
p5.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p5.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F))
# Weak selection.
p5.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "5%", thin_lib = F, p5.s.as),
             data.frame(Selection = "Moderate", Prevalence = "5%", thin_lib = F, p5.m.as),
             data.frame(Selection = "Weak", Prevalence = "5%", thin_lib = F, p5.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p5.s.as)
plot_sens_lod(p5.m.as)
plot_sens_lod(p5.w.as)

Library thinning

# Strong selection.
p5.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p5.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p5.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "5%", thin_lib = T, p5.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "5%", thin_lib = T, p5.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "5%", thin_lib = T, p5.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p5.s.tl.as)
plot_sens_lod(p5.m.tl.as)
plot_sens_lod(p5.w.tl.as)

Prevalence: 2.5%

No library thinning

# Strong selection.
p2.5.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p2.5.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F))
# Weak selection.
p2.5.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "2.5%", thin_lib = F, p2.5.s.as),
             data.frame(Selection = "Moderate", Prevalence = "2.5%", thin_lib = F, p2.5.m.as),
             data.frame(Selection = "Weak", Prevalence = "2.5%", thin_lib = F, p2.5.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p2.5.s.as)
plot_sens_lod(p2.5.m.as)
plot_sens_lod(p2.5.w.as)

Library thinning

# Strong selection.
p2.5.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p2.5.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p2.5.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "2.5%", thin_lib = T, p2.5.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "2.5%", thin_lib = T, p2.5.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "2.5%", thin_lib = T, p2.5.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p2.5.s.tl.as)
plot_sens_lod(p2.5.m.tl.as)
plot_sens_lod(p2.5.w.tl.as)

Prevalence: 1%

No library thinning

# Strong selection.
p1.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p1.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F))
# Weak selection.
p1.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "1%", thin_lib = F, p1.s.as),
             data.frame(Selection = "Moderate", Prevalence = "1%", thin_lib = F, p1.m.as),
             data.frame(Selection = "Weak", Prevalence = "1%", thin_lib = F, p1.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p1.s.as)
plot_sens_lod(p1.m.as)
plot_sens_lod(p1.w.as)

Library thinning

# Strong selection.
p1.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p1.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p1.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "1%", thin_lib = T, p1.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "1%", thin_lib = T, p1.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "1%", thin_lib = T, p1.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p1.s.tl.as)
plot_sens_lod(p1.m.tl.as)
plot_sens_lod(p1.w.tl.as)

Prevalence: 0.1%

No library thinning

# Strong selection.
p01.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5)), stringsAsFactors = F))

# Moderate selection.
p01.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1)), stringsAsFactors = F))
# Weak selection.
p01.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5)), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "0.1%", thin_lib = F, p01.s.as),
             data.frame(Selection = "Moderate", Prevalence = "0.1%", thin_lib = F, p01.m.as),
             data.frame(Selection = "Weak", Prevalence = "0.1%", thin_lib = F, p01.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p01.s.as)
plot_sens_lod(p01.m.as)
plot_sens_lod(p01.w.as)

Library thinning

# Strong selection.
p01.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1.5), thin_library = T), stringsAsFactors = F))

# Moderate selection.
p01.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 1), thin_library = T), stringsAsFactors = F))
# Weak selection.
p01.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rnorm, 
                         signal_params = list(mean = gauss_mean, sd = 0.5), thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "0.1%", thin_lib = T, p01.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "0.1%", thin_lib = T, p01.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "0.1%", thin_lib = T, p01.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p01.s.tl.as)
plot_sens_lod(p01.m.tl.as)
plot_sens_lod(p01.w.tl.as)

Pareto

pareto_params.strong <- list(location = 1.5, shape = 5)
pareto_params.moderate <- list(location = 1.5, shape = 7)
pareto_params.weak <- list(location = 1.5, shape = 9)

Prevalence: 10%

No library thinning

# Strong selection.
p10.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong), stringsAsFactors = F))

# Moderate selection.
p10.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate), stringsAsFactors = F))
# Weak selection.
p10.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "10%", thin_lib = F, p10.s.as),
             data.frame(Selection = "Moderate", Prevalence = "10%", thin_lib = F, p10.m.as),
             data.frame(Selection = "Weak", Prevalence = "10%", thin_lib = F, p10.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p10.s.as)
plot_sens_lod(p10.m.as)
plot_sens_lod(p10.w.as)

Library thinning

# Strong selection.
p10.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong, thin_library = T), stringsAsFactors = F))

# Moderate selection.
p10.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rpareto, 
                         signal_params = papareto_params.weak, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate, thin_library = T), stringsAsFactors = F))
# Weak selection.
p10.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 1800, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak, thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "10%", thin_lib = T, p10.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "10%", thin_lib = T, p10.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "10%", thin_lib = T, p10.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p10.s.tl.as)
plot_sens_lod(p10.m.tl.as)
plot_sens_lod(p10.w.tl.as)

Prevalence: 5%

No library thinning

# Strong selection.
p5.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong), stringsAsFactors = F))

# Moderate selection.
p5.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate), stringsAsFactors = F))
# Weak selection.
p5.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "5%", thin_lib = F, p5.s.as),
             data.frame(Selection = "Moderate", Prevalence = "5%", thin_lib = F, p5.m.as),
             data.frame(Selection = "Weak", Prevalence = "5%", thin_lib = F, p5.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p5.s.as)
plot_sens_lod(p5.m.as)
plot_sens_lod(p5.w.as)

Library thinning

# Strong selection.
p5.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong, thin_library = T), stringsAsFactors = F))

# Moderate selection.
p5.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate, thin_library = T), stringsAsFactors = F))
# Weak selection.
p5.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 900, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak, thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "5%", thin_lib = T, p5.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "5%", thin_lib = T, p5.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "5%", thin_lib = T, p5.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p5.s.tl.as)
plot_sens_lod(p5.m.tl.as)
plot_sens_lod(p5.w.tl.as)

Prevalence: 2.5%

No library thinning

# Strong selection.
p2.5.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong), stringsAsFactors = F))

# Moderate selection.
p2.5.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate), stringsAsFactors = F))
# Weak selection.
p2.5.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "2.5%", thin_lib = F, p2.5.s.as),
             data.frame(Selection = "Moderate", Prevalence = "2.5%", thin_lib = F, p2.5.m.as),
             data.frame(Selection = "Weak", Prevalence = "2.5%", thin_lib = F, p2.5.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p2.5.s.as)
plot_sens_lod(p2.5.m.as)
plot_sens_lod(p2.5.w.as)

Library thinning

# Strong selection.
p2.5.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong, thin_library = T), stringsAsFactors = F))

# Moderate selection.
p2.5.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate, thin_library = T), stringsAsFactors = F))
# Weak selection.
p2.5.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 450, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak, thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "2.5%", thin_lib = T, p2.5.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "2.5%", thin_lib = T, p2.5.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "2.5%", thin_lib = T, p2.5.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p2.5.s.tl.as)
plot_sens_lod(p2.5.m.tl.as)
plot_sens_lod(p2.5.w.tl.as)

Prevalence: 1%

No library thinning

# Strong selection.
p1.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong), stringsAsFactors = F))

# Moderate selection.
p1.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate), stringsAsFactors = F))
# Weak selection.
p1.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "1%", thin_lib = F, p1.s.as),
             data.frame(Selection = "Moderate", Prevalence = "1%", thin_lib = F, p1.m.as),
             data.frame(Selection = "Weak", Prevalence = "1%", thin_lib = F, p1.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p1.s.as)
plot_sens_lod(p1.m.as)
plot_sens_lod(p1.w.as)

Library thinning

# Strong selection.
p1.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong, thin_library = T), stringsAsFactors = F))

# Moderate selection.
p1.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate, thin_library = T), stringsAsFactors = F))
# Weak selection.
p1.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 180, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak, thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "1%", thin_lib = T, p1.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "1%", thin_lib = T, p1.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "1%", thin_lib = T, p1.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p1.s.tl.as)
plot_sens_lod(p1.m.tl.as)
plot_sens_lod(p1.w.tl.as)

Prevalence: 0.1%

No library thinning

# Strong selection.
p01.s.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong), stringsAsFactors = F))

# Moderate selection.
p01.m.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate), stringsAsFactors = F))
# Weak selection.
p01.w.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "0.1%", thin_lib = F, p01.s.as),
             data.frame(Selection = "Moderate", Prevalence = "0.1%", thin_lib = F, p01.m.as),
             data.frame(Selection = "Weak", Prevalence = "0.1%", thin_lib = F, p01.w.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p01.s.as)
plot_sens_lod(p01.m.as)
plot_sens_lod(p01.w.as)

Library thinning

# Strong selection.
p01.s.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.strong, thin_library = T), stringsAsFactors = F))

# Moderate selection.
p01.m.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.moderate, thin_library = T), stringsAsFactors = F))
# Weak selection.
p01.w.tl.as <- rbind(data.frame(dataset = "c903", add_signal_repl(c903, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak, thin_library = T), stringsAsFactors = F),
               data.frame(dataset = "cFour", add_signal_repl(cFour, num_reps = num_reps, num_genes = 18, 
                         signal_fun = rpareto, 
                         signal_params = pareto_params.weak, thin_library = T), stringsAsFactors = F))

pall <- rbind(pall,
             data.frame(Selection = "Strong", Prevalence = "0.1%", thin_lib = T, p01.s.tl.as),
             data.frame(Selection = "Moderate", Prevalence = "0.1%", thin_lib = T, p01.m.tl.as),
             data.frame(Selection = "Weak", Prevalence = "0.1%", thin_lib = T, p01.w.tl.as))

saveRDS(pall, file="pall.rds", compress="xz")

plot_sens_lod(p01.s.tl.as)
plot_sens_lod(p01.m.tl.as)
plot_sens_lod(p01.w.tl.as)

Analysis

Relative performance taken as median of log ratio of sensitivity across all logFC.LOD - drugZ/Mageck.

pall <- rbind(readRDS(file="../Rmd/pall.rds") %>%
                mutate(fc_distr = "Normal"),
              readRDS(file="../Rmd/pall-gamma.rds") %>%
                mutate(fc_distr = "Gamma"),
              readRDS(file="../Rmd/pall-bimod.rds") %>%
                dplyr::filter(row_number() < 1921) %>%
                mutate(fc_distr = "Normal-bimodal")) %>%
  mutate(replicate = rep(1:(n()/16),each=16))

rp <- pall %>%
  dplyr::filter(!is.na(Sensitivity)) %>%
  group_by(replicate, logFC.LOD) %>%
  summarise(fc_distr = fc_distr[1], Selection = Selection[1],
            Prevalence = Prevalence[1],
            log_ratio_sens = log(Sensitivity[Algorithm=="drugZ"]/Sensitivity[Algorithm=="Mageck"])) %>%
  ungroup() %>%
  group_by(replicate) %>%
  summarise(fc_distr = fc_distr[1], Selection = Selection[1],
            Prevalence = Prevalence[1],
            log_ratio_sensitivity.median = median(log_ratio_sens)) %>%
  ungroup()

ggplot(rp, aes(Prevalence, log_ratio_sensitivity.median, color = Selection)) +
  geom_point(position = position_dodge(width=0.75)) +
  geom_boxplot(fill=NA) +
  geom_hline(yintercept = 0, linetype="dashed") +
  facet_wrap(~fc_distr) +
  ylim(-1,1) +
  ggtitle("Sensitivity Ratio: -ve Mageck better; +ve drugZ better")

alex-kalinka-cruk/crispRutils documentation built on March 13, 2021, 7:52 p.m.