R包封面
该R包由国外友人Enblacar完成，目前处于预印本阶段，旨在提供一种简化的方式，为已知的单细胞可视化生成可发布的图形。与“审美愉悦”一词一样主观，这是一组跨不同情节类型实施的主题修改。该软件包也可作为个人项目，具有未来的增长前景。
可以使用以下命令安装此软件包：
if(!requireNamespace("devtools", quietly = T)){
install.packages("devtools") # If not installed.
}
devtools::install_github("enblacar/SCpubr")
同时，为了使该包正常执行，需要安装下列包：
colortools
dplyr
enrichR
forcats
ggbeeswarm
ggplot2
ggpubr
ggrastr
ggrepel
Matrix
Nebulosa
patchwork
pbapply
purrr
rlang
scales
Seurat
stringr
svglite
tidyr
viridis
可以使用以下命令安装所有软件包：
cran_packages <- c("colortools",
"dplyr",
"enrichR",
"forcats",
"ggbeeswarm",
"ggplot2",
"ggpubr",
"ggrepel",
"Matrix",
"patchwork",
"purrr",
"rlang",
"scales",
"Seurat",
"stringr",
"svglite",
"tidyr",
"viridis")
install.packages(cran_packages)
bioconductor_packages <- c("Nebulosa")
if (!require("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install(bioconductor_packages)
seura对象的准备质控与降维
#引用R包
library(dplyr)
library(Seurat)
library(patchwork)
library(ggplot2)
library(SingleR)
library(CCA)
library(clustree)
library(cowplot)
library(monocle)
library(tidyverse)
library(SCpubr)
#设置目录并读取data文件夹下的10X文件
data_dir <- paste0(getwd(),"/data") #路径必须中文
samples=list.files(data_dir)
dir=file.path(data_dir,samples)
afdata <- Read10X(data.dir = dir)
# 简单创建一个seurat对象“sample”，每个feature至少在3细胞中表达同时每个细胞中至少200个feature被检测到
sample <- Seurat::CreateSeuratObject(counts = afdata,
project = "SeuratObject",
min.cells = 3,
min.features = 200)
# 计算一下线粒体与核糖体基因的百分比，在sample@meta.data添加列名"percent.mt"与"percent.rb"
sample <- Seurat::PercentageFeatureSet(sample, pattern = "^MT-", col.name = "percent.mt")
sample <- Seurat::PercentageFeatureSet(sample, pattern = "^RP", col.name = "percent.rb")
# 质控
mask1 <- sample$nCount_RNA >= 1000
mask2 <- sample$nFeature_RNA >= 500
mask3 <- sample$percent.mt <= 20
mask <- mask1 & mask2 & mask3
sample <- sample[, mask]
# 标准化
sample <- Seurat::SCTransform(sample)
# 简单降个维
sample <- Seurat::RunPCA(sample)
sample <- Seurat::RunUMAP(sample, dims = 1:30)
sample <- Seurat::RunTSNE(sample, dims = 1:30)
# 分一下cluster
sample <- Seurat::FindNeighbors(sample, dims = 1:30)
sample <- Seurat::FindClusters(sample, resolution = 1.2)
#简单注释一下
refdata=celldex::HumanPrimaryCellAtlasData()
afdata <- GetAssayData(sample, slot="data")
cellpred <- SingleR(test = afdata,
ref = refdata,
labels = refdata$label.main,
method = "cluster",
clusters = sample@meta.data$seurat_clusters)
metadata <- cellpred@metadata
celltype = data.frame(ClusterID = rownames(cellpred),
celltype = cellpred$labels,
stringsAsFactors = F)
#########细胞注释后结果可视化
newLabels=celltype$celltype
names(newLabels)=levels(sample)
sample=RenameIdents(sample, newLabels)
三种算法简单可视化一下：
#简单可视化一下
p1<-DimPlot(sample, reduction = "pca")
p2<-DimPlot(sample, reduction = "tsne")
p3<-DimPlot(sample, reduction = "umap")
p1+p2+p
image.png
增加一列celltype到sample@metada里：
#增加一列为celltype的对象
kk=as.data.frame(sample@active.ident)
colnames(kk)="celltype"
identical(colnames(sample),row.names(kk))
sample$celltype <-kk$celltype
identical(sample@meta.data[,"celltype"],kk[,"celltype"])
afmeta <- sample@meta.data
metada表头
使用SCpubr的do_DimPlot函数降维：
p4<-do_DimPlot(sample = sample,
plot.title = "af object",
reduction = "pca",
dims = c(2, 1)) #选择展示的主成分，这边是PC2与PC1
p5<-do_DimPlot(sample = sample,
plot.title = "af object",
reduction = "tsne",
dims = c(2, 1))
p6<-do_DimPlot(sample = sample,
plot.title = "af object",
reduction = "umap",
dims = c(2, 1))
p4+p5+p6
image.png
图片调整：
#更改图例位置置顶top并修改列数为2（放在左边则为left）
p4<-do_DimPlot(sample = sample,
plot.title = "My object", #标题
reduction = "pca",
legend.position = "top",
legend.ncol = 2,
dims = c(2, 1))
p5<-do_DimPlot(sample = sample,
plot.title = "My object",
reduction = "tsne",
legend.position = "top",
legend.ncol = 2,
dims = c(2, 1))
p6<-do_DimPlot(sample = sample,
plot.title = "My object",
reduction = "umap",
legend.position = "top",
legend.ncol = 2,
dims = c(2, 1))
p4+p5+p6
image.png
换个展示形式（细胞反倒没那么好看了）：
#使用标签而不是图例
p4<-do_DimPlot(sample = sample,
plot.title = "My object",
reduction = "pca",
legend.position = "top",
legend.ncol = 2,
dims = c(2, 1),
label = TRUE,
legend = FALSE)
p5<-do_DimPlot(sample = sample,
plot.title = "My object",
reduction = "tsne",
legend.position = "top",
legend.ncol = 2,
dims = c(2, 1),
label = TRUE,
legend = FALSE)
p6<-do_DimPlot(sample = sample,
plot.title = "My object",
reduction = "umap",
legend.position = "top",
legend.ncol = 2,
dims = c(2, 1),
label = TRUE,
legend = FALSE)
p4+p5+p6
image.png
换个颜色吧（配色鬼才）：
# colors.use换个颜色吧，红黄蓝绿橙紫（配色鬼才）
colors <- c("Chondrocytes" = "red",
"Endothelial_cells" = "yellow",
"Tissue_stem_cells" = "blue",
"Monocyte" = "green",
"T_cells" = "Orange",
"NK_cell" = "Purple"
)
p4<-do_DimPlot(sample = sample,
plot.title = "My object",
reduction = "pca",
legend.position = "top",
legend.ncol = 2,
dims = c(2, 1),
label = TRUE,
legend = FALSE,
colors.use = colors)
p5<-do_DimPlot(sample = sample,
plot.title = "My object",
reduction = "tsne",
legend.position = "top",
legend.ncol = 2,
dims = c(2, 1),
label = TRUE,
legend = FALSE,
colors.use = colors)
p6<-do_DimPlot(sample = sample,
plot.title = "My object",
reduction = "umap",
legend.position = "top",
legend.ncol = 2,
dims = c(2, 1),
label = TRUE,
legend = FALSE,
colors.use = colors)
p4+p5+p6
image.png
看看多个特征（如基因，线粒体含量等）的分布：
#单个特征基因的展示
do_FeaturePlot(sample,
features = "CD14",
plot.title = "CD14",
ncol = 1,
dims = c(2, 1))
#多个特征基因查询
do_FeaturePlot(sample,
features = c("percent.mt","percent.rb", "CD14"),
plot.title = "A collection of features",
ncol = 3,
dims = c(2, 1))
单个基因
多个基因或特征
**弱化某些细胞亚群的展示，即标成灰色，只看我们关注的细胞群之间的差异
#最后一行输入想弱化展示的细胞名
do_FeaturePlot(sample,
features = "CD14",
plot.title = "CD14",
ncol = 3,
dims = c(2, 1),
idents.highlight = levels(sample)[!(levels(sample) %in% c("Monocyte","Tissue_stem_cells","T_cells","NK_cell" ))])
image.png
换个配色：
do_FeaturePlot(sample,
features = "percent.mt",
plot.title = "percent.mt",
ncol = 1,
dims = c(2, 1),
viridis_color_map = "A")
SCpubr包 ****viridis_color_map****参数支持8种连续变量配色：
- A - 岩浆颜色图。
- B - 地狱色图。
- C - 等离子颜色图。
- D - viridis 颜色图。
- E - cividis 彩色地图。
- F——火箭彩图。
- G - mako 彩色地图。
-
H - 涡轮彩色图。image.png
小提琴图
do_VlnPlot(sample = sample,
features = "ACTB",
boxplot_width = 0.5) #箱子的宽度
image.png
换个配色（配色鬼才）
# colors.use换个颜色吧，红黄蓝绿橙紫（配色鬼才）
colors <- c("Chondrocytes" = "red",
"Endothelial_cells" = "yellow",
"Tissue_stem_cells" = "blue",
"Monocyte" = "green",
"T_cells" = "Orange",
"NK_cell" = "Purple"
)
do_VlnPlot(sample = sample,
features = "ACTB",
boxplot_width = 0.5,
colors.use = colors)
image.png
蜂群图
#蜂群图,看看某个基因的含量吧
do_BeeSwarmPlot(sample = sample,
feature_to_rank = "ACTB",
group.by = "celltype", #按什么分组比较
continuous_feature = T)
image.png
同样可以使用viridis_color_map参数更换魔鬼配色
#魔鬼配色
do_BeeSwarmPlot(sample = sample,
feature_to_rank = "ACTB",
group.by = "celltype",
continuous_feature = T, viridis_color_map = "A")
image.png
**画个点图吧家人们
#点图，将需要展示的基因复制为
genesgenes <- c("IL7R", "CCR7", "CD14", "LYZ", "S100A4", "MS4A1", "CD8A", "FCGR3A", "MS4A7", "GNLY", "NKG7", "FCER1A", "CST3", "PPBP")
do_DotPlot(sample = sample,
features = genes,
flip = T) #翻转坐标轴
image.png
当然也可以分区块进行细胞注释
#当然也可以分区块,这个时候genes就得是列表了
genes <- list("Chondrocytes" = c("IL7R", "CCR7"),
"Endothelial_cells" = c("CD14", "LYZ"),
"Monocyte" = c("CD8A"),
"T_cells" = c("FCGR3A", "MS4A7"),
"NK_cell" = c("GNLY", "NKG7"))
do_DotPlot(sample = sample,features = genes)
image.png
换个配色
#换个配色
do_DotPlot(sample = sample,
features = genes,
colors.use = c("blue","red"))
image.png
柱状图展示一下各细胞群的cluster含量
#柱状图
do_BarPlot(sample = sample,
features = "celltype",
group.by = "seurat_clusters",
legend = T,
horizontal = T,
add.summary_labels = T,
add.subgroup_labels = T,
repel.summary_labels = T,
repel.subgroup_labels = T)
image.png
只有一个样本展示不了样本中细胞亚群的比例了，这里放一下官方文档的示例：
image.png
单细胞的基因集富集展示：
#基因集富集
GS <- list("MAPK single" = c("MAPK1","MAPK2","MAPK4","IL1B","AKT","IL7R", "CCR7"),
"P53 single" = c("BAX","AKR", "BCL2","ACTB","TP53"),
"JAK-STAT single" = c("MAPK2","MAPK4","IL1B","AKT","IL7R", "CCR7"),
"cell cycle" = c("FCGR3A", "MS4A7","GNLY", "NKG7","CD8A"),
"Death" = c("BAX","AKR", "BCL2","MAPK1","MAPK2"))
do_EnrichmentHeatmap(sample = sample,
list_genes = GS,
group.by = "celltype",
transpose = T, #是否装置
column_names_rot = 90,row_names_rot = 0, #行列的宽度
cell_size = 5 #格子大小
)
image.png
关注”生信碱移“公众号后台回复：newafdata，即可获得示例文件与代码
说这么多了，快去实操一下吧
