2021-08-06
上一篇10x單細胞分析我們介紹了如何用Seurat對單細胞數據進行分群分析，這一篇我們介紹一下多個單細胞樣品的分析方法。
測試數(shu)據(ju)選擇(ze)Seurat提(ti)供(gong)的ifnb數(shu)據(ju)集，其中(zhong)包含(han)兩個(ge)PBMC樣(yang)品，一(yi)個(ge)為干(gan)擾素刺激處(chu)理，另一(yi)個(ge)為對照。打開(kai)Rstudio，在(zai)控制臺中(zhong)輸入
devtools::install_github('satijalab/seurat-data')
library(Seurat)
library(SeuratData)
library(patchwork)
InstallData("ifnb")
LoadData("ifnb")下(xia)載好的ifnb數據是一個Seurat對象，我(wo)們(men)需要按照處(chu)理方(fang)式把它(ta)拆成兩個樣品的對象
ifnb.list <- SplitObject(ifnb, split.by = "stim")
輸入ifnb.list查看拆分(fen)后(hou)的對象(xiang)
Seurat在進行多個單(dan)細胞(bao)樣品分(fen)析時，有兩(liang)種處理方式，直接合并和(he)批(pi)次效應校正后(hou)合并。我們分(fen)別介紹這(zhe)兩(liang)種分(fen)析方法(fa)。
1、直接合并
直接(jie)(jie)使用merge函數(shu)將(jiang)兩個(ge)樣(yang)品合并(bing)，這個(ge)方(fang)法直接(jie)(jie)將(jiang)兩個(ge)對(dui)(dui)象內的數(shu)據(ju)矩陣(zhen)分別合并(bing)到一起，生成一個(ge)新(xin)的Seurat對(dui)(dui)象
merged.ifnb<- merge(x = ifnb.list[[1]],y = ifnb.list[2])
合并后的數據(ju)可以當(dang)成(cheng)一(yi)個(ge)單(dan)樣品的Seurat對象進行分(fen)析
merged.ifnb <- NormalizeData(merged.ifnb, normalization.method = "LogNormalize", scale.factor = 10000) merged.ifnb <- FindVariableFeatures(merged.ifnb, selection.method = "vst", nfeatures = 2000) merged.ifnb <- ScaleData(merged.ifnb)merged.ifnb <- RunPCA(merged.ifnb) merged.ifnb <- FindNeighbors(merged.ifnb, dims = 1:30) merged.ifnb <- FindClusters(merged.ifnb) merged.ifnb <- RunUMAP(merged.ifnb, dims = 1:30)
用DimPlot查看(kan)合(he)并對象的分群結果
p1 <- DimPlot(merged.ifnb, reduction = "umap", group.by = "stim") p2 <- DimPlot(merged.ifnb, reduction = "umap", label = TRUE, repel = TRUE)p1 + p2
左(zuo)邊的(de)為兩個(ge)(ge)樣(yang)品(pin)的(de)細胞在UMAP結(jie)果中(zhong)的(de)分布，右邊的(de)為細胞cluster在UMAP中(zhong)的(de)分布。可以看到兩個(ge)(ge)樣(yang)品(pin)幾(ji)乎沒有重疊的(de)細胞群。
2、批次效應矯正
Seurat提供(gong)了一組Integration方法來去除批次效(xiao)應，這些方法首先識別處(chu)于匹配生物學狀態（anchors）的跨數(shu)據集(ji)細胞對(dui)，然后基于這些anchors校(xiao)正數(shu)據集(ji)之間的批次效(xiao)應
首先提取用(yong)來進(jin)行Integration的基因(yin)，然后找到(dao)anchors，基于(yu)anchors進(jin)行批次(ci)效應矯(jiao)正
ifnb.list <- lapply(X = ifnb.list, FUN = function(x) {
x <- NormalizeData(x)
x <- FindVariableFeatures(x, selection.method = "vst", nfeatures = 2000)
})
features <- SelectIntegrationFeatures(object.list = ifnb.list)ifnb.anchors <- FindIntegrationAnchors(object.list = ifnb.list, anchor.features = features)ifnb.combined <- IntegrateData(anchorset = ifnb.anchors)ifnb.combined是(shi)兩個樣(yang)品經過批次(ci)效應矯正后(hou)合并的Seurat對(dui)象，對(dui)這個對(dui)象進行分(fen)群(qun)分(fen)析
DefaultAssay(ifnb.combined) <- "integrated" ifnb.combined <- ScaleData(ifnb.combined, verbose = FALSE) ifnb.combined <- RunPCA(ifnb.combined, npcs = 30, verbose = FALSE) ifnb.combined <- RunUMAP(ifnb.combined, reduction = "pca", dims = 1:30) ifnb.combined <- FindNeighbors(ifnb.combined, reduction = "pca", dims = 1:30) ifnb.combined <- FindClusters(ifnb.combined)
用DimPlot查看合并(bing)對象的分群結果
可(ke)以看到(dao)經過(guo)批(pi)次矯正后，兩個樣品的大部分細胞在UMAP上都是重疊的。
使用Seurat提供的細胞類型(xing)信息對(dui)兩種分析的分群(qun)結果進行注(zhu)釋。
測試樣(yang)品(pin)經過批次(ci)(ci)矯正(zheng)后，相同(tong)類(lei)型(xing)的(de)細(xi)胞都聚到了一起，這樣(yang)更有利于對細(xi)胞類(lei)型(xing)進行(xing)鑒定。實際分析(xi)時(shi)，我們可以對自(zi)己的(de)樣(yang)品(pin)同(tong)時(shi)做這兩種(zhong)分析(xi)，根(gen)據結(jie)果來判(pan)斷是否需要進行(xing)批次(ci)(ci)矯正(zheng)。
3、差異分析
可(ke)以對相同細(xi)胞(bao)類型在(zai)不同的(de)樣品間的(de)表達差(cha)異分析，這里以B細(xi)胞(bao)為例。
首先，從(cong)對象中提取出B細胞
DefaultAssay(ifnb.combined) <- "RNA" B.ifnb <- subset(ifnb.combined, idents = "B")
將B細胞按(an)照樣品類型分開
Idents(B.ifnb) <- "stim"
使用FindMarkers比較兩個樣品間的差(cha)異
b.interferon.response <- FindMarkers(B.ifnb, ident.1 = "STIM", ident.2 = "CTRL", verbose = FALSE)
可(ke)以用(yong)VlnPlot可(ke)視化一下差異(yi)基因在不同(tong)細胞類型(xing)中的(de)表達
VlnPlot(ifnb.combined, features = c("LYZ", "ISG15", "CXCL10"), split.by = "stim", group.by = "celltype", pt.size = 0,ncol=1)
獲得的差異基因還可以進一步的挖掘，如GO/KEGG富集分析等。這里就不展開了，后續的單細胞分析再一一分析吧！
更多資訊(xun)，請關注派森諾生物官網
