Commit 393ed35d authored by Gervaise Henry's avatar Gervaise Henry 🤠
Browse files

Working initial version: works with Pd

parent 7f4b6797
gc()
library(methods)
library(optparse)
library(Seurat)
library(readr)
library(readxl)
library(fBasics)
library(pastecs)
library(qusage)
library(RColorBrewer)
library(monocle)
library(dplyr)
library(viridis)
library(gridExtra)
options(bitmapType="cairo")
setwd("../")
source("./r.scripts/sc-TissueMapper_functions.R")
source("./r.scripts/sc-TissueMapper_process.R")
scFolders()
sc10x <- scLoad(p="Pd")
lg=500
hg=2500
hm=0.1
results <- scQC(sc10x,lg=lg,hg=hg,hm=hm,sub=FALSE)
sc10x <- results[[1]]
counts.cell.raw <- results[[2]]
counts.gene.raw <- results[[3]]
counts.cell.filtered <- results[[4]]
counts.gene.filtered <- results[[5]]
rm(lg)
rm(hg)
rm(hm)
rm(results)
results <- scCellCycle(sc10x,sub=FALSE)
sc10x <- results[[1]]
genes.s <- results[[2]]
genes.g2m <- results[[3]]
rm(results)
save.image(file="./analysis/Data.cc.RData")
sc10x.l <- list()
sc10x.l[["D17"]] <- subset(sc10x,subset= D17=="D17")
sc10x.l[["D27"]] <- subset(sc10x,subset= D27=="D27")
sc10x.l[["D35"]] <- subset(sc10x,subset= D35=="D35")
sc10x <- scCCA(sc10x.l)
rm(sc10x.l)
save.image(file="./analysis/Data.cca.RData")
#gc()
sc10x <- ScaleData(object=sc10x,vars.to.regress=c("nFeature_RNA","percent.mito","S.Score","G2M.Score"),do.par=TRUE,num.cores=45,verbose=FALSE)
#sc10x <- ScaleData(object=sc10x,vars.to.regress=c("nFeature_RNA","percent.mito"),do.par=TRUE,num.cores=45,verbose=FALSE)
#gc()
results <- scPC(sc10x,pc=50,hpc=0.9,file="pre.stress",print="2",cca=TRUE)
sc10x <- results[[1]]
pc.use.prestress <- results[[2]]
rm(results)
save.image(file="./analysis/Data.prestresspc.RData")
sc10x <- scCluster(sc10x,res=0.5,red="pca",dim=pc.use.prestress,print="2",folder="pre.stress")
save.image(file="./analysis/Data.prestresscluster.RData")
genes.stress <- read_delim("./genesets/genes.deg.Stress.csv",",",escape_double=FALSE,trim_ws=TRUE,col_names=TRUE)
genes.stress <- genes.stress[1]
colnames(genes.stress) <- "scDWS.Stress"
results <- scScore(sc10x,score="Stress",geneset=as.list(genes.stress),cut.pt=0.9,anchor=c("EGR1","FOS","JUN"))
sc10x.preStress <- results[[1]]
sc10x <- results[[2]]
rm(results)
save.image(file="./analysis/Data.stress.RData")
#gc()
sc10x <- ScaleData(object=sc10x,vars.to.regress=c("nFeature_RNA","percent.mito","S.Score","G2M.Score"),do.par=TRUE,num.cores=45,verbose=FALSE)
#sc10x <- ScaleData(object=sc10x,vars.to.regress=c("nFeature_RNA","percent.mito"),do.par=TRUE,num.cores=45,verbose=FALSE)
#gc()
results <- scPC(sc10x,pc=50,hpc=0.9,file="post.stress",print="2",cca=TRUE)
sc10x <- results[[1]]
pc.use.poststress <- results[[2]]
rm(results)
save.image(file="./analysis/Data.poststresspc.RData")
res <- c(seq(0.1,0.5,0.1),0.75,seq(1,5,1))
sc10x <- scCluster(sc10x,res=res,red="pca",dim=pc.use.poststress,print="2",folder="ALL")
save.image(file="./analysis/Data.poststresscluster.RData")
# 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)
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.set,gene.set1)
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) <- "Club"
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) <- "Hillock"
gene.set <- c(gene.set,gene.set1)
gene.set1 <- read_delim("./genesets/genes.deg.NE.csv",",",escape_double=FALSE,trim_ws=TRUE,col_names=TRUE)
gene.set1 <- gene.set1[1]
gene.set1 <- as.list(gene.set1)
names(gene.set1) <- "NE"
gene.set <- c(gene.set,gene.set1)
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.set,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)
# genes.leu <- read_excel("./genesets/40425_2017_215_MOESM1_ESM.xlsx",sheet="S3. Candidate markers")
# leu <- as.list(unique(genes.leu[,2]))$Cell
# leu <- leu[-c(9,17,20,24,26,28)]
# #leu.l <- leu[-c(1,3:4,7:8,13:18,21,20:30)]
# #leu.lin <- c("Myeloid","Lymphoid","Lymphoid","Lymphoid","Myeloid","Myeloid","Lymphoid","Myeloid","Myeloid","Myeloid","Myeloid","Myeloid","Lymphoid","Lymphoid","Lymphoid","Myeloid","Lymphoid","Lymphoid","Lymphoid","Lymphoid","Lymphoid","Lymphoid","Lymphoid","Lymphoid","Lymphoid","Lymphoid","Lymphoid")
# #leu.lin <- c("Lymphoid","Myeloid","Myeloid","Myeloid","Myeloid","Myeloid","Myeloid","Lymphoid")
# genes.leu <- genes.leu[unlist(genes.leu[,2]) %in% unlist(leu),]
# #genes.leu.l <- genes.leu[unlist(genes.leu[,2]) %in% unlist(leu.l),]
# #genes.leu.l[nrow(genes.leu.l)+1,] <- list("MS4A7","Macrophages",0)
# #genes.leu.l[nrow(genes.leu.l)+1,] <- list("CD14","Macrophages",0)
# # genes.leu.l[nrow(genes.leu.l)+1,] <- list("CD86","Macrophages",0)
# # genes.leu.l[nrow(genes.leu.l)+1,] <- list("S100A9","Neutrophils",0)
# # genes.leu.l[nrow(genes.leu.l)+1,] <- list("EREG","Neutrophils",0)
# # genes.leu.l[nrow(genes.leu.l)+1,] <- list("S100A8","Neutrophils",0)
# # genes.leu.l[nrow(genes.leu.l)+1,] <- list("FCN1","Neutrophils",0)
# gene.set1 <- split(genes.leu[,1], genes.leu[,2])
# #gene.set1 <- split(genes.leu.l[,1], genes.leu.l[,2])
# gene.set1 <- lapply(gene.set1,unname)
# gene.set1 <- lapply(gene.set1,unlist)
# #names(gene.set1) <- leu.lin
# names(gene.set1) <- leu
# gene.set <- c(gene.set,gene.set1)
rm(gene.set1)
min <- min(table(sc10x$integrated_snn_res.0.5))
results <- scQuSAGE(sc10x,gs=gene.set,save=TRUE,type="lg",id="integrated_snn_res.0.5",ds=0,nm="all.pops",print="2")
sc10x <- results[[1]]
results.cor.all.pops <- results[[2]]
results.clust.all.pops.id <- results[[3]]
rm(results)
rm(gene.set)
\ No newline at end of file
This diff is collapsed.
#sc-TissueMapper
#Author: Gervaise H. Henry
#Email: gervaise.henry@utsouthwestern.edu
#Lab: Strand Lab, Deparment of Urology, University of Texas Southwestern Medical Center
density <- function(x,y,n=100){
dens <- MASS::kde2d(x=x,y=y,n=n)
ix <- findInterval(x,dens$x)
iy <- findInterval(y,dens$y)
ii <- cbind(ix,iy)
return(dens$z[ii]*max(x)*max(y))
}
\ No newline at end of file
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment