Change PC's to be top 80% of 25 calculated

r.scripts/sc_Cluster.R

@@ -9,12 +9,11 @@ library(readr)
 option_list=list(
   make_option("--p",action="store",default="Pr",type='character',help="Project Name"),
   make_option("--g",action="store",default="ALL",type='character',help="Group To Analyze"),
-  make_option("--pc",action="store",default=10,type='integer',help="Number Of Principal Components")
+  make_option("--pc",action="store",default=0,type='integer',help="Number Of Principal Components")
 )
 opt=parse_args(OptionParser(option_list=option_list))
 rm(option_list)
 Project.Name <- opt$p
-PC.Max <- opt$pc
 
 #load data
 setwd("../")
@@ -23,6 +22,12 @@ load(paste0("./analysis/sc10x.",Project.Name,".",opt$g,".qc.Rda"))
 sc10x.Group <- get(paste0("sc10x.",opt$g,".qc"))
 rm(list=paste0("sc10x.",opt$g,".qc"))
 
+if (opt$pc>0){
+  PC.Max <- opt$pc
+} else {
+  PC.Max <- PC.var.8
+}
+
 #create folder structure
 if (!dir.exists(paste0("./analysis/",opt$g))){
   dir.create(paste0("./analysis/",opt$g))

r.scripts/sc_LineageSubClust.R

@@ -14,7 +14,6 @@ option_list=list(
 opt=parse_args(OptionParser(option_list=option_list))
 rm(option_list)
 Project.Name <- opt$p
-PC.Max <- opt$pc
 if (opt$st==TRUE){
   sub.folder <- "/NoStress"
   sub.file <- ".NOStress"
@@ -58,23 +57,39 @@ sc10x.Group.Epi <- FindVariableGenes(object=sc10x.Group.Epi,mean.function=ExpMea
 sc10x.Group.St <- FindVariableGenes(object=sc10x.Group.St,mean.function=ExpMean,dispersion.function=LogVMR,x.low.cutoff=0.0125,x.high.cutoff=3,y.cutoff=0.5,do.plot=FALSE)
 
 #re-PCA on new variable genes
-sc10x.Group.Epi <- RunPCA(object=sc10x.Group.Epi,pc.genes=sc10x.Group.Epi@var.genes,do.print=FALSE,pcs.compute=50)
-sc10x.Group.St <- RunPCA(object=sc10x.Group.St,pc.genes=sc10x.Group.St@var.genes,do.print=FALSE,pcs.compute=50)
+sc10x.Group.Epi <- RunPCA(object=sc10x.Group.Epi,pc.genes=sc10x.Group.Epi@var.genes,do.print=FALSE,pcs.compute=25)
+sc10x.Group.St <- RunPCA(object=sc10x.Group.St,pc.genes=sc10x.Group.St@var.genes,do.print=FALSE,pcs.compute=25)
 
 #generate QC figures
 postscript(paste0("./analysis/",opt$g,"/qc/Plot_PCElbow.Epi.eps"))
-plot <- PCElbowPlot(object=sc10x.Group.Epi,num.pc=50)
+plot <- PCElbowPlot(object=sc10x.Group.Epi,num.pc=25)
 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))
 plot(plot)
 dev.off()
-
+PC.var.epi <- (sc10x.Group.Epi@dr$pca@sdev)^2
+PC.var.epi <- PC.var.epi/sum(PC.var.epi)
+PC.var.epi <- cumsum(PC.var.epi)
+PC.var.epi.8 <- min(which(PC.var.epi>=0.8))
 postscript(paste0("./analysis/",opt$g,"/qc/Plot_PCElbow.St.eps"))
-plot <- PCElbowPlot(object=sc10x.Group.St,num.pc=50)
+plot <- PCElbowPlot(object=sc10x.Group.St,num.pc=25)
 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))
 plot(plot)
 dev.off()
+PC.var.st <- (sc10x.Group.St@dr$pca@sdev)^2
+PC.var.st <- PC.var.st/sum(PC.var.st)
+PC.var.st <- cumsum(PC.var.st)
+PC.var.st.8 <- min(which(PC.var.st>=0.8))
+
+if (opt$pc>0){
+  PC.Max.epi <- opt$pc
+  PC.Max.st <- opt$pc
+} else {
+  PC.Max.epi <- PC.var.epi.8
+  PC.Max.st <- PC.var.st.8
+}
 
 #re-generate tSNE dimensions and generage sample figures
+PC.Max <- PC.Max.epi
 sc10x.Group.Epi <- RunTSNE(object=sc10x.Group.Epi,dims.use=1:PC.Max,do.fast=TRUE)
 postscript(paste0("./analysis/",opt$g,"/global",sub.folder,"/tSNE_Sample.Epi.eps"))
 plot <- TSNEPlot(object=sc10x.Group.Epi,group.by="samples",pt.size=2.5,do.return=TRUE,vector.friendly=FALSE)
@@ -83,6 +98,7 @@ plot <- plot+guides(colour=guide_legend(override.aes=list(size=10)))
 plot(plot)
 dev.off()
 
+PC.Max <- PC.Max.st
 sc10x.Group.St <- RunTSNE(object=sc10x.Group.St,dims.use=1:PC.Max,do.fast=TRUE)
 postscript(paste0("./analysis/",opt$g,"/global",sub.folder,"/tSNE_Sample.St.eps"))
 plot <- TSNEPlot(object=sc10x.Group.St,group.by="samples",pt.size=2.5,do.return=TRUE,vector.friendly=FALSE)
@@ -96,8 +112,10 @@ pops <- c("Epi","St")
 for (j in pops){
   if (j=="Epi"){
     sc10x.Group <- sc10x.Group.Epi
+    PC.Max <- PC.Max.epi
   } else if (j=="St") {
     sc10x.Group <- sc10x.Group.St
+    PC.Max <- PC.Max.st
   }
   for (i  in c(5,2,1,0.5,0.2,0.1)){
     gc()  
@@ -124,6 +142,8 @@ rm(j)
 
 #make generic variables of run specific variables
 rm(PC.Max)
+rm(PC.Max.epi)
+rm(PC.Max.st)
 opt.subclust <- opt
 rm(opt)
 

r.scripts/sc_QC.R

@@ -92,52 +92,35 @@ if (opt$cc==TRUE){
   gc()
 }
 sc10x.Group <- SetAllIdent(object=sc10x.Group,id="ALL")
-sc10x.Group <- RunPCA(object=sc10x.Group,pc.genes=sc10x.Group@var.genes,do.print=FALSE,pcs.compute=50)
-sc10x.Group <- ProjectPCA(object=sc10x.Group,do.print=FALSE,pcs.store=50)
+sc10x.Group <- RunPCA(object=sc10x.Group,pc.genes=sc10x.Group@var.genes,do.print=FALSE,pcs.compute=25)
+sc10x.Group <- ProjectPCA(object=sc10x.Group,do.print=FALSE,pcs.store=25)
 rm(mito.genes)
 rm(percent.mito)
 
 #generate qc figures
 postscript(paste0("./analysis/",opt$g,"/qc/Plot_PCElbow.eps"))
-plot <- PCElbowPlot(object=sc10x.Group,num.pc=50)
+plot <- PCElbowPlot(object=sc10x.Group,num.pc=25)
 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))
 plot(plot)
 dev.off()
-postscript(paste0("./analysis/",opt$g,"/qc/Plot_VizPCA.01.09.eps"))
-VizPCA(object=sc10x.Group,pcs.use=1:9)
+PC.var <- (sc10x.Group@dr$pca@sdev)^2
+PC.var <- PC.var/sum(PC.var)
+PC.var <- cumsum(PC.var)
+PC.var.8 <- min(which(PC.var>=0.8))
+postscript(paste0("./analysis/",opt$g,"/qc/Plot_VizPCA.01.05.eps"))
+VizPCA(object=sc10x.Group,pcs.use=1:5,num.genes=50,nCol=5)
 dev.off()
-postscript(paste0("./analysis/",opt$g,"/qc/Plot_VizPCA.10.18.eps"))
-VizPCA(object=sc10x.Group,pcs.use=10:18)
+postscript(paste0("./analysis/",opt$g,"/qc/Plot_VizPCA.06.10.eps"))
+VizPCA(object=sc10x.Group,pcs.use=6:10,num.genes=50,nCol=5)
 dev.off()
-postscript(paste0("./analysis/",opt$g,"/qc/Plot_VizPCA.19.27.eps"))
-VizPCA(object=sc10x.Group,pcs.use=19:27)
+postscript(paste0("./analysis/",opt$g,"/qc/Plot_VizPCA.11.15.eps"))
+VizPCA(object=sc10x.Group,pcs.use=11:15,num.genes=50,nCol=5)
 dev.off()
-postscript(paste0("./analysis/",opt$g,"/qc/Plot_VizPCA.28.36.eps"))
-VizPCA(object=sc10x.Group,pcs.use=28:36)
+postscript(paste0("./analysis/",opt$g,"/qc/Plot_VizPCA.16.20.eps"))
+VizPCA(object=sc10x.Group,pcs.use=16:20,num.genes=50,nCol=5)
 dev.off()
-postscript(paste0("./analysis/",opt$g,"/qc/Plot_VizPCA.37.45.eps"))
-VizPCA(object=sc10x.Group,pcs.use=37:45)
-dev.off()
-postscript(paste0("./analysis/",opt$g,"/qc/Plot_VizPCA.46.50.eps"))
-VizPCA(object=sc10x.Group,pcs.use=46:50)
-dev.off()
-postscript(paste0("./analysis/",opt$g,"/qc/Heatmap_PC.01.09.eps"))
-PCHeatmap(object=sc10x.Group,pc.use=1:9,cells.use=500,do.balanced=TRUE,label.columns=FALSE,use.full=FALSE)
-dev.off()
-postscript(paste0("./analysis/",opt$g,"/qc/Heatmap_PC.10.18.eps"))
-PCHeatmap(object=sc10x.Group,pc.use=10:18,cells.use=500,do.balanced=TRUE,label.columns=FALSE,use.full=FALSE)
-dev.off()
-postscript(paste0("./analysis/",opt$g,"/qc/Heatmap_PC.19.27.eps"))
-PCHeatmap(object=sc10x.Group,pc.use=19:27,cells.use=500,do.balanced=TRUE,label.columns=FALSE,use.full=FALSE)
-dev.off()
-postscript(paste0("./analysis/",opt$g,"/qc/Heatmap_PC.28.36.eps"))
-PCHeatmap(object=sc10x.Group,pc.use=28:36,cells.use=500,do.balanced=TRUE,label.columns=FALSE,use.full=FALSE)
-dev.off()
-postscript(paste0("./analysis/",opt$g,"/qc/Heatmap_PC.37.45.eps"))
-PCHeatmap(object=sc10x.Group,pc.use=37:45,cells.use=500,do.balanced=TRUE,label.columns=FALSE,use.full=FALSE)
-dev.off()
-postscript(paste0("./analysis/",opt$g,"/qc/Heatmap_PC.46.50.eps"))
-PCHeatmap(object=sc10x.Group,pc.use=46:50,cells.use=500,do.balanced=TRUE,label.columns=FALSE,use.full=FALSE)
+postscript(paste0("./analysis/",opt$g,"/qc/Plot_VizPCA.21.25.eps"))
+VizPCA(object=sc10x.Group,pcs.use=21:25,num.genes=50,nCol=5)
 dev.off()
 
 #make generic variables of run specific variables