書名：Microarray Data Analysis
編輯：Giuseppe Agapito
年份：2022
出版社：HUMAN PRESS
ISBN：ISBN 978-1-0716-1839-4

Pathway Enrichment Analysis of Microarray Data

作者

• Chiara Pastrello, Yun Niu, and Igor Jurisica

摘要

• 微陣列分析通常產(chǎn)生一個(gè)差異基因列表，這些基因需要進(jìn)行注釋以：

  • 鏈接到所研究的表型
  • 幫助規(guī)劃驗(yàn)證實(shí)驗(yàn)
  • 解釋結(jié)果

• 通路富集分析經(jīng)常用于上述目的乌妒，其中通路是人工創(chuàng)建的分子活動(dòng)和過程模型候味。

• 雖然有不同類型的通路富集分析可用垃你，本協(xié)議重點(diǎn)介紹最常見的類型——過度表示分析。

• 數(shù)據(jù)庫多樣性：

  • 許多數(shù)據(jù)庫收集不同的通路集
  • 為相同的通路策劃不同的基因集
  • 因此私杜，選擇最合適的通路來源進(jìn)行富集分析非常重要

• 綜合通路分析工具：

  • 本協(xié)議將使用pathDIP
  • pathDIP整合了22個(gè)主要通路數(shù)據(jù)庫
  • 支持全面的富集分析

• 可視化步驟：

  • 描述使用GSOAP進(jìn)行富集通路可視化的步驟

關(guān)鍵詞

• Pathway enrichment analysis, Pathway orphans, Pathway consolidation, Database coverage, Database overlap, Enrichment visualization

引言

• 微陣列分析的目標(biāo)是識(shí)別與研究條件相關(guān)的基因遭商。
• 完成預(yù)處理和歸一化后，常見步驟是識(shí)別不同條件下差異顯著的基因（如疾病與正常降狠、處理與未處理等）。
• 即使使用嚴(yán)格的統(tǒng)計(jì)過濾庇楞，研究人員通常也會(huì)得到一個(gè)較長的基因列表榜配，需要對(duì)其進(jìn)行優(yōu)先級(jí)排序以進(jìn)行后續(xù)驗(yàn)證，并結(jié)合表型進(jìn)行解釋吕晌。

基因注釋的多樣性

• 基因有多種注釋類型蛋褥，尤其在人類和模式生物中，這些注釋廣泛且豐富睛驳。

• 注釋類型包括：

  • 組織注釋（例如烙心，TissueEnrich [1]）
  • 疾病注釋（例如，Disease Ontology [2]）
  • 突變注釋（例如乏沸，COSMIC [3]）
  • 基因本體論（Gene Ontology [4]）
  • 通路注釋（最常用）

通路數(shù)據(jù)庫的發(fā)展

• 通路是研究人員創(chuàng)建的復(fù)雜生物事件的簡(jiǎn)化表示淫茵，用于總結(jié)和存儲(chǔ)當(dāng)前知識(shí)。
• 1990年代中期蹬跃，首兩個(gè)通路數(shù)據(jù)庫EcoCyc [5]和KEGG [6]被創(chuàng)建匙瘪。
• 機(jī)器可讀形式收集通路數(shù)據(jù)的重要性顯現(xiàn)。
• 根據(jù) pathguide.org 數(shù)據(jù)蝶缀，目前有超過兩百個(gè)數(shù)據(jù)庫用于策劃丹喻、收集、存儲(chǔ)和/或可視化特定或全面的通路集合翁都。
• 然而碍论，研究人員通常只使用少數(shù)幾個(gè)通路數(shù)據(jù)庫進(jìn)行注釋，選擇不總是基于“最佳”數(shù)據(jù)庫柄慰，因?yàn)椤白罴选痹诓煌瑘?chǎng)景和需求下定義不同鳍悠。

通路數(shù)據(jù)庫的異質(zhì)性與重疊

• 通路不是靜態(tài)實(shí)體，缺乏“真實(shí)值”作為參考先煎。
• 策劃自文獻(xiàn)的通路數(shù)據(jù)庫之間預(yù)期會(huì)有一定的重疊贼涩。
• 某些數(shù)據(jù)庫專注于特定反應(yīng)和信號(hào)級(jí)聯(lián)，如代謝通路（Biocyc [8]）或信號(hào)傳導(dǎo)通路（SignaLink [9]）薯蝎，因此不與其他類型數(shù)據(jù)庫重疊遥倦。
• 一般通路數(shù)據(jù)庫涵蓋廣泛類別，預(yù)期有較高的重疊率占锯，但實(shí)際最高重疊率僅約50% [7]袒哥。
• 數(shù)據(jù)庫在包含的通路數(shù)量、通路名稱和每個(gè)通路包含的基因上非常異質(zhì)消略，導(dǎo)致選擇合適的通路數(shù)據(jù)庫具有挑戰(zhàn)性堡称。

[圖片上傳失敗...(image-9b0764-1727822732738)]

<figcaption>image</figcaption>

基因覆蓋率的差異

• 通路數(shù)據(jù)庫在基因覆蓋率上存在差異：

  • Reactome [10]注釋約一萬個(gè)人類基因
  • KEGG注釋約七千個(gè)
  • WikiPathways [11]注釋約六千個(gè)

• Reactome和KEGG之間的最高基因重疊率為六千多個(gè)基因。

• 注釋缺失可能導(dǎo)致通路富集分析中的偏差和相關(guān)通路排名的差異 [12]艺演。

• 例如却紧，WNT信號(hào)通路在不同數(shù)據(jù)庫中的注釋數(shù)量和基因內(nèi)容差異顯著桐臊，導(dǎo)致富集分析結(jié)果差異較大。

[圖片上傳失敗...(image-5adce1-1727822732738)]

<figcaption>image</figcaption>

解決方案：通路數(shù)據(jù)庫的整合

• 整合不同的數(shù)據(jù)庫可以擴(kuò)展注釋基因數(shù)量和涵蓋的通路類型晓殊。
• 例如断凶，pathDIP [13]集成了22個(gè)通路數(shù)據(jù)庫，涵蓋17個(gè)物種巫俺，注釋人類超過13,000個(gè)基因到5,380個(gè)通路认烁。
• 雖然整合有助于解決覆蓋率問題，但也增加了返回的富集通路數(shù)量介汹，導(dǎo)致多重比較調(diào)整后顯著性p值急劇下降却嗡。
• 許多“相似”通路可能被返回，難以確定其具體來源和內(nèi)容相似性嘹承。

解決方案：通路整合與合并

• 通路整合不僅包括數(shù)據(jù)庫的整合窗价，還包括相似通路的合并，以減少通路數(shù)量并擴(kuò)展通路特定覆蓋叹卷。
• 每個(gè)原始數(shù)據(jù)庫中舌镶，相同通路可能有不同名稱，相似名稱可能指不同通路豪娜。
• 需要通過結(jié)合通路內(nèi)容和名稱相似性來分組功能相似的通路（例如餐胀，Pathcards [14]）。
• 最佳策略可能包括文本和內(nèi)容相似性的結(jié)合瘤载。

通路孤立基因的預(yù)測(cè)

• 即使經(jīng)過整合否灾，經(jīng)過策劃的數(shù)據(jù)庫最多覆蓋基因組的三分之二，許多基因可能無注釋鸣奔，成為“通路孤立基因”墨技。
• 預(yù)測(cè)通路注釋可以彌補(bǔ)這一空白。
• 在pathDIP中挎狸，使用策劃的通路進(jìn)行預(yù)測(cè)扣汪，增加了人類超過18,000個(gè)基因的覆蓋率，以及17個(gè)物種中超過120,000個(gè)基因的覆蓋率锨匆。

計(jì)算需求

• 本協(xié)議執(zhí)行的第一組分析將使用基于網(wǎng)絡(luò)的軟件崭别。
• 用戶需要互聯(lián)網(wǎng)連接和網(wǎng)頁瀏覽器。
• 對(duì)于可視化恐锣，我們將使用 R 中的 GSOAP茅主。

獲取差異表達(dá)基因

1. 為了從微陣列數(shù)據(jù)中獲取差異表達(dá)的基因，我們將使用 GEO (https://www.ncbi.nlm.nih.gov/geo)土榴，這是一個(gè)收集了數(shù)千個(gè)針對(duì)許多不同疾病的微陣列數(shù)據(jù)集的數(shù)據(jù)庫诀姚。
2. 擁有自己微陣列數(shù)據(jù)和差異表達(dá)基因列表的研究人員可以跳過此步驟。

執(zhí)行通路富集分析

1. 為了執(zhí)行通路富集分析玷禽，我們將使用 pathDIP 4.0 (http://ophid.utoronto.ca/pathDIP)赫段。
2. 有興趣的讀者可參考 http://ophid.utoronto.ca/pathDIP/API.jsp呀打，將此類分析集成到使用 Java、R 或 Python 的生物信息學(xué)工作流程中糯笙。

在 R 中使用富集分析結(jié)果

• 富集分析的結(jié)果將用于 R聚磺。

• 本協(xié)議中的示例在 R 4.0.3 (https://www.r-project.org) 上運(yùn)行。

• 需要安裝 GSOAP [15] 包炬丸。要安裝它，請(qǐng)?jiān)?R 中運(yùn)行以下命令：

  require(devtools)install_github("tomastokar/gsoap", dependencies=TRUE)
  

方法

3.1 從微陣列數(shù)據(jù)獲取差異基因

1. 連接到 GEO 并搜索疾病：

   • 研究人員在 GEO 網(wǎng)站（https://www.ncbi.nlm.nih.gov/geo）的搜索框中輸入所選疾病蜒蕾。在我們的案例中稠炬，輸入“骨關(guān)節(jié)炎”（“Osteoarthritis”），點(diǎn)擊搜索咪啡，選擇“GEO DataSets 數(shù)據(jù)庫中共有 1292 個(gè)“骨關(guān)節(jié)炎”結(jié)果首启。點(diǎn)擊數(shù)字將打開搜索結(jié)果頁面。

2. 篩選搜索結(jié)果：

   • 根據(jù)研究興趣篩選結(jié)果撤摸。在左側(cè)毅桃，選擇“Series”作為條目類型（見注1）和“Expression profiling by array”作為研究類型；在右側(cè)准夷，選擇“人類”（“Homo sapiens”）作為主要生物體钥飞。對(duì)于本協(xié)議，選擇數(shù)據(jù)集“GSE55584”衫嵌。

3. 使用 GEO2R 進(jìn)行差異基因分析：

   • 在數(shù)據(jù)集頁面读宙，點(diǎn)擊“Analyze with GEO2R”。這將打開用于獲取差異基因的頁面楔绞。

   • 設(shè)置分析中使用的組：

     • 點(diǎn)擊“define groups”结闸，在自由文本字段中輸入“OA”并按回車，然后輸入“RA”并按回車（見注2）酒朵。

     • 選擇適當(dāng)?shù)臉颖静⑵滏溄拥綄?duì)應(yīng)組：

       • 選擇所有“臨床狀態(tài)”為“骨關(guān)節(jié)炎”的樣本桦锄，點(diǎn)擊“OA”。
       • 選擇所有“臨床狀態(tài)”為“類風(fēng)濕性關(guān)節(jié)炎”（“Rheumatoid arthritis”）的樣本蔫耽，點(diǎn)擊“RA”结耀。

     • 定義組將自動(dòng)顯示每個(gè)標(biāo)簽的樣本數(shù)量（本例中為 6 和 10）。

   • 點(diǎn)擊頁面底部的“Analyze”按鈕匙铡，GEO2R 將運(yùn)行差異基因表達(dá)分析饼记。

4. 下載和篩選差異基因：

   • 分析結(jié)果顯示在頁面底部的表格中。
   • 選擇“Download full table”以獲取平臺(tái)上所有基因的結(jié)果（包含在數(shù)據(jù) S2 中）慰枕。
   • 篩選最顯著的基因具则，選擇調(diào)整后的 p 值 < 0.01，折疊變化大于 1 或小于 -1具帮，并具有基因符號(hào)注釋的基因（見注3）博肋。

3.2 通路富集分析

1. 在 pathDIP 中進(jìn)行富集分析：

   • 將前一步驟中獲得的 56 個(gè)基因符號(hào)列表粘貼到 pathDIP 的搜索框中低斋。列表可以由基因符號(hào)、Entrez 基因 ID 或 UniProt ID 組成匪凡。
   • 保持默認(rèn)選項(xiàng)或使用頁面上的篩選器膊畴。注意通路來源，通路數(shù)據(jù)庫根據(jù)存儲(chǔ)的通路類型進(jìn)行顏色編碼（例如病游，信號(hào)傳導(dǎo)唇跨、代謝、多種）（見注4）衬衬。
   • 在本分析中买猖，選擇“文獻(xiàn)策劃集”（literature curated set）。用戶也可以為注釋不足的基因選擇“擴(kuò)展通路關(guān)聯(lián)”（Extended pathway association）滋尉，并選擇是否僅使用實(shí)驗(yàn)蛋白質(zhì)相互作用或?qū)嶒?yàn)和高置信度預(yù)測(cè)蛋白質(zhì)相互作用生成的集玉控，以及預(yù)測(cè)的置信度水平。
   • 選擇是下載還是顯示結(jié)果狮惜。

2. 查看和下載富集分析結(jié)果：

   • pathDIP 將使用所選數(shù)據(jù)庫執(zhí)行過度表示分析（見注5）高诺，并返回三組結(jié)果：

     1. 注釋基因表：列出在 pathDIP 中具有注釋的基因，幫助用戶確認(rèn)所用集和數(shù)據(jù)庫的適用性（見表2）碾篡。
     2. 富集通路表：包含富集通路及其 p 值虱而、調(diào)整后的 p 值，以及指向原始數(shù)據(jù)庫的鏈接开泽。
     3. 基因-通路注釋表薛窥。

   • 如果選擇下載結(jié)果，三個(gè)表格將包含在一個(gè) txt 文件中眼姐。富集通路的列表通常用于生物學(xué)解釋（見注6）并作為出版物中的表格诅迷。

   • 頁面還允許用戶獲取通路/基因矩陣，并執(zhí)行詞匯富集众旗。詞匯富集結(jié)果以每個(gè)有意義的通路詞的 p 值和調(diào)整后的 p 值的表格形式展示罢杉，并提供用于詞云軟件（例如 Wordle—https://www.wordle.net）的詞列表。數(shù)據(jù) S4 顯示了本示例的詞匯富集結(jié)果贡歧。

[圖片上傳失敗...(image-a60ad8-1727822732738)]

<figcaption>image</figcaption>

3.3 GSOAP 繪圖

[圖片上傳失敗...(image-1d96b8-1727822732738)]

<figcaption>image</figcaption>

1. 準(zhǔn)備數(shù)據(jù)文件：

   • 從 pathDIP 獲得的通路列表將在 R 中使用滩租。

   • 為簡(jiǎn)化操作，將下載的文件分成兩部分：

     • 富集通路集：包含“Pathway Source”利朵、“Pathway Name”律想、“p-value”、“q-value (FDR: BH-method)”绍弟、“q-value (Bonferroni)”技即、“External Source ID”標(biāo)題的部分，命名為“pathways_to_be_plotted.txt”樟遣。
     • 基因-通路注釋：包含“UniProt”而叼、“Gene Symbol”身笤、“Entrez Gene”、“Pathway Source”葵陵、“Pathway Name”標(biāo)題的部分液荸，命名為“genes_annotation.txt”。避免使用 MS Excel 處理此部分（見注7）脱篙。

   • 在 R 中讀取文件：

     pathways_to_plot = read.delim("pathways_to_be_plotted.txt",stringsAsFactors = FALSE, header = TRUE)genes_annotations = read.delim("genes_annotation.txt",stringsAsFactors = FALSE, header = TRUE)
     

2. 格式化數(shù)據(jù)以符合 GSOAP 要求：

   • GSOAP 需要一個(gè)至少包含“Pathway”娇钱、“p-value”或“adjusted p-value”以及每個(gè)通路中識(shí)別出的基因符號(hào)以斜杠分隔的文件。

   • 操作步驟：

     1. 聚合注釋文件绊困，使基因按通路分組并以斜杠分隔：

        genes_annotations = aggregate(Gene.Symbol ~ Pathway.Source +Pathway.Name, genes_annotations, paste, collapse = "/")
        

     2. 合并兩個(gè)文件：

        pathways_to_plot = merge(x = pathways_to_plot, y = genes_annotations, by = c("Pathway.Source", "Pathway.Name"), all.x = TRUE)
        

     3. 篩選 p 值經(jīng)過 Bonferroni 校正 < 0.01 的通路文搂，并設(shè)置數(shù)據(jù)框的行名為通路名稱：

        pathways_to_plot = pathways_to_plot[which(pathways_to_plot$q.value..Bonferroni. < 0.01),]rownames(pathways_to_plot) = make.names(pathways_to_plot[,2], unique = TRUE)
        

3. 使用 GSOAP 計(jì)算和繪圖：

   • GSOAP 將使用 Jaccard 距離計(jì)算通路之間的基因重疊，并使用通路顯著性（計(jì)算為 -log10(provided p-value)）來計(jì)算通路之間的接近度考抄。顯著性還用于對(duì)通路進(jìn)行排序，因此在圖中只顯示最顯著的標(biāo)簽（本例中為 5 個(gè)）蔗彤。

   • 使用以下命令進(jìn)行計(jì)算川梅、排序和繪圖：

     layout_path = gsoap_layout(pathways_to_plot, 'Gene.Symbol', 'q.value..FDR..BH.method.')layout_path = layout_path[order(layout_path$significance, decreasing = TRUE),]gsoap_plot(layout_path, as.color = 'cluster', as.alpha = 'significance', which.label = 1:5)
     

   • 圖2顯示了本示例的結(jié)果。

注釋

1. 注1：Entry type 中的 “Series” 指的是 GEO 數(shù)據(jù)集中的系列類型然遏，用于描述一組相關(guān)的實(shí)驗(yàn)數(shù)據(jù)贫途。
2. 注2：定義組時(shí)，確保正確標(biāo)注每個(gè)樣本所屬的組待侵，以確保差異分析的準(zhǔn)確性丢早。
3. 注3：Gene Symbol 注釋確保基因列表中的每個(gè)基因都有對(duì)應(yīng)的標(biāo)準(zhǔn)符號(hào)秧倾，便于后續(xù)分析和解讀怨酝。
4. 注4：通路數(shù)據(jù)庫的顏色編碼有助于快速識(shí)別不同類型的通路，提升分析效率那先。
5. 注5：過度表示分析（Over-Representation Analysis, ORA）是一種常用的通路富集方法农猬，通過檢測(cè)目標(biāo)基因集中某些通路中的基因是否顯著多于預(yù)期來識(shí)別富集通路。
6. 注6：富集通路的生物學(xué)解釋有助于理解差異基因在生物過程中的潛在作用和機(jī)制售淡。
7. 注7：避免使用 MS Excel 處理基因-通路注釋文件斤葱，以防止數(shù)據(jù)格式和內(nèi)容的意外更改。
8. 注8：設(shè)置數(shù)據(jù)框的行名為通路名稱有助于后續(xù)的可視化和數(shù)據(jù)操作揖闸，確保每個(gè)通路都有唯一標(biāo)識(shí)揍堕。

筆記

1. 注1：GEO 系列是由原始提交者提供的數(shù)據(jù)集。當(dāng) GEO 員工策劃該數(shù)據(jù)時(shí)汤纸，它就成為了一個(gè)策劃過的 GEO 數(shù)據(jù)集衩茸。每個(gè)系列和數(shù)據(jù)集都包括在特定平臺(tái)（本例中為 Illumina HumanHT-12 V3.0 表達(dá)珠芯片）上運(yùn)行的一組樣本（其 ID 以 GSM 開頭）。
2. 注2：設(shè)置組時(shí)使用的順序很重要贮泞，因?yàn)?GEO2R 對(duì)所選組不敏感递瑰。分析將始終獨(dú)立于組的名稱執(zhí)行“第一組輸入” vs “第二組輸入”祟牲。
3. 注3：在微陣列中，一個(gè)探針集應(yīng)設(shè)計(jì)為靶向一個(gè)基因抖部，但當(dāng)一個(gè)基因家族包含保守區(qū)域時(shí)说贝，同一個(gè)探針集可能會(huì)靶向不同的相關(guān)基因。另一方面慎颗，按設(shè)計(jì)乡恕，一個(gè)基因被多個(gè)探針集靶向，這些探針集對(duì)基因的不同區(qū)域?qū)R俯萎。這可能導(dǎo)致不同的探針集靶向基因的不同轉(zhuǎn)錄本傲宜，可能具有不同的表達(dá)水平 [14]。在某些情況下夫啊，探針集不映射到已知基因函卒，需要從基因注釋分析等基于基因的分析中排除。尤其重要的是撇眯，在出版物中跟蹤并提供生成感興趣基因列表的探針集报嵌。
4. 注4：必須注意用戶正在處理的數(shù)據(jù)類型。與通路數(shù)據(jù)庫類似熊榛，不同的微陣列具有不同的基因覆蓋率锚国。一些旨在靶向整個(gè)基因組（例如 Affymetrix Human Genome U133 Plus 2.0 Array），而另一些則專門針對(duì)一組基因（例如 CustomArray Signaling Pathway platform version 3 是針對(duì)關(guān)注信號(hào)傳導(dǎo)通路的約 2000 個(gè)基因的陣列 [16]）玄坦。在后一種情況下血筑，使用僅限于信號(hào)傳導(dǎo)特定數(shù)據(jù)庫進(jìn)行富集分析會(huì)更合適。同樣煎楣，取決于實(shí)驗(yàn)設(shè)計(jì)豺总，某個(gè)數(shù)據(jù)庫可能比其他數(shù)據(jù)庫更合適。例如择懂，專注于代謝的研究將受益于使用代謝特定數(shù)據(jù)庫园欣，而研究與疾病相關(guān)的基因則可能受益于使用通用數(shù)據(jù)庫。
5. 注5：通路富集分析可以使用不同的技術(shù)進(jìn)行休蟹，這些技術(shù)可以分為三大類：過度表示分析（ORA）沸枯、基因集富集分析（GSEA）和基于拓?fù)涞耐犯患治觯═PEA），如 [17] 中所述并進(jìn)行比較赂弓。簡(jiǎn)而言之绑榴，ORA 是最常用的方法，因?yàn)樗?jiǎn)單盈魁；它只需要一個(gè)感興趣基因列表翔怎，計(jì)算列表中屬于某個(gè)通路的基因比例（通常使用超幾何檢驗(yàn)或 Fisher 檢驗(yàn)）。GSEA 考慮通路內(nèi)的所有表達(dá)變化，以提供該通路受到影響的概率赤套。提供的基因列表需要進(jìn)行排序飘痛，以及排序值（例如，折疊變化）容握。TPEA 涵蓋使用通路內(nèi)基因之間的連接來計(jì)算通路富集和激活狀態(tài)的方法宣脉。這些方法特別容易受到數(shù)據(jù)庫選擇的影響，因?yàn)橥負(fù)鋵釉黾恿水愘|(zhì)性 [12]剔氏。
6. 注6：獲得的通路數(shù)量可以是可管理的塑猖，如本例，或可能過多谈跛，如 [18]羊苟。在后一種情況下，將單個(gè)通路映射到更高層次的通路本體（例如感憾，Pathway Ontology [19] 或 Reactome）可以幫助減少復(fù)雜性并突出顯示感興趣的過程蜡励。
7. 注7：在第一個(gè)表格中，因?yàn)槲覀儍H處理通路名稱阻桅，任何工具/軟件都可以使用凉倚，包括 MS Excel。然而鳍刷，由于第二個(gè)表格包含基因-通路對(duì)注釋占遥，MS Excel 不應(yīng)被使用俯抖，因?yàn)樗鼤?huì)默認(rèn)將特定的基因符號(hào)轉(zhuǎn)換為日期 [20]输瓜。
8. 注8：設(shè)置行名的常用命令是 rownames。然而芬萍，在來自不同數(shù)據(jù)庫的通路名稱的情況下尤揣，可能有少數(shù)不同的通路具有相同的名稱，這會(huì)導(dǎo)致 R 中出現(xiàn)錯(cuò)誤柬祠。因此北戏，我們使用 make.names 命令，它將為行賦予與通路相同的名稱漫蛔，除非有兩個(gè)（或更多）通路具有相同的名稱（如本例中的 ‘Chemokine signaling’）嗜愈。在后一種情況下，具有相同名稱的行名將帶有索引以進(jìn)行區(qū)分（在本例中莽龟，一個(gè)行名為 ‘Chemokine.signaling’蠕嫁，另一個(gè)為 ‘Chemokine.signaling.1’）。