Nature级别的单细胞细胞比例柱状图咋画？一文拿捏！/SCI论文/科研/研究生/生信分析热

今天小云想为大家分享一下如何对单细胞数据进行处理和分型并绘制高颜值的细胞比例柱状图并做富集检验，又是一种新的分析思路，有兴趣的小伙伴可以跟着小果开始今天的学习，绝对是干货满满奥！

1.如何绘制高颜值的细胞比例柱状图并做富集检验？

首先利用Seurat包进行单细胞数据分析，对分析后的数据进行样本分类信息提取，然后对各组样本进行 one vs other的fisher检验，进行多重性校正，得到各组的p-adj，最终绘制细胞比例堆叠柱状图并添加富集检验结果；接下来跟着小果开始今天的实操内容吧！如果觉得推文不错，点赞加关注奥。

2.准备需要的R包

#安装需要的R包

install.packages(“ggplot2”)

install.packages(“rstatix”)

install.packages(“Seurat”)

remotes::install_github("mojaveazure/seurat-disk")

install.packages("hdf5r")

install.packages(“Matrix”)

#载入需要的R包

library(ggplot2)

library(rstatix)

library(Seurat)

library(SeuratDisk)

library(Matrix)

3.数据读取

# 读取count矩阵、细胞Barcode、基因名和细胞对应的样本信息四个文件。

#mtx格式的count矩阵文件

mtx <- readMM("GSE151914_expression_matrix.mtx.gz")

#细胞Barcode文件

cellID <- read.table("GSE151914_cellIDs.txt.gz")

#基因名文件

geneID <- read.table("GSE151914_genes.txt.gz")

#细胞对应的样本信息文件,通过逗号分割，对应细胞Barcode和样本信息。

metadata <- read.table("GSE151914_metadata.txt.gz", header = T, sep = ",")

4. scRNA常规分析流程

# 修改count行名

rownames(mtx) <- geneID$V1

#修改count列名

colnames(mtx) = cellID$V1

# 制作seurat对象

seu <- CreateSeuratObject(mtx)

#添加样本信息到创建的seurat对象中

seu$Sample = metadata$Sample

#筛选需要的样本数据

seu<-subset(seu,Sample%in%c("Tum_963_WT","Tum_650_KO","Tum_877_WT", "Tum_685_KO")) # 对数化表达值

seu <- NormalizeData(seu)

# 寻找高变异基因                           

seu <- FindVariableFeatures(seu, nfeatures = 1000)   

# 标准化

seu <- ScaleData(seu)

# PCA降维                               

seu <- RunPCA(seu)                                    

# 计算UMAP降维坐标

dim.to.use = 1:12

# 进行UMAP非监督聚类

seu <- RunUMAP(seu, dims = dim.to.use)               

# 构建最近邻图

seu <- FindNeighbors(seu, dims = dim.to.use)          

# 进行聚类

seu <- FindClusters(seu)                              

#保存单细胞分析结果

saveRDS(seu, "seu.rds")

4. 绘制堆叠比例柱状图

#提取样本分类信息

cellinfo <- FetchData(seu, vars = c("Sample", "seurat_clusters"))

#利用正则表达式将sample列进行分割，取出第三个

cellinfo$Run <- gsub("(\\w+)(_)(\\d+)(_)(\\w+)", "\\3", cellinfo$Sample)

#新增一列Run列，通过ifelse判断

cellinfo$Run <- ifelse(test = cellinfo$Run %in% c("963", "650"),

                       yes = "R1", no = "R2")

##利用正则表达式将sample列进行分割，取出第五个

cellinfo$Group <- gsub("(\\w+)(_)(\\d+)(_)(\\w+)", "\\5", cellinfo$Sample)

# 设置分组颜色

cell.col <- setNames(object = c("#1A63A8", "#FC6910"),

                     nm = c("WT", "KO"))

# 制作列联表

plot.data <- as.data.frame(table(cellinfo$seurat_clusters, cellinfo$Run, cellinfo$Group))

p1 <- ggplot(plot.data, aes(x = Var2, y = Freq, fill = Var3)) +

  # 绘制堆叠比例柱状图，将WT和KO的顺序倒过来

geom_bar(stat = "identity", position = position_fill(reverse = T)) +                    

 # 设置不同组对应的颜色

scale_fill_manual(values=cell.col)+                                                  # 设置不同组对应的颜色

  facet_wrap(~Var1,nrow=1)+                                                           # 设置柱状图按seurat_cluster分别显示

  theme_classic()+                                                                       # 修改x和y轴列名

  xlab("Replicate")+                                                                    

  ylab("Fraction of Cells") +                           

 #修改主题 

theme(strip.background = element_rect(fill="grey", color = "white", size = 1),      

        strip.text=element_text(size=12,colour="black"),                                  legend.title=element_blank(),                                                          axis.title=element_text(size=15))                                             

p1

#保存图片

ggsave(file="scbarplot.pdf", width = 12, height = 4)

# 构建seurat_clusters和Group的列联表

tbl <- table(cellinfo$seurat_clusters, cellinfo$Group)

# 进行fisher精确检验和post hoc检验

# post hoc检验：对各组样本进行 one vs other的fisher检验，进行多重性校正，得到各组的p-adj

fisher.test(tbl, simulate.p.value = T) # fisher精确检验

post.hoc <- row_wise_fisher_test(tbl) # post hoc检验

post.hoc$sig.label <- ifelse(test = post.hoc$p.adj.signif == "ns", # 调整显著性显示标签

                             yes = "", no = post.hoc$p.adj.signif) # 不显示NS (No Significant)

# 绘制柱状图

plot.data <- as.data.frame(tbl)

p2 <- ggplot(plot.data, aes(x = Var1, y = Freq, fill = Var2)) +

  # 绘制堆叠比例柱状图，将WT和KO的顺序倒过来

geom_bar(stat = "identity", position = position_fill(reverse = T)) +  

    scale_fill_manual(values = cell.col) +                               

  # 设置y坐标轴的刻度

scale_y_continuous(breaks = seq(0, 1, by = 0.25)) +                   

  #标签的位置和具体内容

geom_text(data = post.hoc,                                            

            aes(x = group, y = 1.1, label = sig.label),                

            inherit.aes = F) +                                          

# 修改x和y轴列名  

xlab("Cluster") +                                                     

  ylab("Fraction of Cells") +

#修改主题  

theme_classic() +                                                     

  theme(legend.title = element_blank(),                                 

        axis.text = element_text(size = 10),                            

        axis.title = element_text(size = 15))                           

p2

#保存图片

ggsave(filename = "scbarplotsig.pdf", width = 8, height = 4)

4. 结果文件

1. scbarplot.pdf

该结果图片为绘制的是细胞比例堆叠柱状图

1. scbarplotsig.pdf

该结果图片为细胞比例柱状图并做了富集检验，横坐标为Cluster,纵坐标为细胞所占的比例，并添加了显著性标签。

今天小云对细胞比例柱状图绘制的分享就到这里啦，进行该分析需要熟练掌握单细胞分析流程，单细胞相关分析内容可以尝试使用本公司新开发的云平台生物信息分析小工具，零代码完成分析，云平台网址：http://www.biocloudservice.com/home.html，包括单细胞分析(http://www.biocloudservice.com/366/366.php),单细胞数据绘制小提琴图(http://www.biocloudservice.com/788/788.php),绘制单细胞tSNE图(http://www.biocloudservice.com/229/229.php)等单细胞分析相关小工具，欢迎小伙伴来尝试奥。