哈佛大學單細胞課程|筆記匯總 （七）

哈佛大學劉小樂教授講授的計算生物學和生物信息學導論 (2020 視頻+資料）

（八）Single-cell RNA-seq clustering analysis

得到高質量的細胞后，我們可以探究細胞群中的不同細胞類型。

Exploration of quality control metrics

為了確定我們的細胞亞群是否可能是由于諸如細胞周期階段或線粒體表達之類的假象引起的，可以對不同cluster之間進行多種可視化分析判斷影響因素。這里我們將調用DimPlot()和FeaturePlot()函數來分析并可視化多種質控指標。

按照樣本對細胞亞群進行分離

從探索每個樣品中不同細胞亞群開始分析：

# 從Seurat對象中提取身份和樣本信息，以確定每個樣本的每個簇的細胞數 n_cells <- FetchData(seurat_integrated, vars = c("ident", "orig.ident")) %>%dplyr::count(ident, orig.ident) %>%tidyr::spread(ident, n)# View table View(n_cells)

使用UMAP對每個樣品不同細胞亞群進行可視化：

# UMAP of cells in each cluster by sample DimPlot(seurat_integrated, label = TRUE, split.by = "sample") + NoLegend()

從上圖沒有看到單獨由樣品來源引起的聚類簇。

按照細胞周期對細胞亞群進行分離

我們可以探討細胞是否在不同的細胞周期階段聚類。當我們利用SCTransform執行歸一化和無意義的變異源回歸（regression of uninteresting sources of variation）時，沒有考慮細胞周期的影響。

如果我們的細胞簇在不同細胞周期上顯示出很大差異，則表明我們要重新運行SCTransform并將S.Score和G2M.Score (細胞周期得分)添加到變量中以進行回歸，然后重新運行其余步驟 。

# Explore whether clusters segregate by cell cycle phase DimPlot(seurat_integrated,label = TRUE, split.by = "Phase") + NoLegend()

從上圖未能看到單個細胞周期主導的聚類簇。

可視化其它指標對細胞亞群分群的影響

我們還可以探索其他指標，例如每個細胞的UMI和基因數量，S期和G2M期的標記，以及通過UMAP進行的線粒體基因表達。其中各個S和G2M分數可提供更多的周期信息。

# Determine metrics to plot present in seurat_integrated@meta.data metrics <- c("nUMI", "nGene", "S.Score", "G2M.Score", "mitoRatio")FeaturePlot(seurat_integrated, reduction = "umap", features = metrics,pt.size = 0.4, sort.cell = TRUE,min.cutoff = 'q10',label = TRUE)

NOTE: The sort.cell argument will plot the positive cells above the negative cells, while the min.cutoff argument will determine the threshold for shading. A min.cutoff of q10 translates to the 10% of cells with the lowest expression of the gene will not exhibit any purple shading (completely gray).

這些指標在不同的cluster中還是比較均勻的，但nUMI和nGene在cluster 3、9、14和15以及也許還有cluster 17中顯示出了更高的值。

探索導致不同分群clusters的PC

我們可以探索不同PC分群clusters的情況，以確定定義的PC是否能夠很好地區分細胞類型。為了可視化此信息，我們需要提取細胞的UMAP坐標信息以及每個PC的相應分數，并通過UMAP查看。

首先識別Seurat對象中的待提取信息，然后使用FetchData（）函數提取。

# Defining the information in the seurat object of interest columns <- c(paste0("PC_", 1:16),"ident","UMAP_1", "UMAP_2")# Extracting this data from the seurat object pc_data <- FetchData(seurat_integrated, vars = columns)

NOTE: How did we know in the FetchData() function to include UMAP_1 to obtain the UMAP coordinates? The Seurat cheatsheet（https://satijalab.org/seurat/essential_commands.html） describes the function as being able to pull any data from the expression matrices, cell embeddings, or metadata.

For instance, if you explore the seurat_integrated@reductions list object, the first component is for PCA, and includes a slot for cell.embeddings. We can use the column names (PC_1, PC_2, PC_3, etc.) to pull out the coordinates or PC scores corresponding to each cell for each of the PCs.

We could do the same thing for UMAP:

# Extract the UMAP coordinates for the first 10 cells seurat_integrated@reductions$umap@cell.embeddings[1:10, 1:2] The FetchData() function just allows us to extract the data more easily.

利用FetchData（）函數能更輕松地提取數據。

在以下UMAP圖中，細胞將根據PC score進行上色。下面的是top16 PCs：

# Adding cluster label to center of cluster on UMAP umap_label <- FetchData(seurat_integrated, vars = c("ident", "UMAP_1", "UMAP_2")) %>%group_by(ident) %>%summarise(x=mean(UMAP_1), y=mean(UMAP_2))# Plotting a UMAP plot for each of the PCs map(paste0("PC_", 1:16), function(pc){ggplot(pc_data, aes(UMAP_1, UMAP_2)) +geom_point(aes_string(color=pc), alpha = 0.7) +scale_color_gradient(guide = FALSE, low = "grey90", high = "blue") +geom_text(data=umap_label, aes(label=ident, x, y)) +ggtitle(pc) }) %>% plot_grid(plotlist = .)

我們可以看到clusters是如何由不同的PC表示。例如，驅動PC_2的基因在cluster 6、11和17中貢獻更大（在15中也可能更高），因此我們可以回顧一下驅動此PC的基因，以了解細胞類型可能是什么：

# Examine PCA results print(seurat_integrated[["pca"]], dims = 1:5, nfeatures = 5)

使用CD79A基因和HLA基因作為PC_2的陽性標記，我們可以假設cluster 6、11和17對應于B細胞。這一點只能為我們細胞分群提供線索，cluster的真正身份還需要通過PC組合來確定。

探究已知的細胞類型marker

細胞聚類之后，我們可以通過尋找已知的marker來探索細胞類型。下圖顯示了帶有marker的簇的UMAP圖，然后顯示了預期的不同細胞類型。

DimPlot(object = seurat_integrated, reduction = "umap", label = TRUE) + NoLegend()

可以利用FeaturePlot()對細胞marker進行展示：

利用Seurat的FeaturePlot（）函數可以輕松地在UMAP可視化上探索已知的標記，進而確定簇的身份。同時我們可以使用RNA assay slot中的標準化計數數據，來訪問所有基因的表達水平（而不僅僅是3000個高度可變的基因）。

# Select the RNA counts slot to be the default assay DefaultAssay(seurat_integrated) <- "RNA"# Normalize RNA data for visualization purposes seurat_integrated <- NormalizeData(seurat_integrated, verbose = FALSE)

同時還需要考慮簇與簇之間marker表達的一致性。例如，如果一個細胞類型有兩個標記基因，而簇中僅表達了其中一個，那么我們就不能絕對地將該簇分配給該細胞類型。

CD14+ monocyte markers

FeaturePlot(seurat_integrated, reduction = "umap", features = c("CD14", "LYZ"), sort.cell = TRUE,min.cutoff = 'q10', label = TRUE)

從Marker表達來看，CD14+單核細胞似乎與cluster 1、3和14相對應。

FCGR3A+ monocyte markers

FeaturePlot(seurat_integrated, reduction = "umap", features = c("FCGR3A", "MS4A7"), sort.cell = TRUE,min.cutoff = 'q10', label = TRUE)

從Marker表達來看，FCGR3A+單核細胞與cluster 9對應。

Macrophages

FeaturePlot(seurat_integrated, reduction = "umap", features = c("MARCO", "ITGAM", "ADGRE1"), sort.cell = TRUE,min.cutoff = 'q10', label = TRUE)

似乎沒有cluster對應于巨噬細胞。

Conventional dendritic cell markers

FeaturePlot(seurat_integrated, reduction = "umap", features = c("FCER1A", "CST3"), sort.cell = TRUE,min.cutoff = 'q10', label = TRUE)

conventional dendritic cells（cDC）對應于cluster 15。

Plasmacytoid dendritic cell markers

FeaturePlot(seurat_integrated, reduction = "umap", features = c("IL3RA", "GZMB", "SERPINF1", "ITM2C"), sort.cell = TRUE,min.cutoff = 'q10', label = TRUE)

Plasmacytoid dendritic cells （pDC）對應于cluster 19。

NOTE: If any cluster appears to contain two separate cell types, it’s helpful to increase our clustering resolution to properly subset the clusters. Alternatively, if we still can’t separate out the clusters using increased resolution, then it’s possible that we had used too few principal components such that we are just not separating out these cell types of interest. To inform our choice of PCs, we could look at our PC gene expression overlapping the UMAP plots and determine whether our cell populations are separating by the PCs included.

上面的分析中還存著待解決的問題：

cluster 7和20的細胞類型標識是什么？

對應于相同細胞類型的簇是否具有生物學上有意義的差異？

這些細胞類型有亞群嗎？

通過為這些簇鑒定其他標記基因，我們能否更相信這些細胞類型身份？

下一步將執行marker識別分析，該分析將輸出簇之間表達差異顯著的基因，能對上面的問題進行解惑。使用這些基因，我們可以確定或提高對簇/亞群身份的信心。

總結

以上是生活随笔為你收集整理的哈佛大学单细胞课程|笔记汇总 （八）的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。