Fix output folders for StressCompare and add to Pd.sh

bash.scripts/sc_TissueMapper-Pd.sh

@@ -13,3 +13,4 @@ module load R/3.4.1-gccmkl_20181025
 sh ./sc_LinkData.sh Pd
 
 Rscript ../r.scripts/sc-TissueMapper_RUN.Pd.R
+Rscript ../r.scripts/sc-TissueMapper_RUN.Pd.StressCompare.R

r.scripts/StressCompare.R → r.scripts/sc-TissueMapper_RUN.Pd.StressCompare.R

@@ -9,6 +9,12 @@ library(RColorBrewer)
 library(VennDiagram)
 library(Cairo)
 
+options(bitmapType="cairo")
+
+if (!dir.exists("../analysis/stress")){
+  dir.create("../analysis/stress")
+}
+
 load("../analysis/sc10x.Stress.Rda")
 sc10x.Stress <- SetAllIdent(object=sc10x.Stress,id="ALL")
 gene.exp <- as.data.frame(as.matrix(sc10x.Stress@raw.data[,colnames(sc10x.Stress@data)]))
@@ -31,7 +37,7 @@ genes.stress.vdb <- data.frame(vdBrink.Stress=genes.stress.vdb)
 genes.homolog <- read_delim("../genesets/Ensemble.mus-hum.txt","\t",escape_double=FALSE,trim_ws=TRUE,col_names=TRUE)
 genes.stress.vdb <- list("vdBrink.Stress"=merge(genes.stress.vdb,genes.homolog[,c(2,4)],by.x="vdBrink.Stress",by.y="Gene name")[,2])
 
-venn.diagram(list(go=unlist(genes.stress.go),dws=unlist(genes.stress.dws),vdb=unlist(genes.stress.vdb)),"../stress/Venn.tiff")
+venn.diagram(list(go=unlist(genes.stress.go),dws=unlist(genes.stress.dws),vdb=unlist(genes.stress.vdb)),"../analysis/stress/Venn.tiff")
 
 
 for (i in c("go","dws","vdb")){
@@ -53,7 +59,7 @@ for (i in c("go","dws","vdb")){
   histo.ntstress <- hist(gene.stress.pct[cell.index.ntstress],breaks=seq(0,ceiling(max(gene.stress.pct)*100)/100,0.005),main="Distribution of Stress Gene Expression in Transcriptome",xlab="% of transcriptome associated to stress genes",plot=FALSE)
   #abline(v=.0575,col="red")
   
-  Cairo(paste0("../stress/Hist.Stress.",i,".png"),width=1000,height=500,bg="white")
+  Cairo(paste0("../analysis/stress/Hist.Stress.",i,".png"),width=1000,height=500,bg="white")
   plot(histo.ntstress$mids,histo.ntstress$density,col=rgb(0,1,0,1/5),xlim=c(0,ceiling(max(gene.stress.pct)*100)/100),ylim=c(0,ceiling(max(histo.stress$density,histo.ntstress$density)/10)*10),type="h",lend="square",lwd=10,xlab="Fraction of transcriptome associated to stress genes",ylab="Density")
   par(new=TRUE)
   plot(histo.stress$mids,histo.stress$density,col=rgb(1,0,0,1/5),xlim=c(0,ceiling(max(gene.stress.pct)*100)/100),ylim=c(0,ceiling(max(histo.stress$density,histo.ntstress$density)/10)*10),type="h",lend="square",lwd=10,axes=FALSE,xlab="",ylab="")
@@ -116,13 +122,13 @@ comp <- cbind(vdb,go[match(rownames(vdb),rownames(go))],dws[match(rownames(vdb),
 colnames(comp) <- c("vdb","go","dws")
 comp <- data.frame(comp)
 
-Cairo(paste0("./Stress.vdb.vs.go.png"),width=1000,height=500,bg="white")
+Cairo(paste0("../analysis/stress/Stress.vdb.vs.go.png"),width=1000,height=500,bg="white")
 ggplot(comp,aes(x=vdb,y=go))+geom_point(size=0.01)+geom_vline(xintercept=cut.vdb.0.95)+geom_hline(yintercept=cut.go.0.95)+annotate("text",x=max(comp$vdb)/4*3,y=max(comp$go)/4*3,col="red",size=10,fontface=2,label=paste0(round(nrow(comp[comp$vdb>cut.vdb.0.95,][comp[comp$vdb>cut.vdb.0.95,"go"]>cut.go.0.95,])/nrow(comp[comp$vdb>cut.vdb.0.95,])*100,1),"%"))+annotate("text",x=max(comp$vdb)/100,y=max(comp$go)/100,col="green",size=10,fontface=2,label=paste0(round(nrow(comp[comp$vdb<=cut.vdb.0.95,][comp[comp$vdb<=cut.vdb.0.95,"go"]<=cut.go.0.95,])/nrow(comp[comp$vdb<=cut.vdb.0.95,])*100,1),"%"))
 dev.off()
-Cairo(paste0("./Stress.vdb.vs.dws.png"),width=1000,height=500,bg="white")
+Cairo(paste0("../analysis/stress/Stress.vdb.vs.dws.png"),width=1000,height=500,bg="white")
 ggplot(comp,aes(x=vdb,y=dws))+geom_point(size=0.01)+geom_vline(xintercept=cut.vdb.0.95)+geom_hline(yintercept=cut.dws.0.95)+annotate("text",x=max(comp$vdb)/4*3,y=max(comp$dws)/4*3,col="red",size=10,fontface=2,label=paste0(round(nrow(comp[comp$vdb>cut.vdb.0.95,][comp[comp$vdb>cut.vdb.0.95,"dws"]>cut.dws.0.95,])/nrow(comp[comp$vdb>cut.vdb.0.95,])*100,1),"%"))+annotate("text",x=max(comp$vdb)/100,y=max(comp$dws)/100,col="green",size=10,fontface=2,label=paste0(round(nrow(comp[comp$vdb<=cut.vdb.0.95,][comp[comp$vdb<=cut.vdb.0.95,"dws"]<=cut.dws.0.95,])/nrow(comp[comp$vdb<=cut.vdb.0.95,])*100,1),"%"))
 dev.off()
-Cairo(paste0("./Stress.go.vs.dws.png"),width=1000,height=500,bg="white")
+Cairo(paste0("../analysis/stress/Stress.go.vs.dws.png"),width=1000,height=500,bg="white")
 ggplot(comp,aes(x=go,y=dws))+geom_point(size=0.01)+geom_vline(xintercept=cut.go.0.95)+geom_hline(yintercept=cut.dws.0.95)+annotate("text",x=max(comp$go)/4*3,y=max(comp$dws)/4*3,col="red",size=10,fontface=2,label=paste0(round(nrow(comp[comp$go>cut.go.0.95,][comp[comp$go>cut.go.0.95,"dws"]>cut.dws.0.95,])/nrow(comp[comp$go>cut.go.0.95,])*100,1),"%"))+annotate("text",x=max(comp$go)/100,y=max(comp$dws)/100,col="green",size=10,fontface=2,label=paste0(round(nrow(comp[comp$go<=cut.go.0.95,][comp[comp$go<=cut.go.0.95,"dws"]<=cut.dws.0.95,])/nrow(comp[comp$go<=cut.go.0.95,])*100,1),"%"))
 dev.off()