In nbenn/singleCellFeatures: Analyze single cell features stored in an openBIS instance

suppressPackageStartupMessages(library(knitr))
suppressPackageStartupMessages(library(singleCellFeatures))
suppressPackageStartupMessages(library(xtable))
options(scipen=1, digits=2, singleCellFeatures.progressBars="none")
opts_chunk$set(cache.extra = rand_seed)

suppressPackageStartupMessages(library(singleCellFeatures))
suppressPackageStartupMessages(library(xtable))
mtor <- findWells(contents="MTOR", experiment="brucella-du-k[12]")
scra <- findWells(well.names=c("E2", "G23", "H2", "J2"),
                  plates=sapply(mtor, getBarcode))

data <- unlist(getSingleCellData(c(mtor, scra)), recursive=FALSE)
data <- lapply(data, cleanData, "lower")
data <- makeFeatureCompatible(data)
reps  <- 100

To find out if good predition/data separation is only given in special situations or occurs readily, stability analysis using glmnet is performed on the MTOR (H6)/SCRAMBLED (H2) pair and on several randomly selected well pairs on the same plates.

First, wells containing siRNA for the gene MTOR are searched for within the kinome-wide Dharmacon unpooled screens (replicates 1 and 2). Then on the plates containing those wells, scrambled control experiments are looked up (in well H2). The data for the resulting 16 wells is loaded, cleaned up (lower 5% quantile in terms of well-level cell counts is discarded) and melted into data frames.

mtor.all1 <- glmBootstrapStability(data[grep(".E2$", names(data))],
                                   data[grep(".H6$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="all", alpha=0.5)
mtor.all2 <- glmBootstrapStability(data[grep(".G23$", names(data))],
                                   data[grep(".H6$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="all", alpha=0.5)
mtor.all3 <- glmBootstrapStability(data[grep(".H2$", names(data))],
                                   data[grep(".H6$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="all", alpha=0.5)
mtor.all4 <- glmBootstrapStability(data[grep(".J2$", names(data))],
                                   data[grep(".H6$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="all", alpha=0.5)

scra.all1 <- glmBootstrapStability(data[grep(".E2$", names(data))],
                                   data[grep(".G23$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="all", alpha=0.5)
scra.all2 <- glmBootstrapStability(data[grep(".E2$", names(data))],
                                   data[grep(".J2$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="all", alpha=0.5)
scra.all3 <- glmBootstrapStability(data[grep(".H2$", names(data))],
                                   data[grep(".J2$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="all", alpha=0.5)
scra.all4 <- glmBootstrapStability(data[grep(".H2$", names(data))],
                                   data[grep(".G23$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="all", alpha=0.5)

\newpage \blandscape

table.all <- Reduce(function(x, y) {
  res <- merge(x, y, all=TRUE, by=0)
  rownames(res) <- res$Row.names
  res$Row.names <- NULL
  return(res)
}, list(mtor.all1, mtor.all2, mtor.all3, mtor.all4,
        scra.all1, scra.all2, scra.all3, scra.all4),
accumulate=FALSE)

colnames(table.all) <- c("E2/H6", "G23/H6", "H2/H6", "J2/H6", 
                         "E2/G23", "E2/J2", "H2/J2", "H2/G23")
table.all[is.na(table.all)] <- 0

table.all$sums <- rowSums(table.all)
table.all <- table.all[order(table.all$sums, decreasing=TRUE),]
print(xtable(head(table.all, n=50), digits=0), size="footnotesize",
      comment=FALSE)

\elandscape \newpage

Stability analysis on MTOR (Well H6) versus SCRAMBLED (Well H2) is performed with 100 resampling iterations, each run on 70% of the data. The top 20 coefficients of each iteration are selected and their frequencies are tabulated. K1 consists of merged *-2C wells, K2 corresponds to *-2D and DU encompasses all 8 well pairs. Models were fit with glmnet, and default parameters (elastic net penalty: $\alpha = 0.5$).

mtor.inf1 <- glmBootstrapStability(data[grep(".E2$", names(data))],
                                   data[grep(".H6$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="infected", alpha=0.5)
mtor.inf2 <- glmBootstrapStability(data[grep(".G23$", names(data))],
                                   data[grep(".H6$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="infected", alpha=0.5)
mtor.inf3 <- glmBootstrapStability(data[grep(".H2$", names(data))],
                                   data[grep(".H6$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="infected", alpha=0.5)
mtor.inf4 <- glmBootstrapStability(data[grep(".J2$", names(data))],
                                   data[grep(".H6$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="infected", alpha=0.5)

scra.inf1 <- glmBootstrapStability(data[grep(".E2$", names(data))],
                                   data[grep(".G23$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="infected", alpha=0.5)
scra.inf2 <- glmBootstrapStability(data[grep(".E2$", names(data))],
                                   data[grep(".J2$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="infected", alpha=0.5)
scra.inf3 <- glmBootstrapStability(data[grep(".H2$", names(data))],
                                   data[grep(".J2$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="infected", alpha=0.5)
scra.inf4 <- glmBootstrapStability(data[grep(".H2$", names(data))],
                                   data[grep(".G23$", names(data))],
                                   n.rep=reps, norm.method="none", 
                                   select.inf="infected", alpha=0.5)

\newpage \blandscape

table.inf <- Reduce(function(x, y) {
  res <- merge(x, y, all=TRUE, by=0)
  rownames(res) <- res$Row.names
  res$Row.names <- NULL
  return(res)
}, list(mtor.inf1, mtor.inf2, mtor.inf3, mtor.inf4,
        scra.inf1, scra.inf2, scra.inf3, scra.inf4),
accumulate=FALSE)

colnames(table.inf) <- c("E2/H6", "G23/H6", "H2/H6", "J2/H6", 
                         "E2/G23", "E2/J2", "H2/J2", "H2/G23")
table.inf[is.na(table.inf)] <- 0

table.inf$sums <- rowSums(table.inf)
table.inf <- table.inf[order(table.inf$sums, decreasing=TRUE),]
print(xtable(head(table.inf, n=50), digits=0), size="footnotesize",
      comment=FALSE)

\elandscape \newpage

Stability analysis on MTOR (Well H6) versus SCRAMBLED (Well H2) is performed with 100 resampling iterations, each run on 70% of the data. The top 20 coefficients of each iteration are selected and their frequencies are tabulated. K1 consists of merged *-2C wells, K2 corresponds to *-2D and DU encompasses all 8 well pairs. Models were fit with glmnet, and default parameters (elastic net penalty: $\alpha = 0.5$).

nbenn/singleCellFeatures documentation built on May 23, 2019, 12:24 p.m.