Update code for function based pipeline

r.scripts/sc-TissueMapper.R

@@ -744,9 +744,9 @@ write.table(round(prop.table(table(sc10x@meta.data[,i.1],sc10x@meta.data[,i.2]),
 }
 
 
-scDEG <- function(sc10x,i="Lin",g.1="Epi",g.2="St",t=2){
+scDEG <- function(sc10x,i="Lin",g.1="Epi",g.2="St",pct=0.25,t=2){
   sc10x <- SetAllIdent(object=sc10x,id=i)
-  genes.deg <- FindMarkers(object=sc10x,ident.1=g.1,ident.2=g.2,min.pct=0.5,only.pos=TRUE,thresh.use=0.25,logfc.threshold=log(t))
+  genes.deg <- FindMarkers(object=sc10x,ident.1=g.1,ident.2=g.2,min.pct=pct,only.pos=TRUE,thresh.use=0.25,logfc.threshold=log(t))
   
   return(genes.deg)
 }

r.scripts/sc-TissueMapper_RUN.D17.R

@@ -265,7 +265,7 @@ sc10x <- scMerge(sc10x,sc10x.Epi,sc10x.St,i.1="Merge",i.2="Merge",nm="Merge_Epi.
 sc10x <- scMerge(sc10x,sc10x,sc10x.Epi.NE,i.1="Merge_Epi.dws_St.go",i.2="NE",nm="Merge_Epi.dws_St.go_NE")
 
 sc10x <- SetAllIdent(object=sc10x,id="Merge_Epi.dws_St.go")
-sc10x@ident <- factor(sc10x@ident,levels=c("BE","LE","OE","NE","Fib","SM","Endo","Leu"))
+sc10x@ident <- factor(sc10x@ident,levels=c("BE","LE","OE1","OE2","Fib","SM","Endo","Leu"))
 postscript("./analysis/tSNE/FINAL/tSNE_FINAL.eps")
 plot <- TSNEPlot(object=sc10x,pt.size=2.5,do.return=TRUE,vector.friendly=FALSE)
 plot <- plot+theme(axis.text.x=element_text(size=20),axis.text.y=element_text(size=20),axis.title.x=element_text(size=20),axis.title.y=element_text(size=20),legend.text=element_text(size=20))
@@ -277,22 +277,22 @@ scTables(sc10x,i.1="samples",i.2="Merge_Epi.dws_St.go")
 scTables(sc10x,i.1="samples",i.2="Merge_Epi.dws_St.go_NE")
 scTables(sc10x,i.1="Merge_Epi.dws_St.go_NE",i.2="Merge_Epi.dws_St.go")
 
-genes.deg.Stress <- scDEG(sc10x.Stress,i="ALL.Stress",g.1="Stress",g.2="ALL",t=5)
+genes.deg.Stress <- scDEG(sc10x.Stress,i="Stress",g.1="Stress",g.2="ALL",pct=0.5,t=2)
 
 genes.deg.Epi <- scDEG(sc10x,i="Lin",g.1="Epi",g.2="St",t=2)
 genes.deg.St <- scDEG(sc10x,i="Lin",g.1="St",g.2="Epi",t=2)
 
-genes.deg.BE <- scDEG(sc10x.Epi.NE,i="Epi.dws.sub",g.1="BE",g.2=c("LE","OE1","OE2"),t=2)
-genes.deg.LE <- scDEG(sc10x.Epi.NE,i="Epi.dws.sub",g.1="LE",g.2=c("BE","LE","OE1"),t=2)
-genes.deg.OE1 <- scDEG(sc10x.Epi.NE,i="Epi.dws.sub",g.1="OE1",g.2=c("BE","LE","OE2"),t=2)
-genes.deg.OE2 <- scDEG(sc10x.Epi.NE,i="Epi.dws.sub",g.1="OE2",g.2=c("BE","LE","OE1"),t=2)
+genes.deg.BE <- scDEG(sc10x.Epi.NE,i="Epi.dws.sub",g.1="BE",g.2=c("LE","OE1","OE2"),pct=0.25,t=2)
+genes.deg.LE <- scDEG(sc10x.Epi.NE,i="Epi.dws.sub",g.1="LE",g.2=c("BE","LE","OE1"),pct=0.25,t=2)
+genes.deg.OE1 <- scDEG(sc10x.Epi.NE,i="Epi.dws.sub",g.1="OE1",g.2=c("BE","LE","OE2"),pct=0.25,t=2)
+genes.deg.OE2 <- scDEG(sc10x.Epi.NE,i="Epi.dws.sub",g.1="OE2",g.2=c("BE","LE","OE1"),pct=0.25,t=2)
 
-genes.deg.NE <- scDEG(sc10x.Epi.NE,i="NE",g.1="NE",g.2="ALL",t=2)
+genes.deg.NE <- scDEG(sc10x.Epi.NE,i="NE",g.1="NE",g.2="ALL",pct=0.25,t=2)
 
-genes.deg.Fib <- scDEG(sc10x.St,i="Merge_Epi.dws_St.go_NE",g.1="Fib",g.2=c("SM","Endo","Leu"),t=2)
-genes.deg.SM <- scDEG(sc10x.St,i="Merge_Epi.dws_St.go_NE",g.1="SM",g.2=c("Fib","Endo","Leu"),t=2)
-genes.deg.Endo <- scDEG(sc10x.St,i="Merge_Epi.dws_St.go_NE",g.1="Endo",g.2=c("Fib","SM","Leu"),t=2)
-genes.deg.Leu <- scDEG(sc10x.St,i="St.go",g.1="Leu",g.2=c("Fib","SM","Endo"),t=2)
+genes.deg.Fib <- scDEG(sc10x.St,i="Merge_Epi.dws_St.go_NE",g.1="Fib",g.2=c("SM","Endo","Leu"),pct=0.25,t=2)
+genes.deg.SM <- scDEG(sc10x.St,i="Merge_Epi.dws_St.go_NE",g.1="SM",g.2=c("Fib","Endo","Leu"),pct=0.25,t=2)
+genes.deg.Endo <- scDEG(sc10x.St,i="Merge_Epi.dws_St.go_NE",g.1="Endo",g.2=c("Fib","SM","Leu"),pct=0.25,t=2)
+genes.deg.Leu <- scDEG(sc10x.St,i="St.go",g.1="Leu",g.2=c("Fib","SM","Endo"),pct=0.25,t=2)
 
 genes.deg.BE.unique <- setdiff(rownames(genes.deg.BE),Reduce(union,list(rownames(genes.deg.LE),rownames(genes.deg.OE1),rownames(genes.deg.OE2),rownames(genes.deg.NE),rownames(genes.deg.Fib),rownames(genes.deg.SM),rownames(genes.deg.Endo),rownames(genes.deg.Leu))))
 genes.deg.LE.unique <- setdiff(rownames(genes.deg.LE),Reduce(union,list(rownames(genes.deg.BE),rownames(genes.deg.OE1),rownames(genes.deg.OE2),rownames(genes.deg.NE),rownames(genes.deg.Fib),rownames(genes.deg.SM),rownames(genes.deg.Endo),rownames(genes.deg.Leu))))
@@ -310,16 +310,17 @@ genes.deg.10 <- c(genes.deg.BE.unique[1:10],genes.deg.LE.unique[1:10],genes.deg.
 genes.deg.10 <- rev(genes.deg.10)
 
 
+sc10x <- SetAllIdent(object=sc10x,id="Merge_Epi.dws_St.go_NE")
 sc10x <- SetAllIdent(object=sc10x,id="Merge_Epi.dws_St.go_NE")
 sc10x@ident <- factor(sc10x@ident,levels=c("BE","LE","OE1","OE2","NE","Fib","SM","Endo","Leu"))
 postscript("./analysis/deg/Dot.DEG.eps",paper="special",width=10,height=5,horizontal=FALSE)
 DotPlot(sc10x,genes.deg.5,x.lab.rot=TRUE,plot.legend=TRUE,dot.scale=4)
 dev.off()
 postscript("./analysis/deg/Heatmap.DEG.eps",paper="special",width=10,height=5,horizontal=FALSE)
-plot(DoHeatmap(sc10x,genes.use=genes.deg.10,slim.col.label=TRUE,group.label.rot=TRUE,group.spacing=0.25,group.order=c("BE","LE","OE1","OE2","NE","Fib","SM","Endo","Leu"),cex.row=2.5))
+plot(DoHeatmap(sc10x,genes.use=genes.deg.10,slim.col.label=TRUE,group.label.rot=TRUE,group.spacing=0.25,cex.row=2.5))
 dev.off()
 
-for (i in ls(pattern="^genes.deg.*.unique")){
+for (i in ls(pattern="^genes.deg*unique")){
   postscript(paste0("./analysis/deg/Violin.",i,".eps"),paper="special",width=10,height=5,horizontal=FALSE)
   plot(VlnPlot(sc10x,features.plot=get(i)[1:10],nCol=5,point.size.use=0.1,size.title.use=15,x.lab.rot=TRUE))
   dev.off()
@@ -327,14 +328,27 @@ for (i in ls(pattern="^genes.deg.*.unique")){
   plot(RidgePlot(sc10x,features.plot=get(i)[1:10],nCol=5,size.x.use=10,size.title.use=15))
   dev.off()
   postscript(paste0("./analysis/deg/Heatmap.",i,".eps"),paper="special",width=10,height=5,horizontal=FALSE)
-  plot(DoHeatmap(sc10x,genes.use=get(i),slim.col.label=TRUE,group.label.rot=TRUE,group.spacing=0.25,group.order=c("BE","LE","OE1","OE2","NE","Fib","SM","Endo","Leu"),cex.row=2.5))
+  plot(DoHeatmap(sc10x,genes.use=get(i),slim.col.label=TRUE,group.label.rot=TRUE,group.spacing=0.25,cex.row=2.5))
   dev.off()
 }
 
-for (i in ls(pattern="^genes.deg.")){
+sc10x.Stress <- SetAllIdent(object=sc10x.Stress,id="Stress")
+postscript("./analysis/deg/Violin.Stress.eps",paper="special",width=10,height=5,horizontal=FALSE)
+plot(VlnPlot(sc10x.Stress,features.plot=rownames(genes.deg.Stress)[1:10],nCol=5,point.size.use=0.1,size.title.use=15,x.lab.rot=TRUE))
+dev.off()
+postscript("./analysis/deg/Ridge.Stress.eps",paper="special",width=10,height=5,horizontal=FALSE)
+plot(RidgePlot(sc10x.Stress,features.plot=rownames(genes.deg.Stress)[1:10],nCol=5,size.x.use=10,size.title.use=15))
+dev.off()
+postscript("./analysis/deg/Heatmap.Stress.eps",paper="special",width=10,height=5,horizontal=FALSE)
+plot(DoHeatmap(sc10x.Stress,genes.use=rownames(genes.deg.Stress),slim.col.label=TRUE,group.label.rot=TRUE,group.spacing=0.25,cex.row=2.5))
+dev.off()
+
+for (i in ls(pattern="^genes.deg")){
   write.table(get(i),file=paste0("./analysis/deg/",i,".csv"),row.names=TRUE,col.names=NA,append=FALSE,sep=",")
 }
 
+save(list=ls(pattern="^genes.deg"),file="./analysis/DEG.Rda")
+rm(ls(pattern="^genes.deg"))
 save(list=ls(pattern="^sc10x"),file="./analysis/sc10x.Rda")
 rm(ls(pattern="^sc10x"))
 save.image(file="./analysis/Data.RData")
\ No newline at end of file