今天小云想為大家分享一下如何對(duì)單細(xì)胞數(shù)據(jù)進(jìn)行處理和分型并繪制高顏值的細(xì)胞比例柱狀圖并做富集檢驗(yàn)，又是一種新的分析思路，有興趣的小伙伴可以跟著小果開(kāi)始今天的學(xué)習(xí)，絕對(duì)是干貨滿滿奧！

1.如何繪制高顏值的細(xì)胞比例柱狀圖并做富集檢驗(yàn)？

首先利用Seurat包進(jìn)行單細(xì)胞數(shù)據(jù)分析，對(duì)分析后的數(shù)據(jù)進(jìn)行樣本分類(lèi)信息提取，然后對(duì)各組樣本進(jìn)行 one vs other的fisher檢驗(yàn)，進(jìn)行多重性校正，得到各組的p-adj，最終繪制細(xì)胞比例堆疊柱狀圖并添加富集檢驗(yàn)結(jié)果；接下來(lái)跟著小果開(kāi)始今天的實(shí)操內(nèi)容吧！如果覺(jué)得推文不錯(cuò)，點(diǎn)贊加關(guān)注奧。

2.準(zhǔn)備需要的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.數(shù)據(jù)讀取

# 讀取count矩陣、細(xì)胞Barcode、基因名和細(xì)胞對(duì)應(yīng)的樣本信息四個(gè)文件。

#mtx格式的count矩陣文件

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

#細(xì)胞Barcode文件

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

#基因名文件

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

#細(xì)胞對(duì)應(yīng)的樣本信息文件,通過(guò)逗號(hào)分割，對(duì)應(yīng)細(xì)胞Barcode和樣本信息。

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

4.?scRNA常規(guī)分析流程

# 修改count行名

rownames(mtx) <- geneID$V1

#修改count列名

colnames(mtx) = cellID$V1

# 制作seurat對(duì)象

seu <- CreateSeuratObject(mtx)

#添加樣本信息到創(chuàng)建的seurat對(duì)象中

seu$Sample = metadata$Sample

#篩選需要的樣本數(shù)據(jù)

seu<-subset(seu,Sample%in%c("Tum_963_WT","Tum_650_KO","Tum_877_WT", "Tum_685_KO")) # 對(duì)數(shù)化表達(dá)值

seu <- NormalizeData(seu)

# 尋找高變異基因 ??????????????????????????

seu <- FindVariableFeatures(seu, nfeatures = 1000) ??

# 標(biāo)準(zhǔn)化

seu <- ScaleData(seu)

# PCA降維? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?

seu <- RunPCA(seu) ???????????????????????????????????

# 計(jì)算UMAP降維坐標(biāo)

dim.to.use = 1:12

# 進(jìn)行UMAP非監(jiān)督聚類(lèi)

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

# 構(gòu)建最近鄰圖

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

# 進(jìn)行聚類(lèi)

seu <- FindClusters(seu) ?????????????????????????????

#保存單細(xì)胞分析結(jié)果

saveRDS(seu, "seu.rds")

4.?繪制堆疊比例柱狀圖

#提取樣本分類(lèi)信息

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

#利用正則表達(dá)式將sample列進(jìn)行分割，取出第三個(gè)

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

#新增一列Run列，通過(guò)ifelse判斷

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

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

##利用正則表達(dá)式將sample列進(jìn)行分割，取出第五個(gè)

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

# 設(shè)置分組顏色

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

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

# 制作列聯(lián)表

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的順序倒過(guò)來(lái)

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

?# 設(shè)置不同組對(duì)應(yīng)的顏色

scale_fill_manual(values=cell.col)+ ?????????????????????????????????????????????????# 設(shè)置不同組對(duì)應(yīng)的顏色

??facet_wrap(~Var1,nrow=1)+ ??????????????????????????????????????????????????????????# 設(shè)置柱狀圖按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)

# 構(gòu)建seurat_clusters和Group的列聯(lián)表

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

# 進(jìn)行fisher精確檢驗(yàn)和post hoc檢驗(yàn)

# post hoc檢驗(yàn)：對(duì)各組樣本進(jìn)行 one vs other的fisher檢驗(yàn)，進(jìn)行多重性校正，得到各組的p-adj

fisher.test(tbl, simulate.p.value = T) # fisher精確檢驗(yàn)

post.hoc <- row_wise_fisher_test(tbl) # post hoc檢驗(yàn)

post.hoc$sig.label <- ifelse(test = post.hoc$p.adj.signif == "ns", # 調(diào)整顯著性顯示標(biāo)簽

?????????????????????????????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的順序倒過(guò)來(lái)

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

????scale_fill_manual(values = cell.col) + ??????????????????????????????

??# 設(shè)置y坐標(biāo)軸的刻度

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

??#標(biāo)簽的位置和具體內(nèi)容

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.?結(jié)果文件

1.?scbarplot.pdf

該結(jié)果圖片為繪制的是細(xì)胞比例堆疊柱狀圖

1.?scbarplotsig.pdf

該結(jié)果圖片為細(xì)胞比例柱狀圖并做了富集檢驗(yàn)，橫坐標(biāo)為Cluster,縱坐標(biāo)為細(xì)胞所占的比例，并添加了顯著性標(biāo)簽。

今天小云對(duì)細(xì)胞比例柱狀圖繪制的分享就到這里啦，進(jìn)行該分析需要熟練掌握單細(xì)胞分析流程，單細(xì)胞相關(guān)分析內(nèi)容可以嘗試使用本公司新開(kāi)發(fā)的云平臺(tái)生物信息分析小工具，零代碼完成分析，云平臺(tái)網(wǎng)址：http://www.biocloudservice.com/home.html，包括單細(xì)胞分析(http://www.biocloudservice.com/366/366.php),單細(xì)胞數(shù)據(jù)繪制小提琴圖(http://www.biocloudservice.com/788/788.php),繪制單細(xì)胞tSNE圖(http://www.biocloudservice.com/229/229.php)等單細(xì)胞分析相關(guān)小工具，歡迎小伙伴來(lái)嘗試奧。