张淑梅，张 彬，刘军厚，刘洪波，苏建忠，王 芳，张 岩

(哈尔滨医科大学生物信息科学与技术学院，黑龙江哈尔滨150081)

目前，癌症是严重威胁人类健康的三大杀手之一，对于这种严重危害人类健康的顽疾，现在的医学界并不十分清楚它的发病机制。同时，人们对基因的本质也渐渐有了更深入地认识。很长一段时间里，人们认为癌症的形成只与基因突变有关［1－4］。但是，越来越多的证据表明，表观遗传修饰对癌症的发生也起着十分重要的作用。基因序列不变，而基因的表型发生了可遗传的变化，称为表观遗传［5］。这是由表观遗传修饰造成的。DNA甲基化是一种重要的表观遗传修饰，它制约着基因的表达。

据报道，人类的基因只是果蝇的2倍多。很难想象DNA的遗传信息可以调控人类这样复杂的生命体发育和生存的全过程［1－4］。维持细胞的功能，决定哪些基因表达、哪些基因不表达，是非常重要的，几个基因的错误表达便会诱发正常细胞发生癌变［2］。

目前普遍认为，DNA甲基化与癌症的发生有密切关系［6］。癌症的甲基化异常表现为总体的甲基化水平降低与启动子区域的甲基化水平升高［7］。例如，抑癌基因与修复基因的高甲基化会导致它们的失活，造成肿瘤抑制丧失与基因损伤增加。

由于涉及基因的“开”与“关”，DNA甲基化对肿瘤的产生起着重大的作用。同时研究表明，某些基因的异常甲基化与多种癌症的产生有着显著的关联［8］。例如，基因P15的甲基化会使基因沉默，并使细胞过度激活与增殖，而这与白血病、淋巴瘤、鳞状细胞癌、肺癌的发生都有重要的联系［9］。是否存在一组甲基化异常的基因，与多种癌症的发生有着重要的关联以及这些基因在不同的癌症中是否起着不同的作用，成为本文关心的问题。通过研究这些问题，会为癌症的预测提供必要的方法，同时也增强了人们对癌症与DNA甲基化关系更进一步的认识。

表观遗传标记可以在被割除的肿瘤和体液中探测到。例如，超甲基化的癌症基因可以在尿斑中探测到，这在膀胱癌的检测中很有意义［10］。DNA甲基化的生物标记物在疾病诊断和预后的领域正在兴起，并且需要在临床实践中广为应用和扩展。

本课题首先通过对不同癌症DNA甲基化数据进行预处理，利用权重基因共表达网络分析方法(WGCNA)筛选出甲基化基因模块，并分析模块向量基因，利用 DAVID(The Database for Annotation，Visualization and Integrated Discovery)进行功能注释，然后对基因模块进行功能分析，得到DNA甲基化与肿瘤间的关系。本课题有助于发现癌症中DNA甲基化的生物标记物，为肿瘤的诊断及治疗提供可能的靶点。

1 数据及方法

1.1 数据获取

本课题所用的癌症甲基化数据来自于GEO数据库，包括乳腺导管癌甲基化数据(编号为GSE14865，平台为 GPL4126，6 个样本)［11］、胃癌甲基化数据(编号为GSE15291，平台为GPL4126，7个样本)、前列腺癌甲基化数据(编号为GSE15298，平台为 GPL4126，20 个样本)［12］、白血病甲基化数据(编号为GSE18400，平台为 GPL4126，60个样本，样本为婴儿期白血病数据和1个对照组)、食管鳞状细胞瘤甲基化数据(编号为 GSE21238，平台为GPL4126，6个样本，其中包括有淋巴结转移食道癌和无淋巴结转移食道癌样本以及转移淋巴结细胞样本)、肺鳞状细胞瘤甲基化数据(编号为GSE9622，平台为GPL4126，5个样本)。

首先进行数据的预处理和标准化，标准化的原则是对同一基因的不同探针对应的数值取平均值，并且只选择对应于启动子的探针。最后获得包含4 029个基因的甲基化数据。

1.2 权重基因共表达网络分析(WGCNA)的简介

网络分析在生物信息学中得到越来越多的应用。WGCNA(Weight Gene Co-express Network Analysis)是一种描述各个样本的基因芯片相关的系统生物学方法。这种方法可以找到高相关的基因模块，可以使用模块特征基因(eigengene)或hub节点间的基因彼此间和外部采样特征来聚类［13］。相关网络促进了基于基因筛选的方法的发展，可以用于识别候选生物标记物或治疗靶点。

1.3 甲基化基因模块筛选的原理

本文通过构建权重基因共表达网络来识别癌症中甲基化基因模块。

为了便于把显著差异的甲基化基因分类成模块，邻接矩阵被转换成拓扑重叠矩阵。拓扑重叠矩阵不仅可以捕捉到xi，xj的直接互作，也可以捕捉到间接互作。这样，定义了一个相似性测度:

其中，ki=代表点的连通性。1－TOMij是层次聚类的距离矩阵。

1.4 基因模块的功能分析

通过WGCNA筛选出甲基化基因模块并量化模块与表型的关系。分析与癌症表型显著相关的基因模块。挖掘出基因模块的向量基因，并利用DAVID生物信息学分析工具对基因模块进行GO功能注释与KEGG通路富集研究。

2 结果

2.1 选择合适的阈值:网络拓扑结构分析

构建一个权重基因网络，选择一个合适的邻接矩阵的阈值β，得到的阈值满足网络接近无于尺度的标准。通过WGCNA，选择一组候选的阈值，并返回被检测的网络参数(见图1)。从图中可看出阈值选择为5最合适，它既保证了网络接近于无尺度网络(模型指数大于0.9，完美无尺度网络的模型适应指数是1)，同时也是使曲线趋于平滑的最小阈值，并且它也使得网络的平均链接程度不会太小，这有利于网络包含足够的信息(例如，挖掘模块)。

图1 阈值分析Fig.1 Threshold value analysis

2.2 使用TOM(拓扑重叠矩阵)聚类

为了减少噪声和伪关联的影响，将邻接矩阵转换为拓扑重叠矩阵(TOM)。通过TOM，利用层次聚类产生一个基因的层次聚类树(见图2)。

在层次聚类树中，通过各个分支的识别(即从树图“剪枝”)得到模块。使用Dynamic Tree Cut的方法［14］，期望获得较大较少的模块，所以设定参数最小模块大小(minModuleSize)为50，这样从树图中剪枝共得到10个模块，标签为1至10，模块大小依次递减，从806至65个基因。模块0保存着所有模块外的基因。

树图中不同的深浅区域代表了不同的模块。找到匹配的模块，并返回各基因模块的宽度(见表1)。

注:树图中每个叶节点代表一个基因，其中密集连接的分支代表了甲基化数值接近的基因;图中不同的颜色代表不同的模块。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 量化模块与表型的关联

分析模块与模型的显著关联。由于已有每个模块的特征基因(eigengene)，使特征基因(eigengene)与表型相关联，并找到最大相关性。由于已有模块与表型，可以可视化这种关联，用颜色标注相关性。

图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，没有显著的相关性。由于选用的基因来自于启动子区域，可以得出结论:在上述基因模块中，白血病对应基因的甲基化程度要比正常细胞的甲基化程度高。

图3 模块与癌症表型关系图Fig.3 Module-trait relationship

食道鳞状细胞瘤(Esophageal Squamous Cell Carcinomas)在模块1、4、10都有显著的相关性。同时还发现，有淋巴结转移食道鳞状细胞瘤(ESCC with metastasis)和转移淋巴结(Metastatic lymph node)比无淋巴结转移食道鳞状细胞瘤在上面的模块中具有更高的相关性。可以得出结论:在上述模块中，有淋巴结转移食道鳞状细胞瘤比无淋巴结转移食道鳞状细胞瘤对应基因的甲基化程度高。

同时，也可看到胃癌与食管鳞状细胞瘤在模块1、4、10都有显著的相关性，说明在这几个基因模块中，两者甲基化程度较高。

2.4 基因模块的功能分析

为了进一步了解上面的基因模块与癌症发生与发展的关系，挖掘上面得到的模块的向量基因。并对这些向量基因进行基因本体功能分类及生物学通路分析。

首先，对与胃癌与食管鳞状细胞瘤显著相关的模块1、4、10的向量基因进行功能注释。这些模块中共得到1 148个基因。

通过DAVID分析，1 148个基因有647个注释到了189类生物学过程，其余为未知功能基因。设定阈值为p≤0.05，则基因注释到96类生物学过程。这些生物学过程主要包括:基因沉默，蛋白质降解过程，己糖降解，Wnt受体信号通路，蛋白激酶活性负调节等(见表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类生物学过程。这些生物学过程主要包括:蛋白激酶活性负调节，细胞增殖调节，调控细胞死亡，参与细胞形态分化等(见表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代谢过程;转录调控;磷酸化蛋白质氨基酸等(见表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表示一系列基因能代表与某些生物学通路相关的生物学功能发生了改变。本次分析中通路发生改变的主要有:产生癌症(Pathways in cancer)，产生肾上皮细胞癌(Renal cell carcinoma)(见表5)。

表5 基因模块1、4、10生物学通路中富集情况Table 5 Gene modules 1、4、10 biological pathway enrichment

通过以上的富集情况，发现基因模块1、4、10的相关基因富集到了产生癌症的通路。由于基因的启动子区域甲基化程度较高，会产生抑制表达的作用。基因表达的缺失导致低氧诱导因子(缺氧诱导因子－α)的积累，从而产生多种生长因子，如血管内皮生长因子和血小板衍生生长因子，使细胞运动、细胞转化、防止细胞凋亡等生物学效应的调节功能缺失，造成了肿瘤的生成。

同时，模块1、4、10的基因也富集到了产生肾上皮细胞癌的通路，这也说明了相关基因启动子区域的甲基化程度较高，影响到多种癌症的发生。

再对与前列腺癌显著相关的模块2的向量基因进行 KEGG通路分析。通过 DAVID分析，采用Fisher精确检验，本次分析中通路发生改变的是:细胞分裂周期(见表6)。

表6 基因模块2生物学通路中富集情况Table 6 Gene module 2 biological pathway enrichment

通过上面的富集，基因模块2的相关基因富集到影响细胞分裂周期的通路上。有丝分裂是一个重复序列的过程，细胞周期蛋白依赖性激酶(CDKs)是关键调节酶，它通过调控细胞基质来控制细胞进程。细胞周期蛋白依赖性激酶抑制剂(CKIs)，如基因CDC7、MAD1L1、CCNB3 等，参与 CDKs的负调控，从而提供了一个通过该细胞周期负调控的通路。而它又反过来激活p53抑癌蛋白。基因CDC7、MAD1L1、CCNB3的高甲基化，抑制p53的表达，同时细胞不能进行正常分裂，从而造成肿瘤细胞的产生。

3 讨论

目前的研究认为DNA甲基化与肿瘤密切相关。肿瘤的DNA甲基化改变表现为总体的甲基化水平降低与启动子区域CpG岛的甲基化水平升高。所筛选的基因模块的向量基因的甲基化水平普遍较高，就是由于基因的启动子区域CpG岛的甲基化异常造成的。

通过对基因模块进行GO功能注释，发现了与各癌症显著相关的甲基化异常的基因模块内的相应基因注释到了诸如基因沉默，Wnt受体信号通路;蛋白激酶活性负调节;细胞增殖调节;调控细胞死亡;参与细胞形态分化等生物过程，而这些生物学过程又与癌症的发生有着显著的关联。说明这些甲基化异常的基因模块对肿瘤的发生与发展起着重大的作用。

同时，对与胃癌与食管鳞状细胞瘤显著相关的模块1、4、10的向量基因进行生物学通路富集分析，得到产生癌症的通路。说明甲基化异常的基因模块确实与肿瘤的生成有着重要的联系。而对于与胃癌与食管鳞状细胞瘤显著相关的模块1、4、10富集到产生肾上皮细胞癌的通路。也说明了甲基化异常的基因模块同时与多种癌症的发生有着千丝万缕的联系。

在本课题中，首先下载了乳腺导管癌、胃癌、前列腺癌、白血病、肺鳞状细胞瘤、食管鳞状细胞瘤等6种癌症及亚型的DNA甲基化数据，经过预处理后利用WGCNA筛选出了甲基化基因模块，通过量化模块与癌症表型的关系发现了与各癌症显著相关的6个基因模块。然后，挖掘这些基因模块的向量基因，对这些基因进行GO功能注释和KEGG通路富集分析。通过GO功能注释发现了基因模块内相应的基因与可能导致肿瘤产生的生物学过程有关;利用KEGG数据库对基因模块的向量基因进行功能聚类，发现模块内的基因富集到产生癌症的通路也说明甲基化异常的基因模块与癌症的发生有着显著的内在关联。同时，也发现某些甲基化异常的基因模块(模块1、4、10)与多种癌症的发生有着显著的关联。基于此，本文有助于发现癌症中的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.