Add plotViolinPlot

workflow/scripts/plotFeature.R

@@ -60,7 +60,7 @@ for (i in chunk$V1){
   if (opt$u == TRUE){
     dr <- dr.umap
     dr$Expression <- as.numeric(exp[which(exp$genes==i),])[-ncol(exp)]
-    png(paste0("./umap/",i,".png"),width=400,height=400,bg="white")
+    png(paste0("./umap/",i,".png"),width=500,height=500,bg="white")
     plot(ggplot(dr,aes(x=UMAP_1,y=UMAP_2,col=Expression))+geom_point(size=2)+scale_colour_gradient(low="grey",high="darkred",limits=c(lo,hi))+theme_cowplot())
     dev.off()
     rm(dr)

workflow/scripts/plotIdentGroup.R

@@ -27,7 +27,7 @@ if (opt$t == TRUE){
   if (!dir.exists("./tsne")){
     dir.create("./tsne")
   }
-  png(paste0("./tsne/ClusterGroup_tSNE_",opt$g,".png"),width=400,height=400,bg="white")
+  png(paste0("./tsne/ClusterGroup_tSNE_",opt$g,".png"),width=500,height=500,bg="white")
   plot(ggplot(dr.tsne,aes_string(x="tSNE_1",y="tSNE_2",col=opt$g))+geom_point(size=2)+labs(col="ID")+theme_cowplot())
   dev.off()
 }
@@ -38,7 +38,7 @@ if (opt$u == TRUE){
   if (!dir.exists("./umap")){
     dir.create("./umap")
   }
-  png(paste0("./umap/ClusterGroup_UMAP_",opt$g,".png"),width=400,height=400,bg="white")
+  png(paste0("./umap/ClusterGroup_UMAP_",opt$g,".png"),width=500,height=500,bg="white")
   plot(ggplot(dr.umap,aes_string(x="UMAP_1",y="UMAP_2",col=opt$g))+geom_point(size=2)+labs(col="ID")+theme_cowplot())
   dev.off()
 }
\ No newline at end of file

workflow/scripts/plotViolinBox.R

+if ("optparse" %in% row.names(installed.packages()) == FALSE){
+  install.packages("optparse",repos="http://cran.us.r-project.org")
+}
+library(optparse)
+library(Matrix)
+library(readr)
+library(ggplot2)
+library(cowplot)
+
+#Get param inputs
+option_list=list(
+  make_option("--g",action="store",type='character',help="Seurat group ident"),
+  make_option("--c",action="store",type='character',help="Gene chunk filename")
+)
+opt=parse_args(OptionParser(option_list=option_list))
+rm(option_list)
+
+#Load group
+id <- read.csv(paste0("./group_",opt$g,".csv"),header=FALSE)
+names(id)[2] <- opt$g
+
+#Load data
+exp <- readMM("./expression.mtx")
+exp <- as.data.frame(as.matrix(exp))
+genes <- read.table("./genes.csv",quote="\"")
+exp$genes <- genes$V1
+
+#Load gene chunk
+chunk <- read.csv(opt$c,header=FALSE)
+
+#Filter expression data
+exp <- exp[exp$genes %in% chunk$V1,]
+
+#Create output folders
+if (!dir.exists(paste0("./",opt$g))){
+  dir.create(paste0("./",opt$g))
+}
+
+#Plot
+for (i in chunk$V1){
+  id.exp <- id
+  id.exp$Expression <- as.numeric(exp[which(exp$genes==i),])[-ncol(exp)]
+  png(paste0("./",opt$g,"/",i,".png"),width=1000,height=500,bg="white")
+  plot(ggplot(id.exp,aes_string(x=opt$g,y="Expression",fill=opt$g))+geom_violin(scale="width",trim=TRUE)+labs(y="Expression")+geom_boxplot(width=0.25,fill="white",outlier.shape=NA)+theme_cowplot())
+  dev.off()
+  rm(id.exp)
+}

workflow/scrna-seq_imgen.nf

@@ -34,6 +34,7 @@ scale = params.scale
 sample = params.sample
 groups_createGroup = file(params.groups).readLines()
 groups_plotIdentGroup = file(params.groups).readLines()
+groups_plotViolinPlot = file(params.groups)
 convert = params.convert
 
 outDir = params.outDir
@@ -95,7 +96,7 @@ process createGroup {
 //Duplicate group IDs
 id.into {
     id_plotIdentGroup
-    id_plotFeature
+    id_plotViolinBox
 }
 
 /*
@@ -103,7 +104,7 @@ id.into {
 */
 process renameGenes {
     tag "NCBI"
-    publishDir "${outDir}/matrix", mode: 'copy', overwrite: false
+    //publishDir "${outDir}/matrix", mode: 'copy', overwrite: false
 
     input:
     file matrixRaw
@@ -137,6 +138,7 @@ process renameGenes {
 matrixFinal.into {
     matrix_plotIdentGroup
     matrix_plotFeature
+    matrix_plotViolinBox
 }
 
 /*
@@ -168,10 +170,11 @@ process plotIdentGroup {
 // Split genes into chunks
 genes.splitText(by: 1000).into {
     genesChunks_plotFeature
+    genesChunks_plotViolinBox
 }
 
 /*
-* plotFeature: plot feature plots in chunks
+* plotFeature: plot feature plots in gene chunks
 */
 process plotFeature {
     publishDir "${outDir}/png/Feature", mode: 'move', overwrite: true
@@ -192,4 +195,32 @@ process plotFeature {
     module load seurat/3.0.0
     seurat-Rscript ${baseDir}/scripts/plotFeature.R --t ${tsne_plotFeature} --u ${umap_plotFeature} --c ${chunk}
     """
+}
+
+/*
+* plotViolinBox: plot violin/box plots in gene chunks
+*/
+process plotViolinBox {
+    publishDir "${outDir}/png/ViolinBox", mode: 'move', overwrite: true
+
+    input:
+    file matrix_plotViolinBox
+    file groups_plotViolinPlot
+    file("*") from id_plotViolinBox.collect()
+    each file(chunk) from genesChunks_plotViolinBox
+
+    output:
+    file("**/*.png") into plotsViolinBox
+
+    script:
+    """
+    hostname
+    ulimit -a
+    module load seurat/3.0.0
+    while IFS="," read group
+    do
+        seurat-Rscript ${baseDir}/scripts/plotViolinBox.R --g \${group} --c ${chunk} &
+    done < groups.txt
+    wait
+    """
 }
\ No newline at end of file