樣品的無監督聚類
檢查基因表達分析最重要的探索性政策之一是多元定标(MDS)圖或類似的圖。該圖以無監督的方式顯示了樣本之間的相似性和不相似性,以便人們可以了解在進行正式測試之前可以檢測差異表達的程度。理想情況下,樣本會在感興趣的主要條件下很好地聚類,并且可以識别任何偏離其群組的樣本,并追蹤誤差來源或額外變化。
lcpm <- cpm(x, log=TRUE)
par(mfrow=c(1,2))
col.group <- group
levels(col.group) <- brewer.pal(nlevels(col.group), "Set1")
col.group <- as.character(col.group)
col.lane <- lane
levels(col.lane) <- brewer.pal(nlevels(col.lane), "Set2")
col.lane <- as.character(col.lane)
plotMDS(lcpm, labels=group, col=col.group)
title(main="A. Sample groups")
plotMDS(lcpm, labels=lane, col=col.lane, dim=c(3,4))
title(main="B. Sequencing lanes")

這種方法可以使用plotMDS函數在limma包中進行。第一個次元解釋了資料中最大的變化比例,其後的次元具有較小的影響并且與之前的次元正交。當實驗設計涉及多個因素時,建議在多個次元上檢查每個因子。如果樣本以任何這些次元中的給定因子聚類,則表明該因子有助于表達差異,并且值得包括線上性模組化中。另一方面,影響很小或沒有影響的因素可能會被排除在下遊分析之外。在這個資料集中,可以看到樣本在次元1和次元2内的實驗組内很好地聚類。雖然所有樣本都是按照群組聚集的,但是觀察到在基礎和LP之間以及基線和ML在次元1上的最大轉錄差異。由于這個原因,預計細胞群體之間的配對比較将導緻用于比較的更多數量的DE基因涉及基礎樣本,并且在比較ML和LP時涉及相對少量的DE基因。
- 另外,Glimma軟體包提供了可以探索多個次元的互動式MDS圖的便利。 glMDSPlot函數在左側面闆中生成一個帶有MDS圖的html頁面(如果launch = TRUE,則在浏覽器中打開該頁面),并顯示右側面闆中每個次元解釋的變化比例的barplot。 單擊條形圖的條形圖會更改MDS圖形中繪制的一對次元,然後懸停在各個點上可以顯示樣本标簽。 顔色方案也可以改變以突出細胞群或測序泳道(批次)。 這個資料集的互動式MDS圖可以在 http://bioinf.wehi.edu.au/folders/limmaWorkflow/glimma-plots/MDS-Plot.html 找到。
library(Glimma)
glMDSPlot(lcpm, labels=paste(group, lane, sep="_"), groups=x$samples[,c(2,5)],
launch=T)
差異表達分析
建立一個設計矩陣和對比
在這項研究中,我們感興趣的是看到哪些基因在三種細胞群體之間的不同水準上表達。 在我們的分析中,假設基礎資料是正态分布的,假設線性模型符合資料。 為了開始,設計矩陣與細胞群體和測序泳道(批次)資訊一起建立。
design <- model.matrix(~0+group+lane)
colnames(design) <- gsub("group", "", colnames(design))
design
contr.matrix <- makeContrasts(
BasalvsLP = Basal-LP,
BasalvsML = Basal - ML,
LPvsML = LP - ML,
levels = colnames(design))
contr.matrix
從計數資料中去除雜波
已經顯示,對于RNA-seq計數資料,方差并不依賴于平均值13 - 原始計數或轉換為對數CPM值時,情況屬實。使用負二項分布模型計數的方法假定二次均方差關系。在limma中,對log-CPM值進行線性模組化,log-CPM值假定為正态分布,并使用由voom函數計算的精确權重來适應均值 - 方差關系。在DGEList對象上運作時,voom通過自動從x本身提取庫大小和規範化因子,将原始計數轉換為log-CPM值。可以使用normalize.method參數在voom中指定對日志CPM值的其他标準化。
v <- voom(x, design, plot=TRUE)
v
vfit <- lmFit(v, design)
vfit <- contrasts.fit(vfit, contrasts=contr.matrix)
efit <- eBayes(vfit)
plotSA(efit)
左圖顯示了該資料集的log-CPM值的均值 - 方差關系。典型地,“voom-plot”顯示由測序實驗中的技術變化和來自不同細胞群的重複樣品之間的生物變異的組合導緻的均值和變化之間的下降趨勢。具有高生物學變異的實驗通常導緻更平坦的趨勢,其中方內插補點在高表達值時平穩。低生物變異的實驗傾向于導緻急劇下降的趨勢。右圖表示基因方差不受基因平均值影響。
檢查DE基因的數量
為了快速了解差異表達水準,可以在表格中總結顯着上調和下調基因的數目。 顯着性是使用預設設定為5%的調整後的p值截止值來定義的。 為了比較Basal和LP中的表達水準總共8,425個DE基因。 在Basal和ML之間發現總共8,510個DE基因(4,338個下調基因和4,172個上調基因),并且在LP和ML之間發現總共5,340個DE基因(2,895個下調和2,445個上調)。 觀察到的涉及基礎群體的較大數量的DE基因與我們對MDS圖的觀察結果一緻。
summary(decideTests(efit))
BasalvsLP BasalvsML LPvsML
Down 4127 4338 2895
NotSig 5740 5655 8825
Up 4298 4172 2445