嵇江淮 趙瀟瀟李乾鵬 安 奕 趙 磊 李冬果*

(1.首都醫(yī)科大學生物醫(yī)學工程學院生物醫(yī)學信息學系，北京 100069；2.首都醫(yī)科大學宣武醫(yī)院麻醉手術(shù)科，北京 100053；3. 國家老年疾病臨床研究中心，北京 100053)

膠質(zhì)母細胞瘤(glioblastoma, GBM)是人類最常見且致死率極高的一種中樞神經(jīng)系統(tǒng)腫瘤，發(fā)生率約占膠質(zhì)瘤的69%[1]。這種腦瘤具有高浸潤性，預后差，患者的中位生存期大約只有一年[2]。近年來，隨著大量的蛋白編碼基因(protein-coding genes, PCGs)被發(fā)現(xiàn)，編碼基因表達的失調(diào)被證明與疾病的進程有密切關(guān)系[3-4]。有研究[5]表明GBM中EGFR基因表達顯著升高，促進GBM患者細胞的凋亡。Meng等[6]發(fā)現(xiàn)TCTN1基因的過表達促進GBM細胞的增生，并且TCTN1的過表達可以作為預測GBM患者的獨立預后因素。這些研究表明PCGs在GBM進程中發(fā)揮非常重要的作用，但是這些研究主要關(guān)注PCGs的表達模式。目前，GBM中絕大部分PCGs的調(diào)控機制(特別是DNA甲基化調(diào)控機制)尚不清楚。

DNA甲基化是表觀遺傳調(diào)控的一種重要模式。在基因啟動子處的DNA甲基化對基因的表達具有重要的作用，并且參與到許多疾病的進程中[7]。一般來說，基因啟動子處的高甲基化往往抑制轉(zhuǎn)錄因子的組合并且下調(diào)甚至沉默癌癥抑制基因，基因啟動子的低甲基化往往激活致癌基因[8]。研究[9]表明GBM中MGMT基因啟動子甲基化和MGMT的表達具有很強的相關(guān)性。然而，這些研究僅僅分析一些特殊的PCGs的DNA甲基化模式，并沒有系統(tǒng)分析腫瘤中PCGs的DNA甲基化的模式以及全面評估腫瘤中DNA甲基化對基因表達的調(diào)控關(guān)系。

隨著高通量測序技術(shù)的發(fā)展，RNA-seq和Infnium 450k等眾多數(shù)據(jù)已經(jīng)應用于癌癥分析。本研究系統(tǒng)整合多組學數(shù)據(jù)，通過重新注釋DNA甲基化陣列，系統(tǒng)分析GBM進程中PCGs的DNA甲基化調(diào)控模式。本研究將會幫助理解GBM中DNA甲基化的調(diào)控機制，為識別GBM生物標志物和潛在的治療靶點提供新的見解。

1 數(shù)據(jù)來源及方法

1.1 數(shù)據(jù)的來源，整合和標準化

GBM的DNA甲基化數(shù)據(jù)來自Infinium HM450k平臺[10]，本研究從TCGA中下載GBM腫瘤樣本的HM450k數(shù)據(jù)，從GEO(GSE41826[11-12],其中包含58個正常膠質(zhì)細胞樣本)中下載正常樣本的HM450k數(shù)據(jù)。正常數(shù)據(jù)是通過同一平臺用同樣方法獲得的。GBM表達譜數(shù)據(jù)來自于TCGA數(shù)據(jù)庫，包括5個正常樣本數(shù)據(jù)和152個臨床樣本信息。為了保證數(shù)據(jù)的質(zhì)量，本研究挑選50個同時具有DNA甲基化數(shù)據(jù)和表達譜數(shù)據(jù)的樣本用作進一步的分析。人類全基因組的注釋數(shù)據(jù)(V19)收集于GENCODE數(shù)據(jù)庫[13]。

對下載的RNA-seq(counts)數(shù)據(jù)進行l(wèi)og2轉(zhuǎn)換并使用R包“edgeR”[14]對數(shù)據(jù)進行標準化，最終獲得PCGs的表達譜數(shù)據(jù)。由于HM450k數(shù)據(jù)存在一定的缺失值，所以，在構(gòu)建甲基化譜之前需要對含有缺失的數(shù)據(jù)進行處理。K最近鄰填充算法 (K-nearest neighbor, KNN) 是用數(shù)據(jù)集中缺失數(shù)據(jù)的k個最近鄰來估計缺失值[15]。本研究計算在所有腫瘤樣本中具有缺失值的探針數(shù)量，并用“DMwR”包中的knnImputation函數(shù)補全剩余的缺失值，以評估甲基化探針的質(zhì)量。最終有89 512個探針被移除，獲得392 867個甲基化位點探針。

1.2 PCGs甲基化譜的構(gòu)建

本研究用映射到PCGs啟動子區(qū)域探針的甲基化水平來構(gòu)建PCGs的甲基化譜。采用Zhi等[16]的重注釋方法將392 867個探針映射到PCGs的啟動子區(qū)域(轉(zhuǎn)錄開始位點到上游10 kb區(qū)域內(nèi))，使用最接近每個轉(zhuǎn)錄開始位點的探針來確定PCGs啟動子的DNA甲基化水平[17]，從而構(gòu)建GBM相關(guān)的PCGs甲基化譜。

1.3 不同甲基化模式PCGs功能的預測

本研究使用基于線性模型設計的R包“l(fā)imma”[18]計算腫瘤和正常樣本之間的差異甲基化和差異表達。通過Benjamini-Hochberg方法校正P值。為了全面分析表達受對應啟動子異常甲基化調(diào)控的PCGs，分別將高甲基化的PCGs和低表達的PCGs及低甲基化的PCGs和高表達的PCGs取交集，得到兩種具有不同生物學意義的情況：高甲基化且低表達的PCGs及低甲基化且高表達的PCGs。

為了預測不同甲基化模式下PCGs的功能，采用富集分析方法，對挑選出的PCGs進行功能和通路的顯著性分析，使用“clusterProfiler”[19]包來預測不同甲基化模式下PCGs的功能。通過Benjamini-Hochberg方法校正P值，如果矯正后的P值≤0.05，該GO項和富集通路就認為是顯著的。

1.4 生存分析

為了鑒別挑選出的PCGs是否具有良好的預后效果，本研究基于152個GBM患者的表達譜信息，依據(jù)中值將患者分為兩組。Kaplan-Meier生存分析和log-rank檢測被用來評估兩組患者的生存差異。

所有的研究均使用R3.5.1完成。

2 結(jié)果

2.1 重注釋Infinium HM 450K陣列并構(gòu)建PCGs的DNA甲基化譜

為了描繪PCGs的DNA甲基化調(diào)控模式，本研究將甲基化數(shù)據(jù)重新注釋到人類PCGs相關(guān)的啟動子區(qū)域，共有125 442個探針落在14 684個PCGs啟動子區(qū)域，通常每一個PCG都有幾個探針落在其啟動子區(qū)域，在此僅保留最接近每個轉(zhuǎn)錄開始位點的探針來確定PCGs啟動子的DNA甲基化狀態(tài)。

2.2 整合差異甲基化和差異表達的PCGs

為了能夠有效地識別GBM相關(guān)的風險標志物，研究差異甲基化可能的生物學意義，基于鑒別出的3 561個差異甲基化的PCGs，其中高甲基化的PCGs有873個(24.5%)，低甲基化的PCGs有2 688個(75.5%)。顯然低甲基化的PCGs的數(shù)量遠遠多于高甲基化的PCGs的數(shù)量。這些結(jié)果顯示，在GBM的發(fā)生發(fā)展中PCGs呈現(xiàn)為更多的低甲基化模式，并且這種全局的低甲基化可能導致致癌基因的激活和影響基因組的穩(wěn)定性。進一步地，本研究鑒別出 6 586 個差異表達的PCGs，其中有2 788個PCGs表達上調(diào)，3 798個PCGs表達下調(diào)。

為了分析異常甲基化對PCGs表達的調(diào)控機制，針對兩種不同生物學意義的情況，識別出240個高甲基化且低表達的PCGs及390個低甲基化且高表達的PCGs。結(jié)果表明，PCGs在腫瘤樣本和正常樣本中表現(xiàn)出差異甲基化和差異表達兩種模式。

2.3 不同甲基化模式下PCGs的富集分析

對不同甲基化模式下的PCGs，采用富集分析方法進行功能和通路的顯著性分析。結(jié)果顯示，高甲基化的PCGs主要富集在神經(jīng)元系統(tǒng)發(fā)展、細胞-細胞信號傳導等生物過程;低甲基化的PCGs富集在許多與腫瘤進程相關(guān)的生物過程，比如細胞黏附、細胞遷移、免疫相關(guān)的細胞增生和血管生成(圖1A、1B)。因此，低甲基化的PCGs可能是影響GBM進程的一個重要模塊。對于KEGG通路富集分析，高甲基化的PCGs主要富集在GABAergic突觸、胰島素分泌等通路中，而低甲基化的PCGs主要富集在金黃色葡萄球菌感染、補體和凝血級聯(lián)等通路中(圖1C、1D)。

圖1 不同甲基化模式下PCGs的富集分析Fig.1 Enrichment analysis of PCGs under different methylation patternsA and B are the GO function analysis of hypermethylated-low expressed PCGs and hypomethylated-highly expressed PCGs,respectively. C and D are KEGG pathway analysis of hypermethylated-low expressed PCGs and hypomethylated-highly expressed PCGs,respectively.The depth of color represents the number of PCGs enriched in the GO item or pathway; PCGs: protein-coding gens.

2.4 GBM中預后良好的PCGs的鑒別

為了評估這些表達受對應啟動子異常甲基化影響的PCGs是否可以作為GBM的預后因子，本研究將表達譜數(shù)據(jù)與臨床信息相結(jié)合，進行基因生存分析(圖2)，基于生存分析的結(jié)果，獲得了51個與GBM整體存活顯著相關(guān)的PCGs(P<0.05)，如表1所示。

圖2 SOCS1和AEBP1基因生存分析Fig.2 Survival analysis of SOCS1 and AEBP1 genesA and B represent Kaplan-Meier curves of SOCS1 and AEBP1 genes in GBM, respectively.The red and blue lines at the bottom represent the number of patients who survived.

表1 預后良好的PCGsTab.1 PCGs with good prognosis

3 討論

DNA甲基化的變化可能破壞癌癥中特定啟動子的調(diào)控[20]。近年來，隨著高通量測序技術(shù)的發(fā)展，表觀遺傳調(diào)控已經(jīng)成為生物醫(yī)學研究的一大熱點。一項研究[21]表明GBM中B3GNT5、FABP7等13個基因啟動子的甲基化和表達水平之間呈現(xiàn)負相關(guān)。

本研究對重注釋得到的PCGs甲基化譜進行差異甲基化分析，結(jié)果顯示,GBM中多呈現(xiàn)一種低甲基化的模式。為了驗證這個結(jié)果，本研究重新注釋了GBM中l(wèi)ncRNA的甲基化譜，結(jié)果發(fā)現(xiàn)在鑒別出的5 567個異常甲基化的lncRNAs中高甲基化的有1 214個(21.8%)，低甲基化的有4 353個(78.2%)。為了評估是否不同甲基化模式的PCGs對應著不同的生物學功能，對這些基因進行功能富集分析。結(jié)果顯示低甲基化的PCGs參與了許多與腫瘤發(fā)生及進程相關(guān)的生物過程。由此推斷這些低甲基化的PCGs更有可能成為GBM診斷和治療的生物標志物。有研究[22-23]表明AEBP1和SOCS1在GBM中過表達,沉默他們的表達可以抑制GBM細胞的增生。在本研究中，AEBP1和SOCS1均屬于低甲基化且高表達的基因，進一步地，AEBP1和SOCS1基因表達水平較低的GBM患者具有更加良好的預后。結(jié)果表明AEBP1和SOCS1在GBM中表達的上調(diào)可能是受它們基因啟動子的低甲基化調(diào)控，而抑制它們的表達可以提高GBM患者的生存概率。

當然，由于GBM不同于一般的腫瘤，樣本的獲取是一個重大的難題。本研究需要進一步搜集大規(guī)模樣本進行重復性計算或者動物實驗來驗證挖掘出的表達受異常甲基化調(diào)控的PCGs以及潛在的治療靶點。

綜上，本研究通過重新注釋DNA甲基化陣列，系統(tǒng)識別GBM中潛在的表達受異常甲基化調(diào)控的PCGs，加深對GBM中PCGs甲基化調(diào)控模式的理解，并對識別GBM風險標志物和潛在的治療靶點提出了新的認識。