SCpubr包：让你的单细胞测序分析图表更加高级好看，代码实操分

2022-07-02 本文已影响0人春秋至

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，即可获得示例文件与代码
说这么多了，快去实操一下吧

SCpubr包：让你的单细胞测序分析图表更加高级好看，代码实操分

seura对象的准备质控与降维

三种算法简单可视化一下：

增加一列celltype到sample@metada里：

使用SCpubr的do_DimPlot函数降维：

图片调整：

换个展示形式（细胞反倒没那么好看了）：

换个颜色吧（配色鬼才）：

看看多个特征（如基因，线粒体含量等）的分布：

**弱化某些细胞亚群的展示，即标成灰色，只看我们关注的细胞群之间的差异

小提琴图

换个配色（配色鬼才）

蜂群图

同样可以使用viridis_color_map参数更换魔鬼配色

**画个点图吧家人们

当然也可以分区块进行细胞注释

换个配色

柱状图展示一下各细胞群的cluster含量

单细胞的基因集富集展示：

猜你喜欢

热点阅读