Add sc_MergeSubClust.R

r.scripts/sc_MergeSubClust.R

+gc()
+.libPaths("/home2/ghenry/R/x86_64-pc-linux-gnu-library/3.3")
+library(optparse)
+library(Seurat)
+
+#retrive command line options
+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("--st",action="store",default=TRUE,type='logical',help="Analyze Data With Stress Pops Removed?"),
+  make_option("--r",action="store",default=0.1,type='numeric',help="Resolution"),
+  make_option("--e",action="store",default="dws",type='character',help="ID for Epi"),
+  make_option("--s",action="store",default="go",type='character',help="ID for St"),
+  make_option("--u",action="store",default="go",type='character',help="ID for UK")  
+)
+opt=parse_args(OptionParser(option_list=option_list))
+rm(option_list)
+Project.Name <- opt$p
+if (opt$st==TRUE){
+  sub.folder <- "/NoStress"
+  sub.file <- ".NOStress"
+} else {
+  sub.folder <- ""
+  sub.file <- ""
+}
+
+#load data
+setwd("../")
+load(paste0("./analysis/sc10x.",Project.Name,".Cumulative.RData"))
+load(paste0("./analysis/sc10x.",Project.Name,".",opt$g,".cluster",sub.file,".IDepi+st.Rda"))
+
+#create folder structure
+if (!dir.exists(paste0("./analysis/",opt$g))){
+  dir.create(paste0("./analysis/",opt$g))
+}
+if (!dir.exists(paste0("./analysis/",opt$g,"/global"))){
+  dir.create(paste0("./analysis/",opt$g,"/global"))
+}
+if (!dir.exists(paste0("./analysis/",opt$g,"/global",sub.folder))){
+  dir.create(paste0("./analysis/",opt$g,"/global",sub.folder))
+}
+
+sc10x.Group <- SetAllIdent(object=sc10x.Group,id=paste0("Epi.",opt$e))
+Clusters.Epi <- as.character( unique(sc10x.Group@ident))
+Clusters.Epi <- Clusters.St[!Clusters.St %in% c("St","?")]
+sc10x.Group <- SetAllIdent(object=sc10x.Group,id=paste0("St.",opt$s))
+Clusters.St <- as.character( unique(sc10x.Group@ident))
+Clusters.St <- Clusters.St[!Clusters.St %in% c("Epi","?")]
+for (i in Clusters.St){
+assign(paste0("clust.st.",i),names(sc10x.Group@ident[sc10x.Group@ident==i]))
+}
+if (file.exists(paste0("./analysis/sc10x.",Project.Name,".",opt$g,".UKSubcluster",sub.file,".IDst.Rda"))){
+  sc10x.Group <- SetAllIdent(object=sc10x.Group,id=paste0("UK.",opt$u))
+  Clusters.UK <- as.character( unique(sc10x.Group@ident))
+  Clusters.UK <- Clusters.UK[!Clusters.UK %in% c("Epi","St")]
+  for (i in Clusters.UK){
+    assign(paste0("clust.uk.",i),names(sc10x.Group@ident[sc10x.Group@ident==i]))
+}}
+sc10x.Group <- SetAllIdent(object=sc10x.Group,id=paste0("Epi.",opt$e))
+for (i in Clusters.St){
+  sc10x.Group <- SetIdent(object=sc10x.Group,cells.use=get(paste0("clust.st.",i)),ident.use=i)
+}
+if (file.exists(paste0("./analysis/sc10x.",Project.Name,".",opt$g,".UKSubcluster",sub.file,".IDst.Rda"))){
+  for (i in Clusters.UK){
+    sc10x.Group <- SetIdent(object=sc10x.Group,cells.use=get(paste0("clust.uk.",i)),ident.use=i)
+  }
+  sc10x.Group@ident <- factor(sc10x.Group@ident,levels=c(Clusters.Epi,Clusters.St,Clusters.UK))
+} else {
+  sc10x.Group@ident <- factor(sc10x.Group@ident,levels=c(Clusters.Epi,Clusters.St))
+}
+sc10x.Group <- StashIdent(object=sc10x.Group,save.name="merge")
+postscript(paste0("./analysis/",opt$g,"/global",sub.folder,"/tSNE_Merge.ID.eps"))
+plot <- TSNEPlot(object=sc10x.Group,pt.size=5,do.label=TRUE,label.size=10,do.return=TRUE)
+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))
+plot <- plot+guides(colour=guide_legend(override.aes=list(size=10)))
+plot(plot)
+dev.off()