1. Preparations before using plot1cell [or packages installing]
# turbo/fasten your R while connecting to network
options(BioC_mirror="https://mirrors.tuna.tsinghua.edu.cn/bioconductor")
options(repos=structure(c(CRAN="https://mirrors.tuna.tsinghua.edu.cn/CRAN/")))
dev.packages <- c("chris-mcginnis-ucsf/DoubletFinder","Novartis/hdf5r","mojaveazure/loomR")
devtools::install_github(dev.packages)
bioc.packages <- c("biomaRt","GenomeInfoDb","EnsDb.Hsapiens.v86","GEOquery","simplifyEnrichment","ComplexHeatmap")
BiocManager::install(bioc.packages)
devtools::install_github("TheHumphreysLab/plot1cell")
## or the development version, devtools::install_github("HaojiaWu/plot1cell")
2. Load packages and example data
library(plot1cell)
iri.integrated <- Install.example() # example data was a seurat object
# too slow
# searching "GSE139107" in geo and download yourself
# then reading data
# the folowed scripts are packaged in function "Install.example"
getGEOSuppFiles(GEO = "GSE139107") # not run/passed while download yourself
setwd("GSE139107/") # not run/passed while download yourself
all_files <- list.files(pattern = "dge")
all_ct <- list()
for (i in 1:length(all_files)) {
ct_data <- read.delim(all_files[i])
ct_data <- Matrix(as.matrix(ct_data), sparse = T)
all_ct[[i]] <- ct_data
print(all_files[i])
}
all.count <- RowMergeSparseMatrices(all_ct[[1]], all_ct[-1])
meta.data <- read.delim("GSE139107_MouseIRI.metadata.txt.gz")
all.count <- all.count[, rownames(meta.data)]
iri <- CreateSeuratObject(counts = all.count, min.cells = 0,
min.features = 0, meta.data = meta.data)
iri.list <- SplitObject(iri, split.by = "orig.ident")
iri.list <- lapply(X = iri.list, FUN = function(x) {
x <- NormalizeData(x, verbose = FALSE)
x <- FindVariableFeatures(x, verbose = FALSE)
})
features <- SelectIntegrationFeatures(object.list = iri.list)
iri.list <- lapply(X = iri.list, FUN = function(x) {
x <- ScaleData(x, features = features, verbose = FALSE)
x <- RunPCA(x, features = features, verbose = FALSE)
})
anchors <- FindIntegrationAnchors(object.list = iri.list,
reference = c(1, 2), reduction = "rpca", dims = 1:50) ## Fast integration using reciprocal PCA (RPCA)
iri.integrated <- IntegrateData(anchorset = anchors, dims = 1:50)
iri.integrated <- ScaleData(iri.integrated, verbose = FALSE)
iri.integrated <- RunPCA(iri.integrated, verbose = FALSE)
iri.integrated <- RunUMAP(iri.integrated, dims = 1:25, min.dist = 0.2)
iri.integrated <- SetIdent(iri.integrated, value = "celltype")
levels(iri.integrated) <- c("PTS1", "PTS2", "PTS3", "NewPT1",
"NewPT2", "DTL-ATL", "MTAL", "CTAL1", "CTAL2", "MD",
"DCT", "DCT-CNT", "CNT", "PC1", "PC2", "ICA", "ICB",
"Uro", "Pod", "PEC", "EC1", "EC2", "Fib", "Per", "Mø",
"Tcell")
levels(iri.integrated@meta.data$Group) <- c("Control", "4hours",
"12hours", "2days", "14days", "6weeks")
DefaultAssay(iri.integrated) <- "RNA"
colnames(iri.integrated@meta.data) # check colnames
head(iri.integrated) # check data
3. Using plot_circlize
to draw a circle plot
circ_data <- prepare_circlize_data(iri.integrated, scale = 0.8 ) # scaled for plot,[to make it bigger or smaller], has no relation to your data
set.seed(666666) # sed 666666 seems wonderfull
# 设置细胞分群信息的颜色
cluster_colors<-rand_color(length(levels(iri.integrated)))
group_colors<-rand_color(length(names(table(iri.integrated$Group))))
rep_colors<-rand_color(length(names(table(iri.integrated$orig.ident))))
### plot and save figures
# 绘制细胞分群圈图
pdf('circlize_plot.pdf', width = 6, height = 6)
plot_circlize(circ_data,do.label = T, pt.size = 0.1, col.use = cluster_colors ,bg.color = 'white', kde2d.n = 200, repel = T, label.cex = 1)
# 添加细胞群注释信息
add_track(circ_data, group = "Group", colors = group_colors, track_num = 2) ## can change it to one of the columns in the meta data of your seurat object
add_track(circ_data, group = "orig.ident",colors = rep_colors, track_num = 3) ## can change it to one of the columns in the meta data of your seurat object
dev.off()
image.png
4. Using complex_dotplot_single
to draw a dot plots for marker genes
png(filename = 'dotplot_single.png', width = 4, height = 6,units = 'in', res = 100) # pdf or png, two formats you shold learn
complex_dotplot_single(seu_obj = iri.integrated, feature = "Havcr1",groups = "Group")
dev.off()
image.png
设置groups和splitby参数对多个分组信息进行分割绘图
iri.integrated@meta.data$Phase<-plyr::mapvalues(iri.integrated@meta.data$Group, from = levels(iri.integrated@meta.data$Group), to = c("Healthy",rep("Injury",3), rep("Recovery",2)))
iri.integrated@meta.data$Phase<-as.character(iri.integrated@meta.data$Phase)
png(filename = 'dotplot_single_split.png', width = 4, height = 6,units = 'in', res = 100)
complex_dotplot_single(iri.integrated, feature = "Havcr1",groups = "Group",splitby = "Phase")
dev.off()
image.png
To visualize the same gene on multiple group factors, simply add more group factor IDs to the groups argument.
png(filename = 'dotplot_more_groups.png', width = 8, height = 6,units = 'in', res = 100)
complex_dotplot_single(seu_obj = iri.integrated, feature = "Havcr1",groups= c("Group","Replicates"))
dev.off()
image.png
Each group factor can be further splitted by its own factor if the splitby argument is provided. Note that in this case, the order of the group factors needs to match the order of splitby factors.
iri.integrated@meta.data$ReplicateID<-plyr::mapvalues(iri.integrated@meta.data$Replicates, from = names(table((iri.integrated@meta.data$Replicates))), to = c(rep("Rep1",3),rep("Rep2",3), rep("Rep3",1)))
iri.integrated@meta.data$ReplicateID<-as.character(iri.integrated@meta.data$ReplicateID)
png(filename = 'dotplot_more_groups_split.png', width = 9, height = 6,units = 'in', res = 200)
complex_dotplot_single(seu_obj = iri.integrated, feature = "Havcr1",groups= c("Group","Replicates"), splitby = c("Phase","ReplicateID"))
dev.off()
### In this example, "Phase" is a splitby factor for "Group" and "ReplicateID" is a splitby factor for "Replicates".
image.png
complex_dotplot_multiple
can be used too
To visualize multiple genes in dotplot format, complex_dotplot_multiple should be used.
png(filename = 'dotplot_multiple.png', width = 10, height = 4,units = 'in', res = 300)
complex_dotplot_multiple(seu_obj = iri.integrated, features = c("Slc34a1","Slc7a13","Havcr1","Krt20","Vcam1"),group = "Group", celltypes = c("PTS1" , "PTS2" , "PTS3" , "NewPT1" , "NewPT2"))
dev.off()
image.png
5. Violin plot to show gene expression across groups
png(filename = 'vlnplot_single.png', width = 4, height = 6,units = 'in', res = 100)
complex_vlnplot_single(iri.integrated, feature = "Havcr1", groups = "Group",celltypes = c("PTS1" , "PTS2" , "PTS3" , "NewPT1" , "NewPT2"))
dev.off()
image.png
Similar to
complex_dotplot_single
, thecomplex_vlnplot_single
function also allows splitting the group factor by another factor with the argument splitby.[i've to say, this is the best fuc for me]
png(filename = 'vlnplot_single_split.png', width = 4, height = 6,units = 'in', res = 100)
complex_vlnplot_single(iri.integrated, feature = "Havcr1", groups = "Group",celltypes = c("PTS1" , "PTS2" , "PTS3" , "NewPT1" , "NewPT2"), splitby = "Phase")
dev.off()
image.png
One gene/multiple group factors violin plot
[one gene per multi groups and for more than one cell types]
:
png(filename = 'vlnplot_multiple.png', width = 6, height = 6,units = 'in', res = 100)
complex_vlnplot_single(iri.integrated, feature = "Havcr1", groups = c("Group","Replicates"),celltypes = c("PTS1" , "PTS2" , "PTS3" , "NewPT1" , "NewPT2"), font.size = 10)
dev.off()
image.png
If we add another gene "PRMT5"
png(filename = 'vlnplot_multiple.png', width = 6, height = 6,units = 'in', res = 100)
complex_vlnplot_single(iri.integrated, feature = c("Havcr1","PRMT5"), groups = c("Group","Replicates"),celltypes = c("PTS1" , "PTS2" , "PTS3" , "NewPT1" , "NewPT2"), font.size = 10)
dev.off()
image.png
[to be continue...]
网友评论