生信統(tǒng)計(jì),利用TCGA數(shù)據(jù)庫(kù)發(fā)5分的套路修煉

1、生信統(tǒng)計(jì)，利用TCGA數(shù)據(jù)庫(kù)發(fā)5分的套路修煉科研經(jīng)驗(yàn)|文獻(xiàn)|實(shí)驗(yàn)|工具| SCI寫作|國(guó)自然作 者：酸菜轉(zhuǎn)載請(qǐng)注明：解螺旋臨床醫(yī)生科研成長(zhǎng)平臺(tái)話說(shuō)天下大勢(shì)，粉久必冷。為了溫暖你們的小心臟，今日祭由“不做實(shí)驗(yàn)卻想發(fā) SCI”的一篇生信+統(tǒng)計(jì)的文章教學(xué)。是不是看到題目就興奮的不要不要的？不做實(shí)驗(yàn)的理由千千萬(wàn)，核心一條就是沒(méi)錢。醫(yī)生做科研， 做實(shí)驗(yàn)，學(xué)統(tǒng)計(jì)，敲代碼，你總得會(huì)一樣。統(tǒng)計(jì)學(xué)深了那是 數(shù)學(xué)問(wèn)題，酸菜自從高考數(shù)學(xué)差一分及格之后就受到了永久 性的心理創(chuàng)傷，用SPSS做做單因素和多因素分析已到極限， 剩下的就只有生信這一條路來(lái)嘗試無(wú)本生意了。今天的教學(xué)案例 2017年新發(fā)表于 Oncotarge

2、t (IF = 5.168 ), 一本讓人愛(ài)恨交織的雜志，一年收稿三千多篇，沒(méi)發(fā)過(guò)的人 都說(shuō)它很水，內(nèi)心卻是弱水三千，吾亦想取一瓢的挽扎糾結(jié)。能發(fā)一篇Oncotarget不錯(cuò)了，要啥自行車。這篇文章的創(chuàng)新在于揭示了EB病毒基因與人類胃癌轉(zhuǎn)錄組的交互作用。EB病毒同胃癌的關(guān)系文獻(xiàn)早有報(bào)道，研究分 子機(jī)制是常規(guī)思路，而用生物信息學(xué)、統(tǒng)計(jì)學(xué)方法從TCGA(The Cancer Genome Atlas) 數(shù)據(jù)庫(kù)挖掘到特征性分子標(biāo)記 進(jìn)行EB病毒與人基因的交互網(wǎng)絡(luò)分析， 則是新穎的切入點(diǎn)， 這個(gè)模式借鑒一下，微生物 +癌癥基因組合，潛在的資源頗豐。文章總體思路是先找到EBV(+)和EBV(-)的胃癌中

3、有差異表達(dá)的基因和通路，隨后分析它們和EBV基因的相關(guān)性，構(gòu)建交互作用網(wǎng)絡(luò)。技術(shù)路線分三大步驟，作者用了一張圖來(lái) 概括：Figure4數(shù)據(jù)分析流程 第一步是數(shù)據(jù)特征提取 (A), 第二步是單因素相關(guān)分析(B),第三步是多因素相關(guān)分析(C)o單因素和多因素分析是臨床研究的入門統(tǒng)計(jì)技術(shù)了，關(guān)鍵性 難度在于測(cè)序數(shù)據(jù)分析方面，基因組數(shù)據(jù)本來(lái)就復(fù)雜，這里 又需要分析人和病毒兩個(gè)轉(zhuǎn)錄組，其中的多重檢驗(yàn)會(huì)導(dǎo)致很 多的假陽(yáng)性率，所以需要降維和校正。第一步：提取特征分子標(biāo)志從TCGA數(shù)據(jù)庫(kù)下載EBV(+)和EBV(-)的原始測(cè)序fastq數(shù) 據(jù)，有25個(gè)EBV(+)和260個(gè)EBV(-)樣本。然后構(gòu)建一個(gè)臨 床

4、特征表格，包括性別、年齡、腫瘤部位、組織學(xué)病理學(xué)診 斷等等，接著用 R語(yǔ)言的dist()函數(shù)計(jì)算兩組樣本間臨床特 征的歐式距離，找由距離最近的樣本進(jìn)行匹配。因?yàn)橛袔讉€(gè)EBV(-)樣本跟EBV(+)樣本的距離是一樣的，所 以最終選到了 20個(gè)EBV(-)樣本，來(lái)跟那25個(gè)EBV(+)樣本 配對(duì)。這樣每個(gè)配對(duì)之間的距離，就比 25個(gè)EBV(+)對(duì)260 個(gè)EBV(-)之間的平均距離小多了，控制了混雜因素，使組間 具有可比性。該文章附件中的表格展示了配對(duì)樣本的一部 分，左邊4列分別為EBV(+)樣本編號(hào)、組織學(xué)診斷、EBV(-) 樣本編號(hào)、組織學(xué)診斷。右邊 2列數(shù)字就是歐式距離，左邊 的是EBV(+

5、)樣本與它所匹配的 EBV(-)樣本的距離，右邊的 是該EBV(+)樣本與原260個(gè)EBV(-)樣本的平均距離，可見 前者是比后者要小多，消除了樣本間差異的偏倚因素。接下來(lái)對(duì)EBV基因進(jìn)行篩選，從原始的 88個(gè)基因中，將在 超過(guò)5個(gè)樣本中原始計(jì)數(shù)(raw count)為0的基因去除，剩 下19個(gè)基因，這就獲得了 EBV相關(guān)的特征分子標(biāo)志。部分EBV基因在EBV(+)樣本中的表達(dá) 然后提取人類基因 組的特征分子標(biāo)志。用Bioconductor的DEseq2包找由兩組樣本間差異表達(dá)的基因(DEx gene),這里就要用到 Bonferroni法對(duì)P值進(jìn)行校正，原始 P值達(dá)到1E-6,即校正 后的P

6、再用 Bioconductor 中的 GAGE (Generally Applicable Gene-set Enrichment for Pathway Analysis)分析法，找到兩組樣本中差 異表達(dá)的通路，共計(jì) 29個(gè)，在EBV(+)樣本中全是上調(diào)的。 同時(shí)用R語(yǔ)言的MEGENA包構(gòu)建那些DEx基因的共表達(dá)網(wǎng) 絡(luò)，找到它們的聚類基因模塊(Gene modules)和其中的關(guān)鍵節(jié)點(diǎn)基因(hub genes)。前者modules是指一組具有較強(qiáng)相 關(guān)性的基因，而hub則是指modules中具有全局性影響的樞 紐基因。這一步是利用人類基因組中本來(lái)就有的交互網(wǎng)絡(luò)數(shù) 據(jù)，來(lái)減少數(shù)據(jù)維度，聚焦到關(guān)

7、鍵信息。這一步找到了91個(gè)節(jié)點(diǎn)基因和27個(gè)基因模塊（一個(gè)節(jié)點(diǎn)基因可由現(xiàn)在不同 的模塊中）。然而，這些數(shù)據(jù)維度還是高，怎么進(jìn)一步抓取核心要素呢？作者又用了主成分分析法（PCA）提取其中起主要作用的信 息，他們只提取了第一個(gè)主成分（ PC1）來(lái)代表該模塊或通 路。這樣，網(wǎng)絡(luò)構(gòu)建和相關(guān)性分析的素材提取完畢，最終獲得19個(gè)EBV基因、91個(gè)人類胃癌相關(guān)的 hub基因、27個(gè) 基因模塊和29個(gè)通路（的PC1）。對(duì)于生信有基礎(chǔ)的同學(xué)，可按照上述方法利用工具自行探索，零基礎(chǔ)的同學(xué)我們還有后續(xù)教學(xué)。第二步，單因素相關(guān)分析單因素相關(guān)分析就是要從19個(gè)EBV基因、91個(gè)人類胃癌相關(guān)的hub基因、27個(gè)基因模塊和2

8、9個(gè)通路中再進(jìn)一步獲取 相關(guān)性最高的組合，采用配對(duì)的Pearson相關(guān)性檢驗(yàn)對(duì)三組數(shù)據(jù)分別進(jìn)行分析：1）基因模塊的PC1和EBV基因表達(dá)數(shù) 據(jù)；2）差異通路的PC1和EBV基因表達(dá)數(shù)據(jù)；3） hub基 因和EBV基因表達(dá)數(shù)據(jù)。用置換檢驗(yàn)法獲取相關(guān)系數(shù)（用C + + ） o至于相關(guān)系數(shù)的顯著性，則沒(méi)有采用通用的p這一步得到了 7個(gè)EBV基因和12個(gè)模塊的PC1呈顯著相關(guān)， 其中相關(guān)性最強(qiáng)的是 BLAF1基因和12號(hào)模塊，相關(guān)系數(shù)達(dá) 到0.662,如下表。3個(gè)EBV基因和4個(gè)人類胃癌的hub基因相關(guān)，其中最強(qiáng)的是 LMP-1和Clorf115,相關(guān)系數(shù)0.754。 又有2個(gè)EBV基因和4個(gè)通路的P

9、C1相關(guān)，其中BLAF4和 4個(gè)通路都相關(guān)，最強(qiáng)的是 BLAF4和磷脂酰肌醇信號(hào)通路， 相關(guān)系數(shù)達(dá)到0.709。好了，單因素分析這步到這里做完。第三步，多因素關(guān)聯(lián)分析這步要用到稀疏典型相關(guān)分析法(sparse CanonicalCorrelation Analysis, sCCA )。典型相關(guān)分析法(CCA )是分 析兩個(gè)數(shù)據(jù)矩陣的經(jīng)典方法，不需要降維提取主成分。但普 通CCA只適用于行比列多的矩陣，即樣本比基因多，這情 況顯然不對(duì)，是樣本少基因多，所以要用 sCCA。這步用R 語(yǔ)言的PMA包來(lái)做。其中用CCA()函數(shù)可得到每個(gè)元素(比如每個(gè)EBV基因)的典型系數(shù)(canonical coef

10、ficient ),它表示該元素對(duì)整個(gè)矩 陣典型相關(guān)系數(shù)的貢獻(xiàn)，如果莫個(gè)元素的典型系數(shù)非零，則 表示它對(duì)全局相關(guān)性有很大貢獻(xiàn)，被選為核心(essential)基因。這步得到22個(gè)基因模塊與19個(gè)EBV基因計(jì)數(shù)矩陣相關(guān)(P 此外，人類hub基因計(jì)數(shù)矩陣和EBV基因計(jì)數(shù)矩陣的典型 相關(guān)系數(shù)是0.806。對(duì)這個(gè)相關(guān)性貢獻(xiàn)最大的核心EBV基因和核心人類hub基因就不一一列由了。所有的差異表達(dá)的通 路都和EBV基因計(jì)數(shù)矩陣達(dá)到了顯著相關(guān)的典型系數(shù)。最 強(qiáng)幾個(gè)的如下表所示。至此，數(shù)據(jù)分析的主要步驟就做完了。但文章沒(méi)完呀，還要做個(gè)總結(jié)，比較一下 Pearson和sCCA 兩種相關(guān)性分析的結(jié)果。咱們接著看。第

11、四步，數(shù)據(jù)匯總先總結(jié)一下EBV基因的情況。如下表，各個(gè)重要EBV基因在Pearson和sCCA分析中，分別與人類胃癌基因模塊、hub基因、通路相關(guān)的情況。人類重要基因模塊在兩種相關(guān)性分 析中的情況則如下表：比較有意思的是模塊 5和模塊12。模塊5在DAVID數(shù)據(jù)庫(kù)的注釋中，有膜、跨膜、離子通道等 關(guān)鍵詞；而模塊12則和饑餓激素通路相關(guān)，在前人的研究 中，它可能會(huì)促進(jìn)胃腸道和胰腺疾病的惡化。由此，文章數(shù) 據(jù)分析的統(tǒng)計(jì)學(xué)結(jié)果終于與生物學(xué)意義產(chǎn)生了交匯！下圖展示了模塊5和EBV基因相互作用的關(guān)系圖，作為上 述相關(guān)性的一個(gè)示例：再來(lái)看看重要的人類hub基因，它們是 C1orf115, CNTD2 和

12、VANGL2 。它們?cè)?Pearson分析中， 和EBV的BALF1和LMP-1相關(guān)，在sCCA分析中也被選為 核心基因。下圖展示了人類hub基因和EBV基因的互作網(wǎng)絡(luò)。人類重要的通路則是凋亡、溶酶體、Jak-STAT信號(hào)通路和磷脂酰肌醇信號(hào)通路，它們?cè)赑earson和sCCA分析中都與EBV的BALF4和/或BALF5基因相關(guān)。值得注意的是， BALF4基因與Jak-STAT信號(hào)通路、磷脂酰肌醇信號(hào)通路中 的大多數(shù)基因都呈現(xiàn)正相關(guān)?？傮w而言，在第一步差異表達(dá)分析中所找到的DEx基因、模塊和通路，與篩選由來(lái)的 EBV基因，大多都有相關(guān)性，這也驗(yàn)證了前人的觀察和實(shí)驗(yàn)的結(jié)果，即EBV和胃癌確實(shí)存在

13、基因?qū)用娴南嗷プ饔谩Ｉ藕徒y(tǒng)計(jì)分析完成后，得到了一堆具有統(tǒng)計(jì)學(xué)相關(guān)性的基 因、通路，這時(shí)候需要我由它們的生物學(xué)相關(guān)性，才是拔高 文章水平的點(diǎn)睛之筆。作者分析了他們找到的每個(gè)具有統(tǒng)計(jì) 學(xué)相關(guān)性的元素在既往文獻(xiàn)報(bào)道中的相關(guān)性，文獻(xiàn)檢索的工 作相當(dāng)細(xì)致。例如LMP-1編碼的潛伏膜蛋白1可導(dǎo)致細(xì)胞生 長(zhǎng)失調(diào)，是一個(gè)致癌蛋白；BALF1編碼抗凋亡蛋白Bcl2同源體，它的表達(dá)也和幾種癌癥的發(fā)生有相關(guān)性，等等等等。當(dāng)然除了驗(yàn)證前人的研究，BALF2可謂一個(gè)新發(fā)現(xiàn)，之前沒(méi) 有報(bào)道過(guò)它和癌癥有什么關(guān)系，本文是第一次提由，如果能 夠經(jīng)過(guò)實(shí)驗(yàn)驗(yàn)證這文章絕非五分這般平庸啊。TCGA的挖掘潛力日益受到科研界的重視，搞腫

14、瘤的這一優(yōu) 勢(shì)得天獨(dú)厚，應(yīng)善加利用。生信的技能學(xué)習(xí)是一個(gè)由淺入深 的漸進(jìn)過(guò)程，解螺旋將推由系列課程，今天放由第一波，由 我們的合作伙伴GCBI制作的TCGA數(shù)據(jù)庫(kù)基礎(chǔ)應(yīng)用教學(xué)， 課件免費(fèi)下載，在服務(wù)號(hào)完成分享任務(wù)后獲取16mins視頻。懶癌晚期患者，請(qǐng)加入解螺旋鉆石會(huì)員，不但能夠?qū)W習(xí)臨床 研究，基金申請(qǐng)和 SCI寫作的深度課程，當(dāng)下還可享受“TCGA基因萬(wàn)能篩”活動(dòng)，1000元，給你用生信方式篩選 一個(gè)靠譜的研究靶標(biāo)基因，實(shí)驗(yàn)驗(yàn)證不成功免費(fèi)重篩，把專 業(yè)的事交給專業(yè)的人辦。下面是“TCGA基因萬(wàn)能篩”活動(dòng)獲得的靶標(biāo)基因示例：服務(wù)特色針對(duì)研究癌癥卻沒(méi)有方向的會(huì)員，整理TCGA數(shù)據(jù)庫(kù)中特定癌癥的癌

15、/癌旁數(shù)據(jù)，包括 miRNA或mRNA數(shù)據(jù)，篩選癌/ 癌旁差異基因，并結(jié)合多種生物信息學(xué)分析內(nèi)容(包括差異 基因功能、信號(hào)通路、生存以及基因網(wǎng)絡(luò)分析等)，為會(huì)員研究特定癌癥的發(fā)生發(fā)展機(jī)制提供有價(jià)值的靶基因及相應(yīng) 的理論基礎(chǔ)。分析結(jié)果展示 差異分析結(jié)果 差異基因功能分析 差異基因信號(hào)通路分析基于 KEGG數(shù)據(jù)庫(kù)，對(duì)差異基因利用 Fisher精確檢驗(yàn)和卡方檢驗(yàn)，把目標(biāo)基因參與的Pathway進(jìn)行顯著性分析，按照P value生存分析利用單基因在所有樣本 中的表達(dá)情況，按照中位數(shù)進(jìn)行分組，利用時(shí)序檢驗(yàn)(log-rank) 檢驗(yàn)進(jìn)行比較兩組的總生存期的差異性，篩選條件為：p調(diào)控關(guān)系每個(gè)基因與基因之間都有非常多的調(diào)控關(guān)系，我們從 基因相關(guān)的轉(zhuǎn)錄因子、miRNA、lncRNA和上下游相關(guān)基因 的四個(gè)角度來(lái)展示基因與基因間的關(guān)系。