其實我一直沒用過WGCNA，因為分析網(wǎng)絡的方法有很多，但是大家好像都更愛用這個。最近幫人分析的幾組數(shù)據(jù)，他們指名要用WGCNA分析，所以就學習一下。
=======WGCNA簡介=========
WGCNA(Weighted Gene Co-Expression NetworkAnalysis, 加權基因共表達網(wǎng)絡分析)，鑒定表達模式相似的基因集合（module）。解析基因集合與樣品表型之間的聯(lián)系，繪制基因集合中基因之間的調(diào)控網(wǎng)絡并鑒定關鍵調(diào)控基因。
?
WGCNA適合于復雜的轉(zhuǎn)錄組數(shù)據(jù)，研究不同器官/組織類型和不同階段的發(fā)育調(diào)控、生物和非生物脅迫的不同時間點響應機制。
======WGCNA中的幾個概念======
共表達網(wǎng)絡：定義為加權基因網(wǎng)絡。點代表基因，邊代表基因表達相關性。加權是指對相關性值進行冥次運算(冥次的值也就是軟閾值 (power,
pickSoftThreshold這個函數(shù)所做的就是確定合適的power))。無向網(wǎng)絡的邊屬性計算方式為abs(cor(genex, geney)) ^ power；有向網(wǎng)絡的邊屬性計算方式為(1+cor(genex, geney)/2) ^ power; sign hybrid的邊屬性計算方式為cor(genex,geney)^power if cor>0 else 0。這種處理方式強化了強相關，弱化了弱相關或負相關，使得相關性數(shù)值更符合無標度網(wǎng)絡特征，更具有生物意義。如果沒有合適的power，一般是由于部分樣品與其它樣品因為某種原因差別太大導致的，可根據(jù)具體問題移除部分樣品或查看后面的經(jīng)驗值。
?
Module(模塊)：高度內(nèi)連的基因集。在無向網(wǎng)絡中，模塊內(nèi)是高度相關的基因。在有向網(wǎng)絡中，模塊內(nèi)是高度正相關的基因。把基因聚類成模塊后，可以對每個模塊進行三個層次的分析：1. 功能富集分析查看其功能特征是否與研究目的相符；2. 模塊與性狀進行關聯(lián)分析，找出與關注性狀相關度最高的模塊；3. 模塊與樣本進行關聯(lián)分析，找到樣品特異高表達的模塊。
?
Connectivity (連接度)：類似于網(wǎng)絡中"度"(degree)的概念。每個基因的連接度是與其相連的基因的邊屬性之和。
?
Module eigengene E: 給定模型的第一主成分，代表整個模型的基因表達譜。
?
Intramodular connectivity: 給定基因與給定模型內(nèi)其他基因的關聯(lián)度，判斷基因所屬關系。
?
Module membership: 給定基因表達譜與給定模型的eigengene的相關性。
?
Hub gene: 關鍵基因 (連接度最多或連接多個模塊的基因)。
?
Adjacency matrix(鄰接矩陣)：基因和基因之間的加權相關性值構成的矩陣。
?
TOM (Topological overlap matrix)：把鄰接矩陣轉(zhuǎn)換為拓撲重疊矩陣，以降低噪音和假相關，獲得的新距離矩陣，這個信息可拿來構建網(wǎng)絡或繪制TOM圖。
基本分析流程如下：
構建基因共表達網(wǎng)絡：使用加權的表達相關性。
識別基因集：基于加權相關性，進行層級聚類分析，并根據(jù)設定標準切分聚類結果，獲得不同的基因模塊，用聚類樹的分枝和不同顏色表示。
如果有表型信息，計算基因模塊與表型的相關性，鑒定性狀相關的模塊。
研究模型之間的關系，從系統(tǒng)層面查看不同模型的互作網(wǎng)絡。
從關鍵模型中選擇感興趣的驅(qū)動基因，或根據(jù)模型中已知基因的功能推測未知基因的功能。導出TOM矩陣，繪制相關性圖。
1. 構建基因關系網(wǎng)絡
1.1 計算基因間相關關系
基因間相似性（similarity）：根據(jù)基因在不同樣品中的表達情況，計算任意兩個基因間的相關關系。用Pearson相關系數(shù)
?
基因共表達矩陣：S=[Sij]
?
Sij 表示基因i和基因j的Pearson相關系數(shù)。
?
軟閾值：通過加權函數(shù)將相關系數(shù)變換，形成鄰接矩陣（Adjacency Matrix）,矩陣中元素連續(xù)化。
?
鄰接函數(shù)：power函數(shù)（冪指數(shù)函數(shù)）
?
aij=power(Sij, β)=|Sij|β
?
需要確定鄰接函數(shù)的參數(shù)β，依據(jù)無尺度網(wǎng)絡原則，即基因表達網(wǎng)絡符合無尺度網(wǎng)絡的冪函數(shù)分布。
1.2 無尺度網(wǎng)絡
網(wǎng)絡圖的點指圖中的每一個節(jié)點，度指與該點的連接數(shù)
?
隨機網(wǎng)絡（Random network），每個節(jié)點的度相對平均
?
無尺度網(wǎng)絡（Scale-free network），少數(shù)節(jié)點具有明顯高于一般點的度，這些點被稱為hub，由少數(shù)hub與其他節(jié)點關聯(lián)，最終構成整個網(wǎng)絡
?
無尺度網(wǎng)絡的冪率分布：節(jié)點連接數(shù)為k的節(jié)點數(shù)h，k與h成反比，負相關
?
尺度：隨機網(wǎng)絡中每個節(jié)點的連接數(shù)符合泊松分布，大部分節(jié)點的連接數(shù)居中，中值稱為隨機網(wǎng)絡的尺度。
?
無尺度網(wǎng)絡符合冪率分布，大多數(shù)點只有很少的連接，少數(shù)點有很多的連接
?
基因相關關系，冪函數(shù)處理后，少數(shù)強相關性不受影響或者影響較小，而相關性弱的取n次冪后，相關性明顯下降。
1.3 確定參數(shù)β
尋找合適的β，使得基因表達關系符合無尺度網(wǎng)絡，度數(shù)高的節(jié)點少，度數(shù)低的節(jié)點多。
節(jié)點度數(shù)k與具有該度數(shù)節(jié)點的個數(shù)h服從冪律分布
具體計算度數(shù)為k的節(jié)點個數(shù)的對數(shù)值log（k），與該節(jié)點出現(xiàn)的概率對數(shù)（log(p(k))）呈現(xiàn)負相關，一般會設置相關系數(shù)大于0.8
為了檢測設置的參數(shù)β是否滿足無尺度網(wǎng)絡，對log10(p(k))和log10(k)作圖，同時為更好評估，對兩者之間的相關系數(shù)做平方，即R2。如果模型R2接近1，則兩者之間為很好的線性關系。
1.4 計算基因間表達關系
評估基因間表達關系：直接關系?
生物體內(nèi)基因間的關系：直接關系+間接關系?
TOM：用拓撲重疊（topologicaloverlap measure，TOM）來計算基因之間關聯(lián)程度，除了分析兩個基因之間的關系，還考慮這兩個基因與其他基因之間的連接。這樣更具有生物學意義。
?
建立TOM矩陣：
TOM公式中，計算i與j之間的關系，不僅考慮了i和j的直接關系，還考慮了第三個基因μ的間接關系。
2 構建基因模塊
2.1 層次聚類樹
基因模塊的劃分基于基因間的連接稀疏性，將TOM矩陣（Similarity）轉(zhuǎn)化為相異度矩陣（Dissimilarity）
利用基于TOM值的相異度
層次聚類建樹
?
建樹方法：動態(tài)剪切樹和靜態(tài)剪切樹
2.2 動態(tài)混合剪切法
第一步：識別滿足設定條件的初級模塊
1.滿足模塊預定義的最低基因數(shù)目
2.距離集群過遠的基因，即使與集群處于同一分支，也去除
3.每個集群與其他周圍的集群顯著不同
4.處在樹分支尖端的每個群集的核心基因緊密相連
第二步：測試步驟
將未分配的基因進行測試，如果足夠接近某個初級群集，則分配進去
通常WGCNA使用動態(tài)混合剪切法建樹
2.3 建樹過程的參數(shù)
模塊最少基因數(shù)目（minModuleSize）
?
合并模塊的最小距離（mincutHeight）計算模塊的特征值，利用模塊特征值建樹，合并距離很近的模塊（如Height小于0.2）
?
模塊特征值（Epigengene）
?
模塊內(nèi)所有基因進行主成分分析（PCA），第一主成分的值即為Epigengene。它代表該模塊內(nèi)基因表達的整體水平。
?
3 篩選基因模塊
3.1 表達模式分析
模塊表達模式分析：模塊在各個樣品中的豐度
?
模塊特征值（Epigengene）：模塊內(nèi)所有基因進行主成分分析（PCA），第一主成分的值即為Epigengene。它代表該模塊內(nèi)基因表達的整體水平。
?
如果某模塊在樣品中特征值正或負表達較高，說明模塊與這個樣品關系緊密。
3.2 模塊與表型性狀關聯(lián)分析
模塊顯著性值（Module significance，MS）：模塊內(nèi)所有基因的基因顯著性值的平均值。
?
基因顯著性值（Gene significance, GS）:基因表達水平與因變量水平的相關系數(shù)。用T檢驗計算每個基因在不同表型樣品組間的差異表達顯著性檢驗P值（Pearson相關系數(shù)），通常將P值取以10底對數(shù)值定義為基因顯著性GS
?
計算各模塊與一表型性狀的MS值，如一個模塊的MS值顯著高于其他模塊，則這一模塊與該性狀存在關聯(lián)關系
?
模塊特征值顯著性（Epigengene significance, ES）:模塊特征值與某一性狀的相關系數(shù)，篩選與性狀關聯(lián)度最高的模塊。
3.3 富集分析
對各個模塊都進行GO和KEGG富集分析，找出與我們研究性狀相關通路相關性最強的模塊進行深入挖掘。
4.4 依據(jù)目標基因篩選模塊
依據(jù)研究目的、前期研究結果和已發(fā)表文獻，有重點關注的目標基因，可直接篩選目標基因所在的基因模塊重點進一步分析。
5 鑒定關鍵基因
5.1 模塊內(nèi)部基因連接度分析
Connectivity(degree)-連接度：與某個基因連接的所有其他基因的總和，即描述一個基因與其他所有基因的關聯(lián)程度，一般用K值表示。
?
Intramodular connectivity KIM-模塊內(nèi)部連接度IC:某個模塊中的基因與該模塊中其他基因的關聯(lián)程度（共表達程度）?？捎脕砗饬磕K身份（module membership,MM).
?
Module Membership MM,or Epigengene-basedconnectivity KME:模塊身份，用一個基因在所有樣本中的表達語與某個模塊特征值的表達譜的相關性，來衡量這個基因在這個模塊中的身份。
?
KME值接近0,說明這個基因不是該模塊的成員：KME接近1或者－1,說明這個基因與該模塊密切相關（正相關或者負相關）。
可以對所有基因計算相對某個模塊的KME值，并不一定要是該模塊的成員。
KME與KIM高度相關。某個模塊中KIM值高的hub基因一定與該模塊的KME也很高。
KME與KIM的區(qū)別：IC衡量基因在特定模塊中的身份，MM衡量基因在全局網(wǎng)絡中的位置。
?
篩選關鍵基因：
TOM值（模塊調(diào)控系表中的weight值）大于閾值（默認是0.15)的兩個基因才認為是相關的，然后計算每個基因的連接度。即先篩選有足夠強度的關系，然后計算連接度。
?
模塊內(nèi)部高連接度的基因，模塊內(nèi)排名前30或者10%(KME或KIM).
篩選關鍵基因：將該基因模塊身份MM相對于基因顯著性GS做散點圖，選擇右上角MM和GS均高的基因進一步分析。
基因顯著性值（Gene significance,GS)因變量水平的相關系數(shù)。衡量基因與表型性狀的關聯(lián)程度，GS越高，說明與表型越相關，越具有生物學意義。GS可以為正值或負值（正相關或負相關）
Cytoscape中一般用weight值（TOM值）來繪制網(wǎng)絡圖。
5.2 特定功能基因分析
高連通性的基因一般位于調(diào)控網(wǎng)絡的上游；低連通性的基因一般位于調(diào)控網(wǎng)絡的下游。
調(diào)控網(wǎng)絡上游一般是調(diào)控因子，如轉(zhuǎn)錄因子；下游一般是功能性的酶或蛋白分子。
重點關注具有調(diào)控功能的基因，典型的為轉(zhuǎn)錄因子，這些基因往往是關鍵基因。
5.3 目標基因關聯(lián)分析
依據(jù)研究目的，選取跟目標基因關系緊密的基因，如篩選與目標基因的TOM值排名前10，或者TOM值大于0.2的基因。
?
可準確篩選與目標基因存在上下游調(diào)控關系的候選基因。
?
當目標基因連接度不高時，可篩選與目標基因TOM值很高，且自身連接度也很高的基因。
===WGCNA安裝===
source("https://bioconductor.org/biocLite.R")
biocLite(c("AnnotationDbi","impute","GO.db", "preprocessCore"))
site="https://mirrors.tuna.tsinghua.edu.cn/CRAN"
install.packages(c("WGCNA","stringr", "reshape2"), repos=site)
本文使用 文章同步助手 同步
