Skip to content

GitLab

Commit e492acbc authored by Gervaise Henry's avatar Gervaise Henry 🤠
Browse files

Remove historical files and save analysis/DATA subdir

parent db82cafc
Expand all Hide whitespace changes
Inline Side-by-side
Showing with 11 additions and 3356 deletions
.gitignore
View file @ e492acbc
analysis/*
analysis/DATA/*
!analysis/DATA/
!.gitkeep
.vscode/
WR/
......
README.md
View file @ e492acbc
Determining cellular heterogeneity in the human prostate with single-cell RNA sequencing
========================================================================================
Generalized Code for the Analysis of scRNA-seq Data
===================================================
* Contact: **Gervaise H. Henry**
* <a href="https://orcid.org/0000-0001-7772-9578" target="orcid.widget" rel="noopener noreferrer" style="vertical-align:top;"><img src="https://orcid.org/sites/default/files/images/orcid_16x16.png" style="width:1em;margin-right:.5em;" alt="ORCID iD icon">orcid.org/0000-0001-7772-9578</a>
......@@ -10,18 +10,21 @@ Determining cellular heterogeneity in the human prostate with single-cell RNA se
* PI: Douglas W. Strand, PhD
* <a href="https://orcid.org/0000-0002-0746-927X" target="orcid.widget" rel="noopener noreferrer" style="vertical-align:top;"><img src="https://orcid.org/sites/default/files/images/orcid_16x16.png" style="width:1em;margin-right:.5em;" alt="ORCID iD icon">orcid.org/0000-0002-0746-927X</a>
* PI Email: [douglas.strand@utsouthwestern.edu](mailto:douglas.strand@utsouthwestern.edu)
<<<<<<< HEAD
=======
* **ANALYZED DATA FOR QUERYING AT: [StrandLab.net](http://strandlab.net/analysis.php)**
* **Raw data at: [GEO](https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE120716) & [GenitoUrinary Development Molecular Anatomy Project (GUDMAP)]("https://doi.org/10.25548/W-R8CM")**
* **Publication at:**
* Cell Reports: PENDING
* [BioRxiv](https://www.biorxiv.org/content/early/2018/10/15/439935)
>>>>>>> 0e5df71a00d5dd5e459450f677ac18bf7c7bfa1d
Data Analysis
-------------
* **Requirements:**
* /analysis/DATA/
* Pd-demultiplex.csv
* **ProjectName**-demultiplex.csv
* 10x/
* aggregation_csv.csv
* GRCh38/
......@@ -45,6 +48,8 @@ Data Analysis
* dplyr (v0.7.6)
* viridis (v0.5.1)
* *and all dependencies*
<<<<<<< HEAD
=======
* **HOW TO RUN**
* 1 Run on 3 patient aggregate
* run bash script [sc\_TissueMapper\-Pr.sh](https://git.biohpc.swmed.edu/StrandLab/sc-TissueMapper_Pr/blob/master/bash.scripts/sc_TissueMapper-Pd.sh)
......@@ -56,6 +61,7 @@ Data Analysis
* run bash script [sc\_TissueMapper\-DS\_D17.sh](https://git.biohpc.swmed.edu/StrandLab/sc-TissueMapper_Pr/blob/master/bash.scripts/sc_TissueMapper-DS_D17.sh)
* 5 Aggregate and compare several downsamples from # 4
* run r script [sc\_TissueMapper\_RUN.DS\_D17.aggr.R](https://git.biohpc.swmed.edu/StrandLab/sc-TissueMapper_Pr/blob/master/r.scripts/sc_TissueMapper_RUN.DS_D17.aggr.R)
>>>>>>> 0e5df71a00d5dd5e459450f677ac18bf7c7bfa1d
* **Pipeline:**
* Link cellranger count/aggr output to analysis
* Create demultiplex file to add custom sample groups
......
bash.scripts/sc_TissueMapper-D17_FACS.sh deleted 100644 → 0
View file @ db82cafc
#!/bin/bash
#SBATCH --job-name R_FullAnalysis.D17FACS
#SBATCH -p 256GB,256GBv1,384GB
#SBATCH -N 1
#SBATCH -t 7-0:0:0
#SBATCH -o job_%j.out
#SBATCH -e job_%j.out
#SBATCH --mail-type ALL
#SBATCH --mail-user gervaise.henry@utsouthwestern.edu
module load R/3.4.1-gccmkl
sh ./sc_LinkData.sh D17_FACS
Rscript ../r.scripts/sc-TissueMapper_RUN.D17_FACS.R
bash.scripts/sc_TissueMapper-D27_FACS.sh deleted 100644 → 0
View file @ db82cafc
#!/bin/bash
#SBATCH --job-name R_FullAnalysis.D27FACS
#SBATCH -p 256GB,256GBv1,384GB
#SBATCH -N 1
#SBATCH -t 7-0:0:0
#SBATCH -o job_%j.out
#SBATCH -e job_%j.out
#SBATCH --mail-type ALL
#SBATCH --mail-user gervaise.henry@utsouthwestern.edu
module load R/3.4.1-gccmkl
sh ./sc_LinkData.sh D27_FACS
Rscript ../r.scripts/sc-TissueMapper_RUN.D27_FACS.R
bash.scripts/sc_TissueMapper-DS_D17.sh deleted 100644 → 0
View file @ db82cafc
#!/bin/bash
#SBATCH --job-name R_FullAnalysis.DS
#SBATCH -p 256GB,256GBv1,384GB
#SBATCH -N 1
#SBATCH -t 7-0:0:0
#SBATCH -o job_%j.out
#SBATCH -e job_%j.out
#SBATCH --mail-type ALL
#SBATCH --mail-user gervaise.henry@utsouthwestern.edu
module load R/3.4.1-gccmkl
Rscript ../r.scripts/sc-TissueMapper_RUN.DS_D17.R
bash.scripts/sc_TissueMapper-Pd.sh deleted 100644 → 0
View file @ db82cafc
#!/bin/bash
#SBATCH --job-name R_FullAnalysis.Pd
#SBATCH -p 256GB,256GBv1,384GB
#SBATCH -N 1
#SBATCH -t 7-0:0:0
#SBATCH -o job_%j.out
#SBATCH -e job_%j.out
#SBATCH --mail-type ALL
#SBATCH --mail-user gervaise.henry@utsouthwestern.edu
module load R/3.4.1-gccmkl
sh ./sc_LinkData.sh Pd
Rscript ../r.scripts/sc-TissueMapper_RUN.Pd.R
Rscript ../r.scripts/sc-TissueMapper_RUN.Pd.StressCompare.R
r.scripts/sc-TissueMapper.R deleted 100644 → 0
View file @ db82cafc
This diff is collapsed.
r.scripts/sc-TissueMapper_CreateTestDataCI.R deleted 100644 → 0
View file @ db82cafc
gc()
.libPaths("/home2/ghenry/R/x86_64-pc-linux-gnu-library/3.4")
library(methods)
library(readr)
library(Matrix)
library(Seurat)
setwd("../")
genes <- read_delim("./analysis/DATA/genes.tsv","\t",escape_double=FALSE,col_names=FALSE,trim_ws=TRUE)
barcodes <- read_delim("./analysis/DATA/barcodes.tsv","\t",escape_double=FALSE,col_names=FALSE,trim_ws=TRUE)
matrix <- readMM("./analysis/DATA/matrix.mtx")
load("./analysis/DATA/sc10x.data.Rda")
sc10x.Group <- SetAllIdent(object=sc10x.Group,id="SubClust.ID+NE")
barcodes.NE <- names(sc10x.Group@ident[sc10x.Group@ident=="NE"])
genes.rnd <- 1:nrow(genes)
#genes.rnd <- sample(genes.rnd,2500)
barcodes.rnd <- 1:nrow(barcodes)
barcodes.NE.index <- c()
for (i in 1:length(barcodes.NE)){
barcodes.NE.index <- c(barcodes.NE.index,which(barcodes==barcodes.NE[i]))
}
barcodes.rnd <- barcodes.rnd[!barcodes.rnd %in% barcodes.NE.index]
barcodes.rnd <- sample(barcodes.rnd,(2500-length(barcodes.NE)))
barcodes.rnd <- c(barcodes.rnd,barcodes.NE.index)
genes.sample <- genes[genes.rnd,]
barcodes.sample <- barcodes[barcodes.rnd,]
matrix.sample <- matrix[genes.rnd,barcodes.rnd]
if (!dir.exists("./analysis/CI.TestData.downsample")){
dir.create("./analysis/CI.TestData.downsample")
}
if (!dir.exists("./analysis/CI.TestData.downsample/10x")){
dir.create("./analysis/CI.TestData.downsample/10x")
}
if (!dir.exists("./analysis/CI.TestData.downsample/10x/filtered_gene_bc_matrices_mex")){
dir.create("./analysis/CI.TestData.downsample/10x/filtered_gene_bc_matrices_mex")
}
if (!dir.exists("./analysis/CI.TestData.downsample/10x/filtered_gene_bc_matrices_mex/GRCh38")){
dir.create("./analysis/CI.TestData.downsample/10x/filtered_gene_bc_matrices_mex/GRCh38")
}
write.table(genes.sample,file="./analysis/CI.TestData.downsample/10x/filtered_gene_bc_matrices_mex/GRCh38/genes.tsv",row.names=FALSE,col.names=FALSE,append=FALSE,sep="\t",quote=FALSE)
write.table(barcodes.sample,file="./analysis/CI.TestData.downsample/10x/filtered_gene_bc_matrices_mex/GRCh38/barcodes.tsv",row.names=FALSE,col.names=FALSE,append=FALSE,sep="\t",quote=FALSE)
writeMM(matrix.sample,file="./analysis/CI.TestData.downsample/10x/filtered_gene_bc_matrices_mex/GRCh38/matrix.mtx")
r.scripts/sc-TissueMapper_RUN.D17_FACS.R deleted 100644 → 0
View file @ db82cafc
gc()
library(methods)
library(optparse)
library(Seurat)
library(readr)
library(fBasics)
library(pastecs)
library(qusage)
library(RColorBrewer)
library(monocle)
library(dplyr)
library(viridis)
source("../r.scripts/sc-TissueMapper.R")
#Create folder structure
setwd("../")
if (!dir.exists("./analysis")){
dir.create("./analysis")
}
if (!dir.exists("./analysis/qc")){
dir.create("./analysis/qc")
}
if (!dir.exists("./analysis/qc/cc")){
dir.create("./analysis/qc/cc")
}
if (!dir.exists("./analysis/tSNE")){
dir.create("./analysis/tSNE")
}
if (!dir.exists("./analysis/tSNE/pre.stress")){
dir.create("./analysis/tSNE/pre.stress")
}
if (!dir.exists("./analysis/pca")){
dir.create("./analysis/pca")
}
if (!dir.exists("./analysis/pca/stress")){
dir.create("./analysis/pca/stress")
}
if (!dir.exists("./analysis/violin")){
dir.create("./analysis/violin")
}
if (!dir.exists("./analysis/violin/stress")){
dir.create("./analysis/violin/stress")
}
if (!dir.exists("./analysis/table")){
dir.create("./analysis/table")
}
if (!dir.exists("./analysis/tSNE/post.stress")){
dir.create("./analysis/tSNE/post.stress")
}
if (!dir.exists("./analysis/cor")){
dir.create("./analysis/cor")
}
if (!dir.exists("./analysis/tSNE/lin")){
dir.create("./analysis/tSNE/lin")
}
if (!dir.exists("./analysis/tSNE/epi")){
dir.create("./analysis/tSNE/epi")
}
if (!dir.exists("./analysis/tSNE/st")){
dir.create("./analysis/tSNE/st")
}
if (!dir.exists("./analysis/tSNE/merge")){
dir.create("./analysis/tSNE/merge")
}
if (!dir.exists("./analysis/pca/ne")){
dir.create("./analysis/pca/ne")
}
if (!dir.exists("./analysis/tSNE/ne")){
dir.create("./analysis/tSNE/ne")
}
if (!dir.exists("./analysis/violin/ne")){
dir.create("./analysis/violin/ne")
}
if (!dir.exists("./analysis/tSNE/FINAL")){
dir.create("./analysis/tSNE/FINAL")
}
if (!dir.exists("./analysis/deg")){
dir.create("./analysis/deg")
}
if (!dir.exists("./analysis/cca")){
dir.create("./analysis/cca")
}
if (!dir.exists("./analysis/diy")){
dir.create("./analysis/diy")
}
if (!dir.exists("./analysis/pseudotime")){
dir.create("./analysis/pseudotime")
}
#Retrieve command-line options
option_list=list(
make_option("--p",action="store",default="DPrF",type='character',help="Project Name"),
make_option("--g",action="store",default="ALL",type='character',help="Group To analyze"),
make_option("--lg",action="store",default=500,type='integer',help="Threshold for cells with minimum genes"),
make_option("--hg",action="store",default=3000,type='integer',help="Threshold for cells with maximum genes"),
make_option("--lm",action="store",default=0,type='numeric',help="Threshold for cells with minimum %mito genes"),
make_option("--hm",action="store",default=0.1,type='numeric',help="Threshold for cells with maximum %mito genes"),
make_option("--lx",action="store",default=0.2,type='numeric',help="x low threshold for hvg selection"),
make_option("--hx",action="store",default=5,type='numeric',help="x high threshold for hvg selection"),
make_option("--ly",action="store",default=1,type='numeric',help="y low threshold for hvg selection"),
make_option("--cc",action="store",default=TRUE,type='logical',help="Scale cell cycle?"),
make_option("--cca",action="store",default=50,type='integer',help="Number of CCAs to cacluate"),
make_option("--acca",action="store",default=30,type='integer',help="Number of CCAs to align"),
make_option("--pc",action="store",default=50,type='integer',help="Number of PCs to cacluate"),
make_option("--res.prestress",action="store",default=1,type='numeric',help="Resolution to cluster, pre-stress"),
make_option("--st",action="store",default=TRUE,type='logical',help="Remove stressed cells?"),
make_option("--stg",action="store",default="dws",type='character',help="Geneset to use for stress ID"),
make_option("--cut.stress",action="store",default=0.9,type='numeric',help="Cutoff for stress score"),
make_option("--res.poststress",action="store",default=0.5,type='numeric',help="Resolution to cluster, post-stress"),
make_option("--cut.ne",action="store",default=0.999,type='numeric',help="Cutoff for NE score")
)
opt=parse_args(OptionParser(option_list=option_list))
rm(option_list)
if (opt$lm==0){opt$lm=-Inf}
sc10x <- scLoad("D17_FACS")
if (opt$cc==TRUE){
results <- scCellCycle(sc10x)
sc10x <- results[[1]]
genes.s <- results[[2]]
genes.g2m <- results[[3]]
rm(results)
} else {
genes.s=""
genes.g2m=""
}
results <- scQC(sc10x,lg=opt$lg,hg=opt$hg,lm=opt$lm,hm=opt$hm)
sc10x <- results[[1]]
counts.cell.raw <- results[[2]]
counts.gene.raw <- results[[3]]
counts.cell.filtered <- results[[4]]
counts.gene.filtered <- results[[5]]
rm(results)
gc()
if (opt$cc==TRUE){
sc10x <- ScaleData(object=sc10x,vars.to.regress=c("nUMI","percent.mito","S.Score","G2M.Score"),display.progress=FALSE,do.par=TRUE,num.cores=45)
} else {
sc10x <- ScaleData(object=sc10x,vars.to.regress=c("nUMI","percent.mito"),display.progress=FALSE,do.par=TRUE,num.cores=45)
}
gc()
results <- scPC(sc10x,lx=opt$lx,hx=opt$hx,ly=opt$ly,cc=opt$cc,pc=50,hpc=0.85,file="pre.stress",cca=FALSE)
sc10x <- results[[1]]
genes.hvg.prestress <- results[[2]]
pc.use.prestress <- results[[3]]
rm(results)
sc10x <- scCluster(sc10x,pc.use=pc.use.prestress,res.use=opt$res.prestress,folder="pre.stress",red="pca")
if (opt$st==TRUE){
results <- scStress(sc10x,stg=opt$stg,res.use=opt$res.prestress,cut=opt$cut.stress)
sc10x <- results[[1]]
counts.cell.filtered.stress <- results[[2]]
sc10x.Stress <- results[[3]]
rm(results)
results <- scPC(sc10x,lx=opt$lx,hx=opt$hx,ly=opt$ly,cc=opt$cc,pc=50,hpc=0.85,file="post.stress",cca=FALSE)
sc10x <- results[[1]]
genes.hvg.poststress <- results[[2]]
pc.use.poststress <- results[[3]]
rm(results)
sc10x <- scCluster(sc10x,pc.use=pc.use.poststress,res.use=0.5,folder="post.stress",red="pca")
}
gene.set1 <- read_delim("./genesets/DEG_Epi_5FC.txt","\t",escape_double=FALSE,trim_ws=TRUE,col_names=FALSE)
gene.set1 <- as.list(gene.set1)
names(gene.set1) <- "Epi"
gene.set <- c(gene.set1)
gene.set1 <- read_delim("./genesets/DEG_FMSt_5FC.txt","\t",escape_double=FALSE,trim_ws=TRUE,col_names=FALSE)
gene.set1 <- as.list(gene.set1)
names(gene.set1) <- "St"
gene.set <- c(gene.set,gene.set1)
rm(gene.set1)
gc()
min.all <- min(table(sc10x@meta.data[,paste0("res",opt$res.poststress)]))
results <- scQuSAGE(sc10x,gs=gene.set,res.use=opt$res.poststress,ds=min.all,nm="Lin",folder="lin")
sc10x <- results[[1]]
results.cor.Lin <- results[[2]]
results.clust.Lin.id <- results[[3]]
rm(results)
rm(gene.set)
sc10x <- SetAllIdent(object=sc10x,id="Lin")
sc10x.Epi <- scSubset(sc10x,i="Lin",g="Epi")
if (any(levels(sc10x@ident)=="Unknown")){
sc10x.St <- scSubset(sc10x,i="Lin",g=c("St","Unknown"))
} else {
sc10x.St <- scSubset(sc10x,i="Lin",g="St")
}
sc10x.Epi <- SetAllIdent(object=sc10x.Epi,id=paste0("res",opt$res.poststress))
sc10x.Epi <- BuildClusterTree(sc10x.Epi,do.reorder=TRUE,reorder.numeric=TRUE,do.plot=FALSE)
sc10x.Epi <- StashIdent(object=sc10x.Epi,save.name=paste0("res",opt$res.poststress))
sc10x.St <- SetAllIdent(object=sc10x.St,id=paste0("res",opt$res.poststress))
sc10x.St <- BuildClusterTree(sc10x.St,do.reorder=TRUE,reorder.numeric=TRUE,do.plot=FALSE)
sc10x.St <- StashIdent(object=sc10x.St,save.name=paste0("res",opt$res.poststress))
sc10x.Epi <- RunTSNE(object=sc10x.Epi,reduction.use="pca",dims.use=1:pc.use.poststress,do.fast=TRUE)
postscript(paste0("./analysis/tSNE/epi/tSNE_Sample.eps"))
plot <- TSNEPlot(object=sc10x.Epi,group.by="samples",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))
plot <- plot+guides(colour=guide_legend(override.aes=list(size=10)))
plot(plot)
dev.off()
postscript(paste0("./analysis/tSNE/epi/tSNE_res",opt$res.poststress,".eps"))
plot <- TSNEPlot(object=sc10x.Epi,pt.size=5,do.label=TRUE,label.size=10,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))
plot <- plot+guides(colour=guide_legend(override.aes=list(size=10)))
plot(plot)
dev.off()
rm(plot)
sc10x.St <- RunTSNE(object=sc10x.St,reduction.use="pca",dims.use=1:pc.use.poststress,do.fast=TRUE)
postscript(paste0("./analysis/tSNE/st/tSNE_Sample.eps"))
plot <- TSNEPlot(object=sc10x.St,group.by="samples",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))
plot <- plot+guides(colour=guide_legend(override.aes=list(size=10)))
plot(plot)
dev.off()
postscript(paste0("./analysis/tSNE/st/tSNE_res",opt$res.poststress,".eps"))
plot <- TSNEPlot(object=sc10x.St,pt.size=5,do.label=TRUE,label.size=10,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))
plot <- plot+guides(colour=guide_legend(override.aes=list(size=10)))
plot(plot)
dev.off()
rm(plot)
gene.set1 <- read_delim("./genesets/genes.deg.BE.csv",",",escape_double=FALSE,trim_ws=TRUE,col_names=TRUE)
gene.set1 <- gene.set1[1]
gene.set1 <- as.list(gene.set1)
names(gene.set1) <- "BE"
gene.set <- c(gene.set1)
gene.set1 <- read_delim("./genesets/genes.deg.LE.csv",",",escape_double=FALSE,trim_ws=TRUE,col_names=TRUE)
gene.set1 <- gene.set1[1]
gene.set1 <- as.list(gene.set1)
names(gene.set1) <- "LE"
gene.set <- c(gene.set,gene.set1)
gene.set1 <- read_delim("./genesets/genes.deg.OE1.csv",",",escape_double=FALSE,trim_ws=TRUE,col_names=TRUE)
gene.set1 <- gene.set1[1]
gene.set1 <- as.list(gene.set1)
names(gene.set1) <- "OE_SCGB"
gene.set <- c(gene.set,gene.set1)
gene.set1 <- read_delim("./genesets/genes.deg.OE2.csv",",",escape_double=FALSE,trim_ws=TRUE,col_names=TRUE)
gene.set1 <- gene.set1[1]
gene.set1 <- as.list(gene.set1)
names(gene.set1) <- "OE_KRT13"
gene.set <- c(gene.set,gene.set1)
rm(gene.set1)
gc()
min.epi <- min(table(sc10x.Epi@meta.data[,paste0("res",opt$res.poststress)]))
results <- scQuSAGE(sc10x.Epi,gs=gene.set,res.use=opt$res.poststress,ds=min.epi,nm="Epi.dws.sc",folder="epi")
sc10x.Epi <- results[[1]]
results.cor.Epi.dws <- results[[2]]
results.clust.Epi.dws.id <- results[[3]]
rm(results)
rm(gene.set)
gene.set1 <- read_delim("./genesets/genes.deg.Endo.csv",",",escape_double=FALSE,trim_ws=TRUE,col_names=TRUE)
gene.set1 <- gene.set1[1]
gene.set1 <- as.list(gene.set1)
names(gene.set1) <- "Endo"
gene.set <- c(gene.set1)
gene.set1 <- read_delim("./genesets/genes.deg.SM.csv",",",escape_double=FALSE,trim_ws=TRUE,col_names=TRUE)
gene.set1 <- gene.set1[1]
gene.set1 <- as.list(gene.set1)
names(gene.set1) <- "SM"
gene.set <- c(gene.set,gene.set1)
gene.set1 <- read_delim("./genesets/genes.deg.Fib.csv",",",escape_double=FALSE,trim_ws=TRUE,col_names=TRUE)
gene.set1 <- gene.set1[1]
gene.set1 <- as.list(gene.set1)
names(gene.set1) <- "Fib"
gene.set <- c(gene.set,gene.set1)
gene.set1 <- read_delim("./genesets/genes.deg.Leu.csv",",",escape_double=FALSE,trim_ws=TRUE,col_names=TRUE)
gene.set1 <- gene.set1[1]
gene.set1 <- as.list(gene.set1)
names(gene.set1) <- "Leu"
gene.set <- c(gene.set,gene.set1)
rm(gene.set1)
gc()
min.st <- min(table(sc10x.St@meta.data[,paste0("res",opt$res.poststress)]))
results <- scQuSAGE(sc10x.St,gs=gene.set,res.use=opt$res.poststress,ds=min.st,nm="St.dws.sc",folder="st")
sc10x.St <- results[[1]]
results.cor.St.go <- results[[2]]
results.clust.St.go.id <- results[[3]]
rm(results)
rm(gene.set)
sc10x.Epi.NE <- scNE(sc10x.Epi,neg="dws",cut=opt$cut.ne)
sc10x <- scMerge(sc10x,sc10x.Epi,sc10x.St,i.1="Epi.dws.sc",i.2="St.dws.sc",nm="Merge_Epi.dws.sc_St.dws.sc")
sc10x <- scMerge(sc10x,sc10x,sc10x.Epi.NE,i.1="Merge_Epi.dws.sc_St.dws.sc",i.2="NE",nm="Merge_Epi.dws.sc_St.dws.sc_NE")
sc10x.Epi <- scMerge(sc10x.Epi,sc10x.Epi,sc10x.Epi.NE,i.1="Epi.dws.sc",i.2="NE",nm="Epi.dws.sc_NE")
sc10x <- SetAllIdent(object=sc10x,id="Merge_Epi.dws.sc_St.dws.sc")
sc10x@ident <- factor(sc10x@ident,levels=c("BE","LE","OE_SCGB","OE_KRT13","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))
plot <- plot+guides(colour=guide_legend(override.aes=list(size=10)))
plot(plot)
dev.off()
scTables(sc10x,i.1="samples",i.2="Merge_Epi.dws.sc_St.dws.sc")
scTables(sc10x,i.1="samples",i.2="Merge_Epi.dws.sc_St.dws.sc_NE")
scTables(sc10x,i.1="Merge_Epi.dws.sc_St.dws.sc_NE",i.2="Merge_Epi.dws.sc_St.dws.sc")
sctSNECustCol(sc10x,i="Lin",bl="Epi",rd="St",file="D17")
sctSNECustCol(sc10x,i="Merge_Epi.dws.sc_St.dws.sc",bl=c("BE","LE","OE_SCGB","OE_KRT13"),rd=c("Fib","SM","Endo","Leu"),file="D17")
sctSNECustCol(sc10x.Epi,i="Epi.dws.sc",bl=c("BE","LE","OE_SCGB","OE_KRT13"),rd="",file="D17")
sctSNECustCol(sc10x.St,i="St.dws.sc",bl="",rd=c("Fib","SM","Endo","Leu"),file="D17")
sctSNEbwCol(sc10x,i=paste0("res",opt$res.poststress),file="ALL",files="D17")
sctSNEbwCol(sc10x.Epi,i=paste0("res",opt$res.poststress),file="Epi",files="D17")
sctSNEbwCol(sc10x.St,i=paste0("res",opt$res.poststress),file="St",files="D17")
sctSNEbwCol(sc10x,i="Merge_Epi.dws.sc_St.dws.sc",file="ALL",files="D17")
sctSNEbwCol(sc10x.Epi,i="Epi.dws.sc",file="Epi",files="D17")
sctSNEbwCol(sc10x.St,i="St.dws.sc",file="St",files="D17")
for (g in c("Epi","St","Unknown")){
sctSNEHighlight(sc10x,i="Lin",g=g,file="D17")
}
for (g in c("BE","LE","OE_SCGB","OE_KRT13")){
sctSNEHighlight(sc10x,i="Merge_Epi.dws.sc_St.dws.sc",g=g,file="D17")
sctSNEHighlight(sc10x.Epi,i="Epi.dws.c",g=g,file="D17")
}
sctSNEHighlight(sc10x.Epi.NE,i="NE",g="NE",file="D17")
for (g in c("Fib","SM","Endo","Leu")){
sctSNEHighlight(sc10x,i="Merge_Epi.dws.sc_St.dws.sc",g=g,file="D17")
sctSNEHighlight(sc10x.St,i="St.dws.sc",g=g,file="D17")
}
for (g in c("D17PrPzF_BE","D17PrPzF_LE","D17PrPzF_OE","D17PrPzF_FMSt")){
sctSNEHighlight(sc10x,i="samples",g=g,file="D17")
}
rm(i)
rm(g)
postscript(paste0("./analysis/diy/Feature.eps"))
FeaturePlot(sc10x,features.plot=c("PECAM1","CD200","PDPN","EPCAM","DPP4","PSCA","VIM","KRT5","KLK3"),cols.use=c("grey","darkred"),reduction.use="tsne")
dev.off()
save(list=ls(pattern="sc10x.Stress"),file="./analysis/sc10x.Stress.Rda")
rm(list=ls(pattern="sc10x.Stress"))
save(list=ls(pattern="sc10x.Epi"),file="./analysis/sc10x.Epi.Rda")
rm(list=ls(pattern="^sc10x.Epi"))
save(list=ls(pattern="sc10x.St"),file="./analysis/sc10x.St.Rda")
rm(list=ls(pattern="sc10x.St"))
save(list=ls(pattern="^sc10x"),file="./analysis/sc10x.Rda")
rm(list=ls(pattern="^sc10x"))
save.image(file="./analysis/Data.RData")
r.scripts/sc-TissueMapper_RUN.D27_FACS.R deleted 100644 → 0
View file @ db82cafc
gc()
library(methods)
library(optparse)
library(Seurat)
library(readr)
library(fBasics)
library(pastecs)
library(qusage)
library(RColorBrewer)
library(monocle)
library(dplyr)
library(viridis)
source("../r.scripts/sc-TissueMapper.R")
#Create folder structure
setwd("../")
if (!dir.exists("./analysis")){
dir.create("./analysis")
}
if (!dir.exists("./analysis/qc")){
dir.create("./analysis/qc")
}
if (!dir.exists("./analysis/qc/cc")){
dir.create("./analysis/qc/cc")
}
if (!dir.exists("./analysis/tSNE")){
dir.create("./analysis/tSNE")
}
if (!dir.exists("./analysis/tSNE/pre.stress")){
dir.create("./analysis/tSNE/pre.stress")