moved plotExplanatoryVariables() to the end of the script to be able to take...

moved plotExplanatoryVariables() to the end of the script to be able to take cell cycle phase information into account; improved guess_annot_db() selection; cosmetics

sc_workflow/sc_quality_check.R

@@ -30,7 +30,7 @@ opts = docopt(doc, commandArgs(TRUE))
 #https://www.bioconductor.org/help/workflows/simpleSingleCell/
 devtools::source_url("https://git.mpi-cbg.de/bioinfo/datautils/raw/v1.40/R/core_commons.R")
 devtools::source_url("https://git.mpi-cbg.de/bioinfo/datautils/raw/v1.45/R/ggplot_commons.R")
-devtools::source_url("https://raw.githubusercontent.com/holgerbrandl/datautils/v1.40/R/bio/diffex_commons.R")
+devtools::source_url("https://git.mpi-cbg.de/bioinfo/ngs_tools/raw/v10/dge_workflow/diffex_commons.R")
 devtools::source_url("https://git.mpi-cbg.de/bioinfo/ngs_tools/raw/v6/common/cp_utils.R")
 load_pack(SingleCellExperiment)
 load_pack(knitr)
@@ -89,7 +89,7 @@ rm(countMatrix)
 
 # get gene information if requested
 if(any(str_detect(rownames(counts(sce)), "^FBgn|^ENSCAFG|^ENSMUSG|^ENSDARG|^ENSPTRG|^ENSG|^ENSCGRG"))) {
-    db <- guess_anno_db(rownames(counts(sce))[1])
+    db <- guess_anno_db(rownames(counts(sce))[which(str_detect(rownames(counts(sce)), "^FBgn|^ENSCAFG|^ENSMUSG|^ENSDARG|^ENSPTRG|^ENSG|^ENSCGRG"))])
 } else {
     db <- ""
 }
@@ -265,13 +265,8 @@ print("Number of outliers detected")
 summary(sce$outlier)
 
 
-#' <br><br>
-#'
-#' ## Relationship of experimental factors and expression
-
 #  log-normalization
 sce <- normalize(sce)
-scater::plotExplanatoryVariables(sce, variables = c(colnames(cell_infos)))
 
 
 #' <br><br>
@@ -319,9 +314,16 @@ if (exists("cc.pairs")) {
 
     print(cc_plot)
     datatable(data.frame(G1 = length(sce$phase[sce$phase == "G1"]), S = length(sce$phase[sce$phase == "S"]), G2M = length(sce$phase[sce$phase == "G2M"])))
+
+    cell_infos %<>% left_join(data %>% select(cell_id, phase))
 }
 
 
+#' <br><br>
+#'
+#' ## Relationship of experimental factors and expression
+scater::plotExplanatoryVariables(sce, variables = c(colnames(cell_infos)))
+
 
 colData(sce) %>% as.data.frame() %>% write_tsv("scater_quality_metrics.txt")
 
@@ -330,5 +332,5 @@ colData(sce) %>% as.data.frame() %>% write_tsv("scater_quality_metrics.txt")
 # get R version and package infos
 writeLines(capture.output(devtools::session_info()), ".sessionInfo.txt")
 
-session::save.session(paste(".", prefix, "quality_check_filtered_genes.dat", sep = ""))
+session::save.session(paste(".", prefix, "sc_quality_check.dat", sep = ""))
 #session::restore.session(".quality_check_filtered_genes.dat")