diff --git a/r.scripts/sc-TissueMapper_RUN.D17.R b/r.scripts/sc-TissueMapper_RUN.D17.R index 1b639cf3b848292c60ee9dfd4bf2307ccd118db9..94c32717231b7068cc6a8939a3c5f69a891aea86 100644 --- a/r.scripts/sc-TissueMapper_RUN.D17.R +++ b/r.scripts/sc-TissueMapper_RUN.D17.R @@ -98,7 +98,7 @@ option_list=list( make_option("--stg",action="store",default="go",type='character',help="Geneset to use for stress ID"), #make_option("--hpc.poststress",action="store",default=0.85,type='numeric',help="Max variance cutoff for PCs to use, post-stress"), make_option("--res.poststress",action="store",default=0.2,type='numeric',help="Resolution to cluster, post-stress"), - make_option("--ds",action="store",default=0,type='integer',help="Number of cells to downsample"), + make_option("--ds",action="store",default=10000,type='integer',help="Number of cells to downsample"), #make_option("--hpc.epi",action="store",default=0.85,type='numeric',help="Max variance cutoff for PCs to use, Epi"), make_option("--res.epi",action="store",default=0.2,type='numeric',help="Resolution to cluster, Epi"), #make_option("--hpc.st",action="store",default=0.85,type='numeric',help="Max variance cutoff for PCs to use, St"), @@ -164,6 +164,7 @@ gene.set1 <- as.list(gene.set1) names(gene.set1) <- "Epi" gene.set <- c(gene.set,gene.set1) rm(gene.set1) +gc() results <- scQuSAGE(sc10x,gs=gene.set,res.use=opt$res.poststress,ds=opt$ds,nm="Lin",folder="lin") sc10x <- results[[1]] results.cor.Lin <- results[[2]] @@ -193,6 +194,7 @@ gene.set1 <- as.list(gene.set1) names(gene.set1) <- "OE" gene.set <- c(gene.set,gene.set1) rm(gene.set1) +gc() results <- scQuSAGE(sc10x.Epi,gs=gene.set,res.use=opt$res.epi,ds=opt$ds,nm="Epi.dws",folder="epi") sc10x.Epi <- results[[1]] results.cor.Epi.dws <- results[[2]] @@ -240,6 +242,7 @@ gene.set1 <- as.list(gene.set1[-1,]) names(gene.set1) <- "Leu" gene.set <- c(gene.set,gene.set1) rm(gene.set1) +gc() results <- scQuSAGE(sc10x.St,gs=gene.set,res.use=opt$res.st,ds=opt$ds,nm="St.go",folder="st") sc10x.St <- results[[1]] results.cor.St.go <- results[[2]] @@ -336,6 +339,14 @@ postscript("./analysis/deg/Heatmap.Stress.eps",paper="special",width=10,height=5 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() +sc10x.Epi.NE <- SetAllIdent(object=sc10x.Epi.NE,id="NE") +postscript("./analysis/deg/Violin.NE.eps",paper="special",width=10,height=5,horizontal=FALSE) +plot(VlnPlot(sc10x.Epi.NE,features.plot=rownames(genes.deg.NE)[1:10],nCol=5,point.size.use=0.1,size.title.use=15,x.lab.rot=TRUE)) +dev.off() +postscript("./analysis/deg/Ridge.NE.eps",paper="special",width=10,height=5,horizontal=FALSE) +plot(RidgePlot(sc10x.Epi.NE,features.plot=rownames(genes.deg.NE)[1:10],nCol=5,size.x.use=10,size.title.use=15)) +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=",") } diff --git a/r.scripts/sc-TissueMapper_RUN.D2.R b/r.scripts/sc-TissueMapper_RUN.D2.R index 66538ad01e0fb4eae777dcf7bcf54f0ffaee65ff..60600a5c8993950df4ed56f53f3b721f64a92db4 100644 --- a/r.scripts/sc-TissueMapper_RUN.D2.R +++ b/r.scripts/sc-TissueMapper_RUN.D2.R @@ -101,7 +101,7 @@ option_list=list( make_option("--stg",action="store",default="dws",type='character',help="Geneset to use for stress ID"), #make_option("--hpc.poststress",action="store",default=0.85,type='numeric',help="Max variance cutoff for PCs to use, post-stress"), make_option("--res.poststress",action="store",default=0.2,type='numeric',help="Resolution to cluster, post-stress"), - make_option("--ds",action="store",default=25000,type='integer',help="Number of cells to downsample"), + make_option("--ds",action="store",default=10000,type='integer',help="Number of cells to downsample"), #make_option("--hpc.epi",action="store",default=0.85,type='numeric',help="Max variance cutoff for PCs to use, Epi"), make_option("--res.epi",action="store",default=0.2,type='numeric',help="Resolution to cluster, Epi"), #make_option("--hpc.st",action="store",default=0.85,type='numeric',help="Max variance cutoff for PCs to use, St"), @@ -153,6 +153,7 @@ rm(results) rm(sc10x.D17) rm(sc10x.D27) +save.image(file="./analysis/Aligned.RData") #results <- scPC(sc10x,lx=opt$lx,hx=opt$hx,ly=opt$ly,cc=opt$cc,pc=opt$pc,hpc=opt$hpc,file="pre.stress",cca=TRUE) #sc10x <- results[[1]] @@ -160,11 +161,11 @@ rm(sc10x.D27) #pc.use <- results[[3]] #rm(results) -pc.use <- 25 +pc.use <- 15 sc10x <- scCluster(sc10x,pc.use=pc.use,res.use=opt$res.prestress,folder="pre.stress",red="cca.aligned") if (opt$st==TRUE){ - results <- scStress(sc10x,stg=opt$stg,res.use=opt$res.prestress,pc.use=pc.use) + results <- scStress(sc10x,stg=opt$stg,res.use=opt$res.prestress,pc.use=pc.use,cut=0.95) sc10x <- results[[1]] counts.cell.filtered.stress <- results[[2]] sc10x.Stress <- results[[3]] @@ -190,6 +191,7 @@ gene.set1 <- as.list(gene.set1) names(gene.set1) <- "Epi" gene.set <- c(gene.set,gene.set1) rm(gene.set1) +gc() results <- scQuSAGE(sc10x,gs=gene.set,res.use=opt$res.poststress,ds=opt$ds,nm="Lin",folder="lin") sc10x <- results[[1]] results.cor.Lin <- results[[2]] @@ -227,6 +229,7 @@ gene.set1 <- as.list(gene.set1) names(gene.set1) <- "OE_KRT13" gene.set <- c(gene.set,gene.set1) rm(gene.set1) +gc() results <- scQuSAGE(sc10x.Epi,gs=gene.set,res.use=opt$res.epi,ds=opt$ds,nm="Epi.dws.sc",folder="epi") sc10x.Epi <- results[[1]] results.cor.Epi.dws.sc <- results[[2]] @@ -250,6 +253,7 @@ gene.set1 <- as.list(gene.set1) names(gene.set1) <- "CGC" gene.set<- c(gene.set,gene.set1) rm(gene.set1) +gc() results.cor.Epi.lgea <- scQuSAGEsm(sc10x.Epi,gs=gene.set,ds=opt$ds,nm="Epi.dws.sc",folder="lgea") rm(gene.set) @@ -276,6 +280,7 @@ gene.set1 <- as.list(gene.set1) names(gene.set1) <- "Leu" gene.set <- c(gene.set,gene.set1) rm(gene.set1) +gc() results <- scQuSAGE(sc10x.St,gs=gene.set,res.use=opt$res.st,ds=opt$ds,nm="St.dws.sc",folder="st") sc10x.St <- results[[1]] results.cor.St.dws.sc <- results[[2]] @@ -283,7 +288,7 @@ results.clust.St.dws.sc.id <- results[[3]] rm(results) rm(gene.set) -sc10x.Epi.NE <- scNE(sc10x.Epi,neg="dws") +sc10x.Epi.NE <- scNE(sc10x.Epi,neg="dws",cut=0.95) sc10x.Epi <- scMergeSubClust(sc10x.Epi,i="Epi.dws.sc",g=c("BE","LE","OE_SCGB","OE_KRT13"),nm="Merge") @@ -372,6 +377,14 @@ postscript("./analysis/deg/Heatmap.Stress.eps",paper="special",width=10,height=5 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() +sc10x.Epi.NE <- SetAllIdent(object=sc10x.Epi.NE,id="NE") +postscript("./analysis/deg/Violin.NE.eps",paper="special",width=10,height=5,horizontal=FALSE) +plot(VlnPlot(sc10x.Epi.NE,features.plot=rownames(genes.deg.NE)[1:10],nCol=5,point.size.use=0.1,size.title.use=15,x.lab.rot=TRUE)) +dev.off() +postscript("./analysis/deg/Ridge.NE.eps",paper="special",width=10,height=5,horizontal=FALSE) +plot(RidgePlot(sc10x.Epi.NE,features.plot=rownames(genes.deg.NE)[1:10],nCol=5,size.x.use=10,size.title.use=15)) +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=",") } diff --git a/r.scripts/sc-TissueMapper_RUN.D3.R b/r.scripts/sc-TissueMapper_RUN.D3.R index 6e91a098bf0ecb0c0dba9fb10bb4e7fb7c06ad96..f06a90bc87ee94b20daad5fb7ccd8cfe89f9303b 100644 --- a/r.scripts/sc-TissueMapper_RUN.D3.R +++ b/r.scripts/sc-TissueMapper_RUN.D3.R @@ -101,7 +101,7 @@ option_list=list( make_option("--stg",action="store",default="dws",type='character',help="Geneset to use for stress ID"), #make_option("--hpc.poststress",action="store",default=0.85,type='numeric',help="Max variance cutoff for PCs to use, post-stress"), make_option("--res.poststress",action="store",default=0.2,type='numeric',help="Resolution to cluster, post-stress"), - make_option("--ds",action="store",default=25000,type='integer',help="Number of cells to downsample"), + make_option("--ds",action="store",default=10000,type='integer',help="Number of cells to downsample"), #make_option("--hpc.epi",action="store",default=0.85,type='numeric',help="Max variance cutoff for PCs to use, Epi"), make_option("--res.epi",action="store",default=0.2,type='numeric',help="Resolution to cluster, Epi"), #make_option("--hpc.st",action="store",default=0.85,type='numeric',help="Max variance cutoff for PCs to use, St"), @@ -168,6 +168,7 @@ rm(results) rm(sc10x.D17) rm(sc10x.D27) rm(sc10x.D35) +save.image(file="./analysis/Aligned.RData") #results <- scPC(sc10x,lx=opt$lx,hx=opt$hx,ly=opt$ly,cc=opt$cc,pc=opt$pc,hpc=opt$hpc,file="pre.stress",cca=TRUE) #sc10x <- results[[1]] @@ -175,11 +176,11 @@ rm(sc10x.D35) #pc.use <- results[[3]] #rm(results) -pc.use <- 25 +pc.use <- 15 sc10x <- scCluster(sc10x,pc.use=pc.use,res.use=opt$res.prestress,folder="pre.stress",red="cca.aligned") if (opt$st==TRUE){ - results <- scStress(sc10x,stg=opt$stg,res.use=opt$res.prestress,pc.use=pc.use) + results <- scStress(sc10x,stg=opt$stg,res.use=opt$res.prestress,pc.use=pc.use,cut=0.95) sc10x <- results[[1]] counts.cell.filtered.stress <- results[[2]] sc10x.Stress <- results[[3]] @@ -205,6 +206,7 @@ gene.set1 <- as.list(gene.set1) names(gene.set1) <- "Epi" gene.set <- c(gene.set,gene.set1) rm(gene.set1) +gc() results <- scQuSAGE(sc10x,gs=gene.set,res.use=opt$res.poststress,ds=opt$ds,nm="Lin",folder="lin") sc10x <- results[[1]] results.cor.Lin <- results[[2]] @@ -242,6 +244,7 @@ gene.set1 <- as.list(gene.set1) names(gene.set1) <- "OE_KRT13" gene.set <- c(gene.set,gene.set1) rm(gene.set1) +gc() results <- scQuSAGE(sc10x.Epi,gs=gene.set,res.use=opt$res.epi,ds=opt$ds,nm="Epi.dws.sc",folder="epi") sc10x.Epi <- results[[1]] results.cor.Epi.dws.sc <- results[[2]] @@ -265,6 +268,7 @@ gene.set1 <- as.list(gene.set1) names(gene.set1) <- "CGC" gene.set<- c(gene.set,gene.set1) rm(gene.set1) +gc() results.cor.Epi.lgea <- scQuSAGEsm(sc10x.Epi,gs=gene.set,ds=opt$ds,nm="Epi.dws.sc",folder="lgea") rm(gene.set) @@ -291,6 +295,7 @@ gene.set1 <- as.list(gene.set1) names(gene.set1) <- "Leu" gene.set <- c(gene.set,gene.set1) rm(gene.set1) +gc() results <- scQuSAGE(sc10x.St,gs=gene.set,res.use=opt$res.st,ds=opt$ds,nm="St.dws.sc",folder="st") sc10x.St <- results[[1]] results.cor.St.dws.sc <- results[[2]] @@ -298,7 +303,7 @@ results.clust.St.dws.sc.id <- results[[3]] rm(results) rm(gene.set) -sc10x.Epi.NE <- scNE(sc10x.Epi,neg="dws") +sc10x.Epi.NE <- scNE(sc10x.Epi,neg="dws"cut=0.95) sc10x.Epi <- scMergeSubClust(sc10x.Epi,i="Epi.dws.sc",g=c("BE","LE","OE_SCGB","OE_KRT13"),nm="Merge") @@ -387,6 +392,14 @@ postscript("./analysis/deg/Heatmap.Stress.eps",paper="special",width=10,height=5 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() +sc10x.Epi.NE <- SetAllIdent(object=sc10x.Epi.NE,id="NE") +postscript("./analysis/deg/Violin.NE.eps",paper="special",width=10,height=5,horizontal=FALSE) +plot(VlnPlot(sc10x.Epi.NE,features.plot=rownames(genes.deg.NE)[1:10],nCol=5,point.size.use=0.1,size.title.use=15,x.lab.rot=TRUE)) +dev.off() +postscript("./analysis/deg/Ridge.NE.eps",paper="special",width=10,height=5,horizontal=FALSE) +plot(RidgePlot(sc10x.Epi.NE,features.plot=rownames(genes.deg.NE)[1:10],nCol=5,size.x.use=10,size.title.use=15)) +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=",") }