8.6 Harmony, 3’ vs 5’ 10k PBMC
使用harmony比任何其他方法都要快得多,并且在最近的標準測試中發(fā)現(xiàn)其表現(xiàn)相當好,還可以方便地與Seurat交互。讓我們首先合并對象(不進行整合)。
> pbmc_harmony <- merge(srat_3p,srat_5p)
Warning message:
In CheckDuplicateCellNames(object.list = objects) :
Some cell names are duplicated across objects provided. Renaming to enforce unique cell names.
進行和以前一樣的分析:
> pbmc_harmony <- NormalizeData(pbmc_harmony, verbose = F)
> pbmc_harmony <- FindVariableFeatures(pbmc_harmony, selection.method = "vst", nfeatures = 2000, verbose = F)
> pbmc_harmony <- ScaleData(pbmc_harmony, verbose = F)
> pbmc_harmony <- RunPCA(pbmc_harmony, npcs = 30, verbose = F)
> pbmc_harmony <- RunUMAP(pbmc_harmony, reduction = "pca", dims = 1:30, verbose = F)
再次繪制整合前的UMAP圖:
> DimPlot(pbmc_harmony,reduction = "umap") +
plot_annotation(title = "10k 3' PBMC and 10k 5' PBMC cells, before integration")

在合并的Seurat對象上運行RunHarmony,使用orig.ident作為批處理:
> pbmc_harmony <- pbmc_harmony %>% RunHarmony("orig.ident", plot_convergence = T)
Harmony 1/10
Harmony 2/10
Harmony 3/10
Harmony 4/10
Harmony 5/10
Harmony converged after 5 iterations
Warning: Invalid name supplied, making object name syntactically valid. New object name is Seurat..ProjectDim.RNA.harmony; see ?make.names for more details on syntax validity

檢查生成的降維結(jié)果:
> harmony_embeddings <- Embeddings(pbmc_harmony, 'harmony')
> harmony_embeddings[1:5, 1:5]
harmony_1 harmony_2 harmony_3 harmony_4 harmony_5
AAACCCACATAACTCG-1_1 -9.206607 -2.351619 -2.374652 -1.897186 -1.011885
AAACCCACATGTAACC-1_1 7.124223 21.600131 -0.292039 1.530283 -5.792142
AAACCCAGTGAGTCAG-1_1 -18.134725 3.405369 5.256459 4.220001 3.961466
AAACGAACAGTCAGTT-1_1 -18.103262 15.279955 12.301681 -18.115094 31.785955
AAACGAACATTCGGGC-1_1 11.097966 -2.330278 -2.723953 1.546468 1.552332
查看PCA圖:
> p1 <- DimPlot(object = pbmc_harmony, reduction = "harmony", pt.size = .1, group.by = "orig.ident") + NoLegend()
> p2 <- VlnPlot(object = pbmc_harmony, features = "harmony_1", group.by = "orig.ident", pt.size = .1) + NoLegend()
> plot_grid(p1,p2)

進行UMAP和聚類:
> pbmc_harmony <- pbmc_harmony %>%
RunUMAP(reduction = "harmony", dims = 1:30, verbose = F) %>%
FindNeighbors(reduction = "harmony", k.param = 10, dims = 1:30) %>%
FindClusters() %>%
identity()
繪制與上章相似的UMAP圖:
> pbmc_harmony <- SetIdent(pbmc_harmony,value = "orig.ident")
> DimPlot(pbmc_harmony,reduction = "umap") +
plot_annotation(title = "10k 3' PBMC and 10k 5' PBMC cells, after integration (Harmony)")

> DimPlot(pbmc_harmony, reduction = "umap", group.by = "orig.ident", pt.size = .1, split.by = 'orig.ident') + NoLegend()

harmony的結(jié)果看起來比Seurat的要差一點:
> pbmc_harmony <- SetIdent(pbmc_harmony,value = "seurat_clusters")
> DimPlot(pbmc_harmony,label = T) + NoLegend()

最后來看一下cluster內(nèi)容:
> plot_integrated_clusters(pbmc_harmony)

> rm(pbmc_harmony)
Cluster及其內(nèi)容與我們在Seurat整合后獲得的結(jié)果非常相似。為了進行更詳細的分析,需要進行細胞類型注釋。
8.7 LIGER, 3’ vs 5’ 10k PBMC
與其他方法類似,我們創(chuàng)建一個統(tǒng)一的對象并對其進行標準化/HVG/scale。LIGER在scale時不會將數(shù)據(jù)中心化,因此ScaleData中指定了do.center選項。最后兩個函數(shù)是使用orig.ident作為批處理變量運行rliger的包裝器。
> pbmc_liger <- merge(srat_3p,srat_5p)
> pbmc_liger <- NormalizeData(pbmc_liger)
> pbmc_liger <- FindVariableFeatures(pbmc_liger)
> pbmc_liger <- ScaleData(pbmc_liger, split.by = "orig.ident", do.center = F)
> pbmc_liger <- RunOptimizeALS(pbmc_liger, k = 30, lambda = 5, split.by = "orig.ident")
> pbmc_liger <- RunQuantileNorm(pbmc_liger, split.by = "orig.ident")
可以選擇在RunQuantileNorm之后執(zhí)行Louvain聚類(FindNeighbors和FindClusters),以將結(jié)果與之前的整合方法進行比較。這里使用相同的參數(shù)(k=10,Louvain聚類的默認分辨率)。
> pbmc_liger <- FindNeighbors(pbmc_liger,reduction = "iNMF",k.param = 10,dims = 1:30)
> pbmc_liger <- FindClusters(pbmc_liger)
和前面一樣,我們進行降維并繪圖:
> pbmc_liger <- RunUMAP(pbmc_liger, dims = 1:ncol(pbmc_liger[["iNMF"]]), reduction = "iNMF", verbose = F)
> pbmc_liger <- SetIdent(pbmc_liger,value = "orig.ident")
> DimPlot(pbmc_liger,reduction = "umap") + plot_annotation(title = "10k 3' PBMC and 10k 5' PBMC cells, after integration (LIGER)")

> DimPlot(pbmc_liger, reduction = "umap", group.by = "orig.ident", pt.size = .1, split.by = 'orig.ident') + NoLegend()

使用LIGER整合數(shù)據(jù)進行聚類似乎更加精細:

這些cluster看起來非常不同,并且每個cluster的分布似乎證實了這一點(兩個cluster分別被認為是3'和5'數(shù)據(jù)集所特有的):
> plot_integrated_clusters(pbmc_liger)

> rm(pbmc_liger)
> rm(srat_3p)
> rm(srat_5p)
往期內(nèi)容:
重生之我在劍橋大學(xué)學(xué)習(xí)單細胞RNA-seq分析——8. scRNA-seq數(shù)據(jù)整合(1)
重生之我在劍橋大學(xué)學(xué)習(xí)單細胞RNA-seq分析——8. scRNA-seq數(shù)據(jù)整合(2)