Ⅰ GO富集分析
GO富集是組學數據分析常用的手段，通常用來挖掘差異基因中GO term的富集程度。Fisher's exact test是常用的統計檢驗方法，但這種方法存在明顯的缺點。很多公司提供的測序分析結果都普遍使用這樣的方法，導致很多後續的分析與實驗結果不一致的情況。對於這種情況，目前還有其他演算法來彌補這些缺點。（文中例子來源於《the Gene Ontology handbook》）
fisher's exact test是基於超幾何分布來計算的，單邊檢驗就是超幾何檢驗。通常用來檢驗兩組分類是否有顯著差異。
m：研究物種的基因數；
n：研究的樣本中基因數；
mt：總體中被注釋到term t（GO 詞條t）的基因數；
nt：樣本中被注釋到term t的基因數。
隨機變數Xt表示樣本中被觀察到的term t 的數目，所有根據超幾何分別，觀察到k個term t 的概率P(Xt)是：
零假設H0：樣本中出現的term t 的數目與總體中總的term t數目沒有正關聯。也就是說樣本中的term t數目的比例與總體中term t的接近。
為了拒絕H0，使用單尾檢驗：
一個簡單的例子：假設總體中有18個基因，其中有5個被注釋到binding這個term，轉錄組分析發現有5個差異表達的基因，其中有3個被注釋到binding這個term，為了說明binding這個term是否是overrepresentation，用上面的Fisher's exact test計算p值：
在現實中，我們不可能只對某一個term進行檢驗，而是對很多term進行檢測，即多重檢驗，但這樣就會導致假陽性的term數目非常高。所以我們需要對p值進行校正。通常是使用Benjamini-Hochberg校正方法來控制預期的錯誤發現率（false discoveries rate-FDR）進行校正。( 如何通俗地解釋錯誤發現率（FDR） )。盡管多重檢驗的校正可以減少假陽性，但並不能從根本上解決GO（或KEGG）富集的問題。
GO富集的根本問題在於一個基因對應的GO term有多個，一個term對應多個gene，同時還有層級關系。這樣導致如果一個term顯著富集，那和它共享很多基因的term也會顯著富集。
有很多其他的演算法來試圖解決這個問題，其中包括parent-child approach、topology-based algorithms、model-based approaches和gene set enrichment analysis。下面是對這些演算法的簡單介紹：
該演算法還是基於Fisher's exact test，只不過考慮了term的父節點。在計算概率時，不再是在總體m中取樣，而是從term的父節點中取樣，所以計算公式變成了：
當一個term有多個父節點時計算就變得復雜了，具體方法還得參考原始文獻（improved detection of overrepresentation of Gene-Ontology annotatins with parent child analysis）
這兩種方法原理反正我沒看懂，有興趣的可以看原始文獻：
1、Improved scoring of functional groups from gene expression data by decorrelating GO graph structure.
2、GOing Bayesian: model-based gene set analysis of genome-scale data
該演算法首先根據感興趣的特徵（比如差異基因的表達量）對基因進行排序，形成一個列表。零假設是某個基因集（genes encoding procts in a metabolic pathway, located in the same cytogenetic band, or sharing the sam GO category）里基因順序與這個列表沒有關聯，即排序是隨機的。對應的備擇假設是它們之間有關聯。如果基因集里的基因都聚集在基因列表的前端或底端或者非隨機分布，我們就傾向於相信它們之間有關聯。
S：想研究的基因集；
L：整個排序的基因列表；
統計量：Kolmogorov-Smirnov(KS)
step 1:計算富集得分（Enrichment Score）。按順序從頭到尾逐個比較L中的基因與S中的基因，加和統計量，如果兩者相同就增加KS統計量，反之就減少KS統計量。增加的多少與這個基因和表型的相關性有關。最後ES就是KS的最大方差值。
step 2、檢驗ES的顯著性。重復k次隨機選擇的大小為nt的基因集（Nt1,...,Ntk)，p值計算公式為：
step3、使用FDR進行多重檢驗的校正。
相關軟體：GSEA-P software
1、Benjamini Y, Hochberg Y. Controlling the false discovery rate: a practercal and powerful approch to multiple testing.
2、 https://www.hu.com/question/3560619
3、Gene set enrichment analysis: A knowledge-based approach for interpreting genome-wide expression profiles
4、 《the Gene Ontology handbook》
Ⅱ 「GO富集分析」從原理到實踐 ~ 零基礎掌握
原本，我並無寫這一稿件的想法。主要原因有二：
如果要找合理解釋，那麼針對第一點，就是每天仍然有大量新接觸生信數據分析的朋友；針對第二點，......在前兩天我推的文稿《零基礎快速完成基因功能注釋 / GO / KEGG / PFAM...》中，評論區答應了下，閱讀過5000，那就寫一寫富集分析。於是，如果不寫，總是不對。如果要寫，只能現在寫。畢竟有些事情，現在不做，以後真的不會做。
對於這一塊，完全陌生的朋友，尤其是不少生物學背景朋友，有必要溫習一下數理統計基礎。這一稿件只做原理最簡單的但使用最廣泛其速度最快的Over-Represence Analysis模式的富集分析講演。其他模式，不涉及。
回到主題，先舉個經典的抽球例子：
小紅小綠小藍三個人自稱有超能力，可以用手摸摸球就分辨出黑球白球，於是我們找來黑袋子，放100個球，其中20個白球80個黑球，讓三人分別無放回地抽取。
小紅隨機抽出來10個球，其中2個白球8個黑球，情況即，
抽球中白球比例與背景白球比例完全一致，說明小紅抽球結果隨機。
球放回去，小綠來抽球，抽出來的10個球，其中3個白球7個黑球，情況即，
這是經典的抽球案例，抽取到的白球個數的概率分布為超幾何分布。基於此，我們可以簡單計算抽取到比小綠抽取到球個數（或更多即更極端）的概率如何，在 R語言中計算，即
而對於小藍的情況，那麼概率如何？
在 TBtools 中也可以計算，只是寫法有點區別
可以看到，盡管這只是一次抽球，小綠抽球中白球比例（或更極端情況）出現的概率是31.88%+，還是挺高的，於是我們有較高的把握說，小綠嘛，只是走了狗屎運。相反，小藍抽球中白球比例或更極端情況出現的概率幾乎為 0 ，我們幾乎沒啥把握說，小藍走狗屎運....換句話說，我們有理由相信，或許小藍真有抽白球的超能力.....
說了這么多，那麼跟基因集合富集分析有啥關系？....基因集合功能富集分析。那麼我們就需要有一個基因集合（如差異表達基因集合或ChIP-seq的Peaks或GWAS定位的系列區間），還有一個功能標簽（如 生長素信號轉導相關 ）。於是黑白球案例可以簡單調整一下。假定現在這個物種一共有100個基因，其中20個基因與生長素信號轉導相關，80個沒有注釋到與生長素信號轉導相關（換句話說，約等於無關），我們做了對植株做了處理，和CK分別測定轉錄表達譜，通過差異表達分析，鑒定到10個差異表達基因，其中2個與生長素信號轉導相關，而另外8個則沒注釋到生長素信號轉導相關，簡單畫一下，即
好，剩下的兩個就不替換了。整體上，ORA模式的富集分析，本身就是經典的抽球案例，感興趣的自行替換就可以了。
基本原理，相信都搞清楚了。不過還是有兩三點需要注意：
具體如何做物種所有基因的背景注釋，請參考前述推文《零基礎快速完成基因功能注釋 / GO / KEGG / PFAM...》。
首先，打開 TBtools GO 富集分析界面
整體如上，一共三個文件：
具體示例如下
點擊 Start ，隨後等待即可。完成時會有彈窗提示。查看輸出文件
（寫到這里，突然覺得這些都沒啥意思，不知為何....就不詳細寫了，大夥自己看看列名，猜猜吧）
很多時候，我們會選擇，篩選第一列，只看 Biological Process。一般這些與我們的生物學認知會貼近一些。
基因集合功能富集分析，是一個常常被談起的話題，甚至近期都有不少新方法或演算法被提出。感興趣的朋友可以去了解。這份教程，只與大夥說最簡單，但也是使用最為廣泛的一種富集分析模式。無論是不是 TBtools 用戶，理論上來說，都可以輕松理解並掌握，從原理到實踐。
寫到一半，其實我已經不想寫了。原因非常簡單，這也是為什麼在我之前，並沒有一個人寫出來 TBtools 類似的工具。不是寫不了，而是不想寫。有時候，隨著能力增長和知識積累，往往不再願意做一些簡單的事情。或許這還涉及到年齡的增長，角色的轉變，責任的變化....雲雲。
小時候，我以為寫 TBtools 玩玩；
後來，我以為我會一直寫下去；
現在，，，，，，
Ⅲ 什麼是GO富集分析，常說的GO功能分析、功能分析、Pathway分析是什麼意思
Gene Ontology可分為分子功能（Molecular Function），生物過程（biological process）和細胞組成（cellular component）三個部分。蛋白質或者基因可以通過ID對應或者序列注釋的方法找到與之對應的GO號，而GO號可對於到Term，即功能類別或者細胞定位。
功能富集分析: 功能富集需要有一個參考數據集，通過該項分析可以找出在統計上顯著富集的GO Term。該功能或者定位有可能與研究的目前有關。
GO功能分類是在某一功能層次上統計蛋白或者基因的數目或組成，往往是在GO的第二層次。此外也有研究都挑選一些Term，而後統計直接對應到該Term的基因或蛋白數。結果一般以柱狀圖或者餅圖表示。
1.GO分析
根據挑選出的差異基因，計算這些差異基因同GO 分類中某（幾）個特定的分支的超幾何分布關系，GO 分析會對每個有差異基因存在的GO 返回一個p-value，小的p 值表示差異基因在該GO 中出現了富集。
GO 分析對實驗結果有提示的作用，通過差異基因的GO 分析，可以找到富集差異基因的GO分類條目，尋找不同樣品的差異基因可能和哪些基因功能的改變有關。
2.Pathway分析
根據挑選出的差異基因，計算這些差異基因同Pathway 的超幾何分布關系，Pathway 分析會對每個有差異基因存在的pathway 返回一個p-value，小的p 值表示差異基因在該pathway 中出現了富集。
Pathway 分析對實驗結果有提示的作用，通過差異基因的Pathway 分析，可以找到富集差異基因的Pathway 條目，尋找不同樣品的差異基因可能和哪些細胞通路的改變有關。與GO 分析不同，pathway 分析的結果更顯得間接，這是因為，pathway 是蛋白質之間的相互作用，pathway 的變化可以由參與這條pathway 途徑的蛋白的表達量或者蛋白的活性改變而引起。而通過晶元結果得到的是編碼這些蛋白質的mRNA 表達量的變化。從mRNA 到蛋白表達還要經過microRNA 調控，翻譯調控，翻譯後修飾（如糖基化，磷酸化），蛋白運輸等一系列的調控過程，mRNA 表達量和蛋白表達量之間往往不具有線性關系，因此mRNA 的改變不一定意味著蛋白表達量的改變。同時也應注意到，在某些pathway 中，如EGF/EGFR 通路，細胞可以在維持蛋白量不變的情況下，通過蛋白磷酸化程度的改變（調節蛋白的活性）來調節這條通路。所以晶元數據pathway 分析的結果需要有後期蛋白質功能實驗的支持，如Western blot/ELISA，IHC（免疫組化），over expression（過表達），RNAi（RNA 干擾），knockout（基因敲除），trans gene（轉基因）等。
3.基因網路分析
目的：根據文獻，資料庫和已知的pathway 尋找基因編碼的蛋白之間的相互關系(不超過1000 個基因)。
Ⅳ RNA-Seq(9)：使用GSEA做GO/KEGG富集分析
最廣為人知的富集分析做法是把上調、下調基因分別或者合並，拿來做GO和KEGG富集分析。經常有一些數據集，拿差異基因做得不到結果，那是因為確實富集不到任何通路，是正常的。不妨試試GSEA，不是拿差異基因，而是拿全部基因作為輸入。
GSEA與GO，KEGG分析區別：GO，KEGG分析更加依賴差異基因，實則是對一部分基因的分析 (忽略差異不顯著的基因)，而GSEA是從全體基因的表達矩陣中找出具有協同差異 (concordant differences)的基因集，故能兼顧差異較小的基因
GO，KEGG富集是定性的分析，GSEA考慮到了表達或其它度量水平的值的影響。GSEA分析不需要指定閾值（p值或FDR）來篩選差異基因，在沒有經驗存在的情況下分析我們感興趣的基因集，而這個基因集不一定是顯著差異表達的基因。GSEA分析可以將那些GO/KEGG富集分信息中容易遺漏掉的差異表達不顯著卻有著重要生物學意義的基因包含在內。
另外，對於時間序列數據或樣品有定量屬性時，GSEA的優勢會更明顯，不需要每個分組分別進行富集，直接對整體進行處理。
數據准備,製作geneList
我們現在知道Cytokine-cytokine receptor interaction setSize enrichmentScore是被抑制的，如果還想看一下這個通路裡面的基因是如何變化的，應該怎麼辦呢，pathview 可以幫到我們。