癌癥中DNA甲基化基因模塊篩選

張淑梅，張 彬，劉軍厚，劉洪波，蘇建忠，王 芳，張 巖

(哈爾濱醫(yī)科大學生物信息科學與技術學院，黑龍江哈爾濱150081)

目前，癌癥是嚴重威脅人類健康的三大殺手之一，對于這種嚴重危害人類健康的頑疾，現(xiàn)在的醫(yī)學界并不十分清楚它的發(fā)病機制。同時，人們對基因的本質(zhì)也漸漸有了更深入地認識。很長一段時間里，人們認為癌癥的形成只與基因突變有關［1－4］。但是，越來越多的證據(jù)表明，表觀遺傳修飾對癌癥的發(fā)生也起著十分重要的作用?；蛐蛄胁蛔?，而基因的表型發(fā)生了可遺傳的變化，稱為表觀遺傳［5］。這是由表觀遺傳修飾造成的。DNA甲基化是一種重要的表觀遺傳修飾，它制約著基因的表達。

據(jù)報道，人類的基因只是果蠅的2倍多。很難想象DNA的遺傳信息可以調(diào)控人類這樣復雜的生命體發(fā)育和生存的全過程［1－4］。維持細胞的功能，決定哪些基因表達、哪些基因不表達，是非常重要的，幾個基因的錯誤表達便會誘發(fā)正常細胞發(fā)生癌變［2］。

目前普遍認為，DNA甲基化與癌癥的發(fā)生有密切關系［6］。癌癥的甲基化異常表現(xiàn)為總體的甲基化水平降低與啟動子區(qū)域的甲基化水平升高［7］。例如，抑癌基因與修復基因的高甲基化會導致它們的失活，造成腫瘤抑制喪失與基因損傷增加。

由于涉及基因的“開”與“關”，DNA甲基化對腫瘤的產(chǎn)生起著重大的作用。同時研究表明，某些基因的異常甲基化與多種癌癥的產(chǎn)生有著顯著的關聯(lián)［8］。例如，基因P15的甲基化會使基因沉默，并使細胞過度激活與增殖，而這與白血病、淋巴瘤、鱗狀細胞癌、肺癌的發(fā)生都有重要的聯(lián)系［9］。是否存在一組甲基化異常的基因，與多種癌癥的發(fā)生有著重要的關聯(lián)以及這些基因在不同的癌癥中是否起著不同的作用，成為本文關心的問題。通過研究這些問題，會為癌癥的預測提供必要的方法，同時也增強了人們對癌癥與DNA甲基化關系更進一步的認識。

表觀遺傳標記可以在被割除的腫瘤和體液中探測到。例如，超甲基化的癌癥基因可以在尿斑中探測到，這在膀胱癌的檢測中很有意義［10］。DNA甲基化的生物標記物在疾病診斷和預后的領域正在興起，并且需要在臨床實踐中廣為應用和擴展。

本課題首先通過對不同癌癥DNA甲基化數(shù)據(jù)進行預處理，利用權重基因共表達網(wǎng)絡分析方法(WGCNA)篩選出甲基化基因模塊，并分析模塊向量基因，利用 DAVID(The Database for Annotation，Visualization and Integrated Discovery)進行功能注釋，然后對基因模塊進行功能分析，得到DNA甲基化與腫瘤間的關系。本課題有助于發(fā)現(xiàn)癌癥中DNA甲基化的生物標記物，為腫瘤的診斷及治療提供可能的靶點。

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

1.1 數(shù)據(jù)獲取

本課題所用的癌癥甲基化數(shù)據(jù)來自于GEO數(shù)據(jù)庫，包括乳腺導管癌甲基化數(shù)據(jù)(編號為GSE14865，平臺為 GPL4126，6 個樣本)［11］、胃癌甲基化數(shù)據(jù)(編號為GSE15291，平臺為GPL4126，7個樣本)、前列腺癌甲基化數(shù)據(jù)(編號為GSE15298，平臺為 GPL4126，20 個樣本)［12］、白血病甲基化數(shù)據(jù)(編號為GSE18400，平臺為 GPL4126，60個樣本，樣本為嬰兒期白血病數(shù)據(jù)和1個對照組)、食管鱗狀細胞瘤甲基化數(shù)據(jù)(編號為 GSE21238，平臺為GPL4126，6個樣本，其中包括有淋巴結轉移食道癌和無淋巴結轉移食道癌樣本以及轉移淋巴結細胞樣本)、肺鱗狀細胞瘤甲基化數(shù)據(jù)(編號為GSE9622，平臺為GPL4126，5個樣本)。

首先進行數(shù)據(jù)的預處理和標準化，標準化的原則是對同一基因的不同探針對應的數(shù)值取平均值，并且只選擇對應于啟動子的探針。最后獲得包含4 029個基因的甲基化數(shù)據(jù)。

1.2 權重基因共表達網(wǎng)絡分析(WGCNA)的簡介

網(wǎng)絡分析在生物信息學中得到越來越多的應用。WGCNA(Weight Gene Co-express Network Analysis)是一種描述各個樣本的基因芯片相關的系統(tǒng)生物學方法。這種方法可以找到高相關的基因模塊，可以使用模塊特征基因(eigengene)或hub節(jié)點間的基因彼此間和外部采樣特征來聚類［13］。相關網(wǎng)絡促進了基于基因篩選的方法的發(fā)展，可以用于識別候選生物標記物或治療靶點。

1.3 甲基化基因模塊篩選的原理

本文通過構建權重基因共表達網(wǎng)絡來識別癌癥中甲基化基因模塊。

為了便于把顯著差異的甲基化基因分類成模塊，鄰接矩陣被轉換成拓撲重疊矩陣。拓撲重疊矩陣不僅可以捕捉到xi，xj的直接互作，也可以捕捉到間接互作。這樣，定義了一個相似性測度:

其中，ki=代表點的連通性。1－TOMij是層次聚類的距離矩陣。

1.4 基因模塊的功能分析

通過WGCNA篩選出甲基化基因模塊并量化模塊與表型的關系。分析與癌癥表型顯著相關的基因模塊。挖掘出基因模塊的向量基因，并利用DAVID生物信息學分析工具對基因模塊進行GO功能注釋與KEGG通路富集研究。

2 結果

2.1 選擇合適的閾值:網(wǎng)絡拓撲結構分析

構建一個權重基因網(wǎng)絡，選擇一個合適的鄰接矩陣的閾值β，得到的閾值滿足網(wǎng)絡接近無于尺度的標準。通過WGCNA，選擇一組候選的閾值，并返回被檢測的網(wǎng)絡參數(shù)(見圖1)。從圖中可看出閾值選擇為5最合適，它既保證了網(wǎng)絡接近于無尺度網(wǎng)絡(模型指數(shù)大于0.9，完美無尺度網(wǎng)絡的模型適應指數(shù)是1)，同時也是使曲線趨于平滑的最小閾值，并且它也使得網(wǎng)絡的平均鏈接程度不會太小，這有利于網(wǎng)絡包含足夠的信息(例如，挖掘模塊)。

圖1 閾值分析Fig.1 Threshold value analysis

2.2 使用TOM(拓撲重疊矩陣)聚類

為了減少噪聲和偽關聯(lián)的影響，將鄰接矩陣轉換為拓撲重疊矩陣(TOM)。通過TOM，利用層次聚類產(chǎn)生一個基因的層次聚類樹(見圖2)。

在層次聚類樹中，通過各個分支的識別(即從樹圖“剪枝”)得到模塊。使用Dynamic Tree Cut的方法［14］，期望獲得較大較少的模塊，所以設定參數(shù)最小模塊大小(minModuleSize)為50，這樣從樹圖中剪枝共得到10個模塊，標簽為1至10，模塊大小依次遞減，從806至65個基因。模塊0保存著所有模塊外的基因。

樹圖中不同的深淺區(qū)域代表了不同的模塊。找到匹配的模塊，并返回各基因模塊的寬度(見表1)。

注:樹圖中每個葉節(jié)點代表一個基因，其中密集連接的分支代表了甲基化數(shù)值接近的基因;圖中不同的顏色代表不同的模塊。Notes:In the tree diagram，every leaf node represents a gene，and branches densely connected represent the genes which have the similar methylation values.The different colors represent different modules.

表1 各基因模塊中的向量基因Table1 Vector genes in the modules

2.3 量化模塊與表型的關聯(lián)

分析模塊與模型的顯著關聯(lián)。由于已有每個模塊的特征基因(eigengene)，使特征基因(eigengene)與表型相關聯(lián)，并找到最大相關性。由于已有模塊與表型，可以可視化這種關聯(lián)，用顏色標注相關性。

圖3中可以清晰地看到模塊與癌癥表型的相關性。模塊0的基因是樹圖中剔除的基因，從圖3中也可看出它與各癌癥表型的相關性較差，因此不予考慮。

可以看到，胃癌(gastric cancer)與模塊1、4、10有較強的相關性(p≤0.05);前列腺癌(prostate cancer)與模塊2有較強的相關性(p≤0.05)。說明這兩種癌癥在上述模塊中甲基化程度較高。

同時，可以看到各種嬰兒期白血病(ALL)，如MLL－AF4白血病、MLL－ENL白血病、未擴散白血病(Untranslocated infant ALL)在模塊7、9都有較強的相關性(p≤0.05)，而正常人體細胞在這兩個模塊中p值都大于0.05，沒有顯著的相關性。由于選用的基因來自于啟動子區(qū)域，可以得出結論:在上述基因模塊中，白血病對應基因的甲基化程度要比正常細胞的甲基化程度高。

圖3 模塊與癌癥表型關系圖Fig.3 Module-trait relationship

食道鱗狀細胞瘤(Esophageal Squamous Cell Carcinomas)在模塊1、4、10都有顯著的相關性。同時還發(fā)現(xiàn)，有淋巴結轉移食道鱗狀細胞瘤(ESCC with metastasis)和轉移淋巴結(Metastatic lymph node)比無淋巴結轉移食道鱗狀細胞瘤在上面的模塊中具有更高的相關性?？梢缘贸鼋Y論:在上述模塊中，有淋巴結轉移食道鱗狀細胞瘤比無淋巴結轉移食道鱗狀細胞瘤對應基因的甲基化程度高。

同時，也可看到胃癌與食管鱗狀細胞瘤在模塊1、4、10都有顯著的相關性，說明在這幾個基因模塊中，兩者甲基化程度較高。

2.4 基因模塊的功能分析

為了進一步了解上面的基因模塊與癌癥發(fā)生與發(fā)展的關系，挖掘上面得到的模塊的向量基因。并對這些向量基因進行基因本體功能分類及生物學通路分析。

首先，對與胃癌與食管鱗狀細胞瘤顯著相關的模塊1、4、10的向量基因進行功能注釋。這些模塊中共得到1 148個基因。

通過DAVID分析，1 148個基因有647個注釋到了189類生物學過程，其余為未知功能基因。設定閾值為p≤0.05，則基因注釋到96類生物學過程。這些生物學過程主要包括:基因沉默，蛋白質(zhì)降解過程，己糖降解，Wnt受體信號通路，蛋白激酶活性負調(diào)節(jié)等(見表2)。

表2 與胃癌和食管鱗狀細胞瘤顯著相關的模塊向量基因的功能富集聚簇Table 2 Functional annotation for module vector genes significantly associated with gastric cancer and ESCC

同理，對與前列腺癌顯著相關的模塊2的向量基因進行功能注釋。注釋的372個基因有204個注釋到了127類生物學過程。設定閾值為p≤0.05，則基因注釋到79類生物學過程。這些生物學過程主要包括:蛋白激酶活性負調(diào)節(jié)，細胞增殖調(diào)節(jié)，調(diào)控細胞死亡，參與細胞形態(tài)分化等(見表3)。

表3 與前列腺癌顯著相關的模塊向量基因的功能富集聚簇Table 3 Functional annotation for module vector genes significantly associated with prostate cancer

同理，對與白血病顯著相關的模塊7、9的向量基因進行功能注釋。注釋的148個基因有141個注釋到了27類生物學過程。設定閾值為p≤0.05，則基因注釋到14類生物學過程。這些生物學過程主要包括:磷代謝過程;mRNA代謝過程;轉錄調(diào)控;磷酸化蛋白質(zhì)氨基酸等(見表4)。

表4 與白血病顯著相關的模塊向量基因的功能富集聚簇Table 4 Functional annotation for module vector genes significantly associated with leukemia

接著，對與胃癌與食管鱗狀細胞瘤顯著相關的模塊1、4、10的向量基因進行KEGG通路分析。這些模塊中共得到1 148個基因。

通過DAVID分析，對向量基因進行生物學通路富集分析。采用Fisher精確檢驗，p＜0.05表示一系列基因能代表與某些生物學通路相關的生物學功能發(fā)生了改變。本次分析中通路發(fā)生改變的主要有:產(chǎn)生癌癥(Pathways in cancer)，產(chǎn)生腎上皮細胞癌(Renal cell carcinoma)(見表5)。

表5 基因模塊1、4、10生物學通路中富集情況Table 5 Gene modules 1、4、10 biological pathway enrichment

通過以上的富集情況，發(fā)現(xiàn)基因模塊1、4、10的相關基因富集到了產(chǎn)生癌癥的通路。由于基因的啟動子區(qū)域甲基化程度較高，會產(chǎn)生抑制表達的作用?；虮磉_的缺失導致低氧誘導因子(缺氧誘導因子－α)的積累，從而產(chǎn)生多種生長因子，如血管內(nèi)皮生長因子和血小板衍生生長因子，使細胞運動、細胞轉化、防止細胞凋亡等生物學效應的調(diào)節(jié)功能缺失，造成了腫瘤的生成。

同時，模塊1、4、10的基因也富集到了產(chǎn)生腎上皮細胞癌的通路，這也說明了相關基因啟動子區(qū)域的甲基化程度較高，影響到多種癌癥的發(fā)生。

再對與前列腺癌顯著相關的模塊2的向量基因進行 KEGG通路分析。通過 DAVID分析，采用Fisher精確檢驗，本次分析中通路發(fā)生改變的是:細胞分裂周期(見表6)。

表6 基因模塊2生物學通路中富集情況Table 6 Gene module 2 biological pathway enrichment

通過上面的富集，基因模塊2的相關基因富集到影響細胞分裂周期的通路上。有絲分裂是一個重復序列的過程，細胞周期蛋白依賴性激酶(CDKs)是關鍵調(diào)節(jié)酶，它通過調(diào)控細胞基質(zhì)來控制細胞進程。細胞周期蛋白依賴性激酶抑制劑(CKIs)，如基因CDC7、MAD1L1、CCNB3 等，參與 CDKs的負調(diào)控，從而提供了一個通過該細胞周期負調(diào)控的通路。而它又反過來激活p53抑癌蛋白?；駽DC7、MAD1L1、CCNB3的高甲基化，抑制p53的表達，同時細胞不能進行正常分裂，從而造成腫瘤細胞的產(chǎn)生。

3 討論

目前的研究認為DNA甲基化與腫瘤密切相關。腫瘤的DNA甲基化改變表現(xiàn)為總體的甲基化水平降低與啟動子區(qū)域CpG島的甲基化水平升高。所篩選的基因模塊的向量基因的甲基化水平普遍較高，就是由于基因的啟動子區(qū)域CpG島的甲基化異常造成的。

通過對基因模塊進行GO功能注釋，發(fā)現(xiàn)了與各癌癥顯著相關的甲基化異常的基因模塊內(nèi)的相應基因注釋到了諸如基因沉默，Wnt受體信號通路;蛋白激酶活性負調(diào)節(jié);細胞增殖調(diào)節(jié);調(diào)控細胞死亡;參與細胞形態(tài)分化等生物過程，而這些生物學過程又與癌癥的發(fā)生有著顯著的關聯(lián)。說明這些甲基化異常的基因模塊對腫瘤的發(fā)生與發(fā)展起著重大的作用。

同時，對與胃癌與食管鱗狀細胞瘤顯著相關的模塊1、4、10的向量基因進行生物學通路富集分析，得到產(chǎn)生癌癥的通路。說明甲基化異常的基因模塊確實與腫瘤的生成有著重要的聯(lián)系。而對于與胃癌與食管鱗狀細胞瘤顯著相關的模塊1、4、10富集到產(chǎn)生腎上皮細胞癌的通路。也說明了甲基化異常的基因模塊同時與多種癌癥的發(fā)生有著千絲萬縷的聯(lián)系。

在本課題中，首先下載了乳腺導管癌、胃癌、前列腺癌、白血病、肺鱗狀細胞瘤、食管鱗狀細胞瘤等6種癌癥及亞型的DNA甲基化數(shù)據(jù)，經(jīng)過預處理后利用WGCNA篩選出了甲基化基因模塊，通過量化模塊與癌癥表型的關系發(fā)現(xiàn)了與各癌癥顯著相關的6個基因模塊。然后，挖掘這些基因模塊的向量基因，對這些基因進行GO功能注釋和KEGG通路富集分析。通過GO功能注釋發(fā)現(xiàn)了基因模塊內(nèi)相應的基因與可能導致腫瘤產(chǎn)生的生物學過程有關;利用KEGG數(shù)據(jù)庫對基因模塊的向量基因進行功能聚類，發(fā)現(xiàn)模塊內(nèi)的基因富集到產(chǎn)生癌癥的通路也說明甲基化異常的基因模塊與癌癥的發(fā)生有著顯著的內(nèi)在關聯(lián)。同時，也發(fā)現(xiàn)某些甲基化異常的基因模塊(模塊1、4、10)與多種癌癥的發(fā)生有著顯著的關聯(lián)?；诖耍疚挠兄诎l(fā)現(xiàn)癌癥中的DNA甲基化的生物標記物，為腫瘤的診斷及治療提供可能的靶點。

References)

［1］ WU C，MORRIS J R.Genes，genetics，and epigenetics:A correspondence［J］.Science，2001，293(5532):1103 －5.

［2］ WOLFFE A P.Chromatin remodeling:Why it is important in cancer［J］.Oncogene，2001，20(24):2988 －90.

［3］ PENNISI E.Behind the scenes of gene expression［J］.Science，2001，293(5532):1064－7.

［4］ VALLBOHMER D，BRABENDER J，YANG D，et al.DNA methyltransferases messenger RNA expression and aberrant methylation of CpG islands in non-small-cell lung cancer:association and prognostic value［J］.Clinical Lung Cancer，2006，8(1):39 －44.

［5］ ALVAREZ-VENEGAS R，AVRAMOVA Z.Methylation patterns of histone H3 Lys 4，Lys 9 and Lys 27 in transcriptionally active and inactive Arabidopsis genes and in atx1 mutants［J］.Nucleic Acids Res，2005，33(16):5199－207.

［6］ AHMAD I，RAO，D.N.Chemistry and biology of DNA methyltransferases［J］.Critical Reviews in Biochemistry and Molecular Biology，1996，31(5－6):361－380.

［7］ VERTINO P M，YEN R W，GAO J，et al.De novo methylation of CpG island sequences in human fibroblasts overexpressing DNA(cytosine-5-)-methyltransferase［J］.Mol Cell Biol，1996，16(8):4555 －65.

［8］ AHUJA N，LI Q，MOHAN A L，et al.Aging and DNA methylation in colorectal mucosa and cancer［J］.Cancer Res，1998，58(23):5489 －94.

［9］ WHEELER J M，BECK N E， KIM H C， et al.Mechanisms of inactivation of mismatch repair genes in human colorectal cancer cell lines:the predominant role of hMLH1［J］.Proc Natl Acad Sci U S A，1999，96(18):10296－301.

［10］ BAYLIN S B.Tying it all together:epigenetics，genetics，cell cycle，and cancer［J］.Science，1997，277(5334):1948－9.

［11］ MURAKAMI J，ASAUMI J，MAKI Y，et al.Influence of CpG island methylation status in O6-methylguanine-DNA methyltransferase expression of oral cancer cell lines［J］.Oncol Rep，2004，12(2):339－45.

［12］NEPHEW K P，HUANG T H.Epigenetic gene silencing in cancer initiation and progression ［J］.Cancer Lett，2003，190(2):125－33.

［13］ ZHANG B，HORVATH S.A general framework for weighted gene co-expression network analysis［J］.Stat Appl Genet Mol Biol，2005，4(1):1 －43.

［14］ LANGFELDER P，ZHANG B，HORVATH S.Defining clusters from a hierarchical cluster tree:the Dynamic Tree Cut package for R ［J］.Bioinformatics，2008，24(5):719－20.