冠狀動脈疾病(Coronary artery disease,CAD),是全球最常見的造成死亡的心血管疾病之一。隨著高通量測序技術的發展,許多研究借此對CAD的分子機制展開研究。但是以往大量的研究都聚焦于CAD患者與正常患者間的差異而忽視了CAD患者間的差異。事實是CAD患者存在臨床上的異質性。作者希望像研究腫瘤一樣,根據CAD患者的基因表達譜進行分型,研究引起CAD的分子機制。
開始之前先看一下數據來源,作者一共分析了GEO數據庫中兩個數據集的樣本的外周血芯片數據(nCAD=352,nControl=263)。兩個數據集分別為GSE12288(n=222),GSE20686(n=393)。其中GSE20686包含兩個批次GSE20680,GSE20681)。兩個數據集的兩個芯片平臺共同檢測到的基因為11314個
作者使用SVA包中的ComBat函數消除三組不同批次樣本間的批次效應,并用消除批次前后的樣本主成分分析圖來驗證。
A:消除批次效應前的主成分分析圖,三個批次的樣本分被區分開(不同顏色的點代表不同批次的樣本)
B:消除批次效應后的主成分分析圖,三個批次的樣本混合在一起,說明批次效應被消除
圖1. 對各樣本的主成分分析
在消除了各樣本基因表達譜間的批次效應后,作者想要根據基因表達譜對CAD樣本進行聚類,劃分CAD亞型。這里作者采用的是一致性聚類的方法(ConsensusClusterPlus包)。B圖是不同聚類數目下各亞組的聚類一致性得分。當聚類數目為3時,各個亞組的聚類一致性得分均大于0.8,故作者決定將樣本聚成三類。A圖時聚類數目為3時的一致性矩陣熱圖。
作者通過以上方法將樣本分為了3個亞組(n=68,183,97)
圖2. 對CAD樣本采取一致性聚類
根據一致性聚類將CAD患者分為三個亞組后,作者對三組間在臨床特征上的差異展開研究
A:三組在男性比例上的差異,只在Ⅱ和Ⅲ組間存在顯著差異
B:三組在CAD指數(一種評價CAD患者等級的指標,這里只有GSE12288數據集中樣本的數據),Ⅰ組與其它兩組相比,CAD指數顯著偏高
C:三組在年齡上的差異,Ⅰ組的年齡與其它兩組相比顯著偏大
圖3. 比較不同亞組之間在臨床特征上的差異
此外,作者用雙因素方差分析比較自己得到的CAD分類,年齡(兩個主效應)以及兩者的交互效應對CAD患者CAD指數的影響(表1),發現自己得到的CAD分類是一個獨立于年齡的顯著效應(p=0.0369),說明CAD患者根據基因表達譜進行分類不僅影響CAD患者的CAD指數,也暗示了CAD患者內在的生物學差異
表1. 方差分析的結果
作者在三個CAD亞組間,每個CAD亞組與對照組之間進行基因差異表達分析(表2)。與以往比較某個基因在兩組間的logFC不同,這里認定一個基因是差異表達基因的條件是該基因在在兩個比較組間均值(mean)的差異要大于0.2(可能由于在CAD中基因差異表達倍數并不大所以選用此標準)且adjp<0.05
在三個CAD亞組間的差異表達分析中,Ⅰ組有2410個表達上調的DEGs,Ⅱ組中有40個表達上調的DEGs,Ⅲ組中有2217個表達上調的DEGs(表2第4列)
表2. 不同亞組間基因差異分析結果
上述通過組間比較得到的各組中表達上調的基因,需要排除在對照組中也高表達的可能
A-C:作者根據上文差異分析得到的三個基因集在各CAD亞組與對照組間進行GSEA分析,結果表明這些基因顯著富集在各自的亞組中
D:WGCNA分析得到的6個共表達模塊中的基因(縱軸)在各亞組以及對照組樣本(橫軸)中的表達量熱圖
圖4. CAD亞組中特定的表達上調基因的表達模式
作者根據上文的得到的在三個CAD亞組間差異表達的共4667個基因對所有CAD樣本進行WCGNA分析(基因加權共表達網絡分析)。共得到6個共表達模塊(圖4.D),表一中的第5列標出了各亞組間表達上調的基因對應的共表達模塊
圖4.D 共表達模塊中基因在各樣本中表達量熱圖
在得到了6個共表達模塊后,作者從MsigDB數據庫中獲取了全部的KEGG通路基因集,對共表達模塊中的基因近行KEGG通路富集分析,選出了各模塊富集結果最顯著的通路進行展示(圖5,顏色越深表示富集越顯著)
最后,作者分析了6個共表達模塊的模塊特征基因與CAD患者年齡以及CAD指數的相關性
模塊1,3,6與CAD患者的年齡以及CAD指數呈顯著正相關
模塊2,4與CAD患者的年齡以及CAD指數呈顯著負相關
圖. 6模塊特征基因與臨床特征的相關性熱圖
小結
本文作為一篇非腫瘤生信分析例子有很多地方值得我們學習。作者先對跨平臺的數據集消除批次效應,以根據基因表達譜對CAD患者分型為目的,用一致性聚類的方法將CAD患者分為三個亞組,并比較了三組在臨床特征上的差異。接著作者分析了三個亞組間各自的表達上調基因,據此對所有樣本進行WCGNA分析,得到了六個共表達模塊。最后作者對各模塊基因集進行KEGG通路分析以及研究各模塊特征基因與臨床特征的相關性。