學習目標

1. 瞭解計數資料變換方法的重要性
2. 瞭解 PCA (principal component analysis)
3. 瞭解如何使用 PCA 和層次聚類評估樣本質量

1. 質控

DESeq2 工作流程的下一步是 QC，其中包括樣本和基因程度上，以對計數資料執行 QC 檢查，以幫助我們確保樣本或重複看起來良好。

2. 樣本QC

RNA-seq 分析中一個有用的初始步驟通常是評估樣本之間的整體相似性：

1. 哪些樣本彼此相似，哪些不同？
2. 這是否符合實驗設計的預期？
3. 資料集中的主要變異來源是什麼？

為了探索樣本的相似性，我們將使用主成分分析 (PCA) 和層次聚類方法執行樣本級 QC。這些方法或工具使我們能夠檢查重複彼此之間的相似程度（聚類），並確保實驗條件是資料變化的主要來源。樣品級 QC

這些無監督聚類方法使用 log2 變換的歸一化計數執行。log2 轉換改進了視覺化的距離。我們將不使用普通的 log2 變換，而是使用正則化對數變換 (rlog)，以避免因大量低計數基因而產生的任何偏差；

• 為什麼需要進行資料轉換？

許多用於多維資料探索性分析的常用統計方法，尤其是聚類和排序方法（例如，主成分分析等），最適合（至少近似地）同方差資料；這意味著可觀察量的方差（即，這裡是基因的表達值）不依賴於均值。然而，在 RNA-seq 資料中，方差隨平均值增加。例如，如果直接對歸一化讀取計數矩陣執行 PCA

• 使用 rlog的優勢：

因此，DESeq2 提供了正則化對數轉換，或簡稱為 rlog。對於計數高的基因，rlog 轉換與普通的 log2 轉換差別不大。然而，對於計數較低的基因，這些值會縮小到所有樣本中基因的平均值。這樣做是為了使 rlog 轉換後的資料近似同方差。

DESeq2 建議大型資料集（100 個樣本）使用方差穩定變換 (vst

) 而不是 rlog 來進行計數變換，因為 rlog 函式可能需要執行很長時間，而 vst() 函式在類似情況下更快。

3. PCA

主成分分析 (PCA) 是一種用於強調變化並在資料集中降維的技術。這是一種非常重要的技術，用於質量控制和 Bulk RNA-seq 和單細胞 RNA-seq 資料的分析。

3.1. PCA plots

本質上，如果兩個樣本的基因表達水平相似，這些基因對給定 PC（主成分）表示的變異有顯著貢獻，則它們將在表示該 PC 的軸上靠近繪製。因此，我們期望生物重複具有相似的分數（因為我們的期望是相同的基因正在發生變化）並聚集在一起。通過視覺化一些示例 PCA 圖最容易理解這一點。

我們在下面有一個示例資料集和一些相關的 PCA 圖，以瞭解如何解釋它們。實驗的元資料如下所示。感興趣的主要條件是處理。

在 PC1 和 PC2 上進行視覺化時，我們沒有看到樣本因處理而分開，因此我們決定探索資料中存在的其他變異來源。我們希望我們已經在我們的元資料表中包含了所有可能的已知變異源，並且我們可以使用這些因素來為 PCA 圖著色。

我們從cage因子開始，但cage因子似乎無法解釋 PC1 或 PC2 上的變化。

然後，我們按 sex 因素著色，這似乎在 PC2 上分離樣本。這是需要注意的好資訊，因為我們可以在下游使用它來解釋模型中由於 sex 引起的變化並將其迴歸。

接下來我們探索 strain 因子，發現它可以解釋 PC1 上的變化。

很高興我們能夠確定 PC1 和 PC2 的變異來源。通過在我們的模型中考慮它，我們應該能夠檢測到更多因處理而差異表達的基因。

令人擔憂的是，我們看到兩個樣本沒有與正確的 strain 聚類。這表明可能存在樣本交換，應進行調查以確定這些樣本是否確實是標記的 strain。如果我們發現有一個交換，我們可以交換元資料中的樣本。但是，如果我們認為它們被正確標記或不確定，我們可以從資料集中刪除樣本。

我們仍然沒有發現處理是否是 strain 和 sex 後變異的主要來源。因此，我們探索 PC3 和 PC4 以檢視處理是否正在推動由這兩個 PC 中的任何一個代表的變異。

我們發現樣本在 PC3 上通過處理分離，並且對我們的 DE 分析持樂觀態度，因為我們感興趣的條件，處理，在 PC3 上分離，我們可以迴歸驅動 PC1 和 PC2 的變化。

根據前幾個主要成分解釋了多少變化，您可能想要探索更多（即考慮更多成分並繪製成對組合）。即使您的樣品沒有通過實驗變數清楚地分離，您仍然可以從 DE 分析中獲得生物學相關的結果。如果您期望效果大小非常小，那麼訊號可能會被無關的變異源淹沒。在您可以識別這些來源的情況下，在您的模型中考慮這些來源很重要，因為它為檢測 DE 基因的工具提供了更多功能。

4. 層次聚類

與 PCA 類似，層次聚類是另一種互補的方法，用於識別資料集中的模式和潛在異常值。熱圖顯示資料集中所有成對樣本組合的基因表達相關性。由於大多數基因沒有差異表達，樣本之間通常具有很高的相關性（值高於 0.80）。低於 0.80 的樣本可能表示您的資料和/或樣本汙染中存在異常值。

沿軸的分層樹指示哪些樣本彼此更相似，即聚集在一起。頂部的色塊表示資料中的子結構，您會希望看到您的重複一起作為每個樣本組的一個塊。我們的期望是樣本聚集在一起類似於我們在 PCA 圖中觀察到的分組。

在下圖中， Wt_3 和 KO_3 樣本沒有與其他重複聚類在一起。我們想要探索 PCA 以檢視我們是否看到相同的樣本聚類。

5. Mov10 QC

現在我們已經很好地理解了通常用於 RNA-seq 的 QC 步驟，讓我們為 Mov10 資料集進行 QC。

5.1. 資料轉換

• 轉換 MOV10 資料集的歸一化計數

為了促進 PCA 和層次聚類視覺化方法的距離或聚類，我們需要通過對歸一化計數應用 rlog 變換來調節均值的方差。

歸一化計數的 rlog 轉換僅在該質量評估期間對於這些視覺化方法是必需的。我們不會使用這些轉換後的計數來確定差異表達。

# rlog 轉換
rld <- rlog(dds, blind=TRUE)

blind=TRUE 引數是為了確保 rlog() 函式不考慮我們的樣本組——即以無偏見的方式進行轉換。在執行質量評估時，包含此選項很重要。

rlog() 函式返回一個 DESeqTransform 物件，這是另一種特定於 DESeq 的物件。您不只是獲得轉換值矩陣的原因是因為用於計算 rlog 轉換的所有引數（即大小因子）都儲存在該物件中。我們使用此物件繪製 PCA 和層次聚類圖以進行質量評估。

5.2. PCA

我們現在已準備好進行 QC 步驟，讓我們從 PCA 開始吧！

DESeq2 有一個內建函式，可以在後臺使用 ggplot2生成 PCA 圖。這很棒，因為它使我們不必輸入程式碼行，也不必擺弄不同的 ggplot2 層。此外，它直接將 rlog 物件作為輸入，從而省去了我們從中提取相關資訊的麻煩。

plotPCA() 需要兩個引數作為輸入：DESeqTransform 物件和 intgroup，即元資料中包含有關實驗樣本組資訊列的名稱。

# Plot PCA 
plotPCA(rld, intgroup="sampletype")

預設情況下，plotPCA()使用前 500 個最易變的基因。您可以通過新增 ntop= 引數並指定您希望函式考慮的基因數量來更改此設定。

plotPCA() 函式將只返回 PC1 和 PC2 的值。如果您想探索資料中的其他 PC，或者如果您想識別對 PC 貢獻最大的基因，您可以使用 prcomp() 函式。例如，要繪製任何 PC，我們可以執行以下程式碼：

 # Input is a matrix of log transformed values
 rld <- rlog(dds, blind=T)
 rld_mat <- assay(rld)
 pca <- prcomp(t(rld_mat))

 # Create data frame with metadata and PC3 and PC4 values for input to ggplot
 df <- cbind(meta, pca$x)
 ggplot(df) + geom_point(aes(x=PC3, y=PC4, color = sampletype))

5.3. Hierarchical Clustering

• MOV10 資料集層次聚類

DESeq2中沒有內建函式來繪製熱圖來顯示所有樣本之間的成對相關性和層次聚類資訊；我們將使用 pheatmap 包中的 pheatmap() 函式。此函式不能使用 DESeqTransform 物件作為輸入，但需要矩陣或資料框。因此，要做的第一件事是使用名為 assay() 的函式，從 rld 物件檢索該資訊，該函式將 DESeqTransform 物件中的資料轉換為簡單的二維資料結構。

# Extract the rlog matrix from the object
rld_mat <- assay(rld)    

接下來，我們需要計算所有樣本的成對相關值。我們可以使用 cor() 函式來做到這一點：

# Compute pairwise correlation values
rld_cor <- cor(rld_mat) 

讓我們看一下相關矩陣的列名和行名。

head(rld_cor)  

head(meta)  

您會注意到它們與我們在開始時使用的元資料資料框中為樣本提供的名稱相匹配。這很重要，因此我們可以使用下面的註釋引數在頂部繪製一個色塊。此塊可輕鬆實現層次聚類的視覺化。

# Load pheatmap package
library(pheatmap)

# Plot heatmap using the correlation matrix and the metadata object
pheatmap(rld_cor, annotation = meta)

當您使用 pheatmap() 進行繪圖時，層次聚類資訊用於將相似的樣本放在一起，並且該資訊由沿軸的樹結構表示。註釋引數接受一個數據框作為輸入，在我們的例子中它是元資料框。

總體而言，我們觀察到高相關性 (> 0.999)，表明沒有異常樣本。此外，與 PCA 圖類似，您會看到樣本按樣本組聚集在一起。總之，這些圖向我們表明資料質量很好，我們有信心可以進行差異表達分析。

歡迎Star -> 學習目錄

更多教程 -> 學習目錄

本文由mdnice多平臺釋出