韦云真，刘晓娟，王 芳，苏建忠，张岩∗，刘洪波∗

（1.哈尔滨医科大学生物信息科学与技术学院，哈尔滨150081；2.哈尔滨医科大学附属第一医院康复医学科，哈尔滨150001）

doi：10.3969/j.issn.1672－5565.2015.03.05

癌症DNA甲基化调控位点的识别

韦云真1，刘晓娟2，王 芳1，苏建忠1，张岩1∗，刘洪波1∗

（1.哈尔滨医科大学生物信息科学与技术学院，哈尔滨150081；
2.哈尔滨医科大学附属第一医院康复医学科，哈尔滨150001）

DNA甲基化是一种重要的表观遗传学修饰，在基因的转录调控方面具有重要的作用。异常的DNA甲基化可以导致癌症等复杂疾病发生，癌基因相关的DNA甲基化调控位点的识别对于解析癌症的发生发展机制及识别新的癌症标记具有重要意义。本研究通过整合The Cancer Genome Atlas（TCGA）的泛癌症基因组的高通量甲基化谱和基因表达谱，识别癌基因相关的DNA甲基化调控位点。对于每种癌症分批次计算CpG位点甲基化与相关基因表达之间的相关性，并筛选调控下游基因的CpG位点（包括强调控位点、弱调控位点和不调控位点），结果表明仅有一半的CpG位点对下游基因具有调控作用；对癌症间共享的调控位点的分析发现不同癌症间共享的调控位点不尽相同，表明癌症特异的甲基化调控位点的存在。进一步地，对差异甲基化和差异表达基因的功能富集分析揭示了受甲基化调控的基因确实参与了癌症发生发展相关的功能。本研究的结果是对当前甲基化调控位点集的重要补充，也是识别癌症新型分子标记特征的重要资源。

DNA甲基化；基因表达；转录调控；癌症

DNA甲基化是一种重要的表观遗传学修饰，在CpG岛（DNA的CG序列密集区）上发生，对调控转录基因具有重要的作用［1－3］。甲基化位点可随DNA的复制而遗传，因为DNA复制后，甲基化酶可将新合成的未甲基化的位点进行甲基化［4］。DNA的甲基化可引起基因的失活。如CpG岛位于某基因的启动子区域，CpG岛的甲基化会显著降低甚至完全沉默该基因的转录，继而影响蛋白的表达。CpG岛的甲基化程度越高，基因表达的程度越低。目前对于甲基化调控基因表达，以及进一步的生物学影响的研究很多［5］。然而对具体的基因，尚没有一个完整的甲基化调控区域的图谱。

最近几年来，表观遗传学领域发展十分迅速。DNA甲基化修饰就是一个非常重要的部分，参与基因表达调控、转座子沉默、X染色体失活、基因印记、以及癌症发生等重要生物学过程［6－8］。近年来随着研究技术和方法的进步，全基因组DNA甲基化的研究广泛兴起，很多物种的全基因组甲基化图谱破译了出来，DNA甲基化全局水平的研究不仅有利于宏观层面上了解DNA甲基化的规律和特性，同时也为深入分析DNA甲基化生物学调控及功能奠定了基础。如今，在当前领域已经取得了一些进展，例如发现了DNA甲基化酶的DNMT家族，并且对其作用机制和生理功能进行了一些研究［9－11］。甲基化与癌症的发生有关系［12－14］。研究发现抑癌基因启动子高度甲基化后，可以令这些基因表达受到抑制，同癌症的发生有着十分密切的关系，从而研究DNA甲基化抑制剂的使用，将有助于预防人类肿瘤的发生。

基于高通量的DNA甲基化数据［15］，研究DNA甲基化与基因表达之间的关系，并建立生物学测度筛选对基因表达之间的关系，筛选对基因具有调控作用的DNA甲基化区域，最后绘制基因组范围内参与基因调控的DNA甲基化区域。研究成果将有助于对表观遗传调控机制更深入理解。

1 材料与方法

1.1 材料

研究使用The Cancer Genome Atlas（TCGA）上同时具有甲基化数据和表达数据的所有27K高通量癌症数据［16］。如表1所示，根据以上挑选条件，总计11癌症符合条件，15套数据，81个处理批次（Batch）。不同癌症内包含着不同的样本数与处理批次。其中BRCA，OV与READ癌症数据同时包含有癌症样本与正常样本，如表1所示。基因表达数据与甲基化数据样本数数量一致，并且是一一对应的，基因表达部分使用level2数据，level2数据内容为探针名－取log2后的表达值。甲基化部分使用的是level3数据，level3数据内容为甲基化位点名－甲基化值。

表1 研究采用的数据Table1 Data of the research

1.2 方 法

1.2.1 数据预处理

将Methylation、Expression分批次进行数据的标准化。对于每一个甲基化位点，我们按不同的Batch计算这个甲基化位点一一对应样本的Methylation部分及Expression部分数据的P值及皮尔森相关系数（PCC），进行DNA甲基化－基因调控关系定量。合并各Batch的结果并进行进一步的分析。

1.2.2 癌症甲基化与表达值的PCC分布曲线绘制

保留皮尔森相关系数显著（P＜0.05）的CpG位点－基因对，并利用R语言里的ggplot2包进行分布曲线的绘制。

为了验证在P值取不同阈值的情况下，PCC值的分布是否是稳定的，取不同的显著性P值下进行重复研究，显著性阈值分别取P＜1.0×10－2，P＜1.0× 10－3，P＜1.0×10－4，P＜1.0×10－5，P＜1.0×10－6，P＜1.0×10－7。并分别绘制PCC值分布曲线图。

1.2.3 Batch间相关性分析

使用预处理数据研究癌症Batch间的相关性，我们进行了Batch间的相关性分析，对于每一个Batch，我们认为P小于0.05时的PCC值是有用的，进行保留，而P值大于等于0.05时的PCC值认为是没有用的，将这时的PCC值更改为0。接着，将所有的P值列删除，只保存PCC值。计算两两Batch间的皮尔森相关系数。用cluster进行双向聚类。并将聚类结果用TreeView进行可视化。

1.2.4 甲基化位点与Batch相关性分析

对于每一个Batch，当P值小于0.05时，PCC值有效，保留原值，当P值大于等于0.05时，认为是不显著的，PCC值改为0。这样，做成一张横向为81个癌症Batch名，纵向为25 851个甲基化cg位点的表格。通过这个表格，可以看出cg位点与每一个癌症Batch的关系。使用cluster软件，对这个表格进行筛选，筛选出来的甲基化位点数为2 420个。将这个横向为81个Batch，纵向为2 420个甲基化位点的表格用cluster进行欧式距离的双向聚类。

1.2.5 相关性数据离散化分析

对甲基化位点与Batch相关性分析做进一步分析，把PCC值分为五个区间，这五个区间分别是－1～－0.4，－0.4～－0.1，－0.1～－0.1，0.1～0.4，0.4～1.0。其中，－1～－0.4区间代表强负相关；－0.4～－0.1区间代表弱负相关；－0.1～0.1区间代表不相关；0.1～0.4区间代表弱正相关；0.4～1.0的区间，代表强正相关。按照这五个部分所占百分比做成饼图。

1.2.6 筛选强、弱以及无相关位点－转录本对

筛选各癌症强相关位点数与弱相关位点数。对于每一个位点－转录本对，样本总数为N个，大于等于0.4的样本值有X个，如果X/N大于等于0.5，则认为该位点对癌症是有强调控作用的；如果大于0.1的样本值有Y个，如果Y/N大于等于0.5，则认为该位点对癌症是有弱调控作用的。

1.2.7 绘制综合癌症数据与所有单个数据的韦恩图

将挑选出来的综合癌症数据，以及BRCA癌症数据，GBM癌症数据，KIRC癌症数据，KIRP癌症数据的基因转录本，画韦恩图。

1.2.8 GO注释

将筛选出来的强相关位点，弱相关位点对应的癌症关联的转录本，放入DAVID里进行GO注释，查看其生物学途径，分子功能，细胞组件。

表达数据差异筛选与甲基化数据差异筛选并进行GO注释。在下载下来的11套数据中，其中BRCA癌症数据，OV癌症数据，READ癌症数据中同时含有正常样本与癌症样本。对这三套数据，进行癌症样本表达数据和正常样本表达数据的差异筛选，使用SAM方法进行差异筛选。再分别做这三套癌症的甲基化数据的差异基因筛选。把筛选出来的差异表达数据转录本与差异甲基化数据转录本放入DAVID里进行GO注释，查看与其相关的生物学途径，分子功能，细胞组件，并进行GO分类富集分析。在做各癌症GO注释描述及GO分类富集分析的时候，认为P－Value 及Benjiaminj值小于0.01时是显著的。

2 结果与讨论

2.1 癌症甲基化与表达值的PCC分布

基于TCGA的高通量的泛癌DNA甲基化和基因表达谱数据，我们利用皮尔森相关系数（PCC）对CpG位点对基因表达调控作用进行了定量。如图1所示，每一条曲线代表的是一个癌症Batch的PCC分布情况，图中共有81条重叠曲线。在弱负相关与弱正相关处出现两个峰值。当PCC呈现弱负相关时，基因出现表达，这符合我们所说的，甲基化程度低，表达程度高，然而当PCC值呈现正负相关时，也出现了一个峰值，但右边的峰值略低于左边。另外，图1出现了与其他 Batch不相似的 Batch曲线，粉色∗的曲线为READ＿7＿Bacth＿1758，紫色∗的曲线为OV＿7＿Bacth＿1141数据，蓝色∗的曲线为133＿OV＿1138数据。

根据选取六个不同的P值值域，画出了六个PCC分布曲线图，如图2所示，在P值取不同临界值的情况下，绝大多数癌症Batch的PCC分布曲线没有发生明显的变化，个别的曲线随着P值的变化而发生改变，这证明了PCC的分布情况是比较稳定的。

图1 当P＝0.05时癌症甲基化与表达值的PCC分布Fig.1 When P＝0.05，cancer methylation and exp ression of value distribution of the PCC

图2 取不同P值情况下PCC的分布Fig.2 When P value take different cases，the distribution of the PCC

2.2 Batch间相关性分析

红色∗越深，代表着相似性越显著。红色∗最显著的斜对角线是每个Batch和自身的相似性，因此最为显著。从Batch间相关性分析的聚类可以看到，图3可视图呈现出块状聚集的分布，处于相同癌症中的Batch的聚类效果比较显著。而对于不同癌症间的Batch，聚类效果不明显。不同癌症之间没有明显的联系。

2.3 甲基化位点与Batch相关性分析

使用 cluster进行欧式距离的双向聚类，再用treeview进行可视化，得到图4。红色∗的部分代表着某个位点的甲基化对该癌症Batch有调控的作用。红色∗越深表示调控的作用越显著。绿色∗的部分代表着这个位点的甲基化对该Batch的调控不显著。从横向来看，分析的是这2 420个CpG位点调控着哪些癌症Batch，从纵向来看，分析的是Batch共享哪些CpG位点的调控。由图4可以看到，一些CpG位点显著调控着所有的Batch，为所有癌症所共享；一些CpG位点显著着调控个别 Batch，而对其他的Batch调控是不显著的，是癌症特异的；有些CpG位点在图上显示对所有的Batch都没有显著的调控，这是27K数据一个不足的地方，位点信息仍存在着缺失。这个图绘制了甲基化区域的调控图谱。

图3 Batch间相关性分析Fig.3 Correlation analysis of Batch

2.4 相关性数据离散化分析

根据强负相关，弱负相关，无相关，弱正相关，强正相关这五个部分个数，画出图5这个饼图。浅蓝∗部分为强负相关区域，橙色∗部分为弱负相关区域，灰色∗部分为无相关区域，黄色∗部分为弱正相关区域，深蓝∗部分强正相关区域。从饼图上可以很明确的看出每一个部分所占的百分比。从饼图上可以看到大多数PCC值都是无相关的。弱相关部分远多于强相关部分。

图4 甲基化位点与Batch相关性分析Fig.4 Correlation analysis ofmethylation site and Batch

图5 不同区间PCC值范围所占百分比Fig.5 Different interval PCC percentage value range

2.5 筛选强相关、弱相关以及无相关位点－转录本对

对综合了所有癌症数据的excel表进行筛选之后，如图6所示，在综合了所有癌症批次数据下挑选出来的强相关位点有186个、弱相关位点16 280个，与无相关位点25 280个。其中无相关位点占61％，是绝大多数，弱相关位点占39％，而强相关位点只有186个，只有一小部分。对于单个的癌症而言，无相关位点同样占据了绝大多数，弱相关位点多于强相关位点。

图6 强相关位点、弱相关位点、无相关位点所占百分比Fig.6 Related sites，weak related sites，no relevant sites for percentage

2.6 绘制综合癌症数据与单个癌症数据的韦恩图

进一步我们研究了各癌症间共享的强相关位点的数量（见图7），可见，有些位点和其他癌症都有关联，并不是局限于某个癌症，为这几个癌症共享；有些转录本被若干个癌症所共享；而有些转录本是癌症特异的，只与该癌症相关，不调控其他的癌症。例如基因Ddx43的转录本NM＿018665，被这所有的五个数据集合共享，Ddx43与个体死亡有关。基因Dynlrb2的转录本NM＿130897，同时被 BRCA与 KIRP共享，Dynlrb2调控动力蛋白。基因LAPTM5的转录本NM＿006762，同时被GBM与KIRC共享，LAPTM5和溶酶体multispanning膜蛋白5有关。基因kazald1的转录本NM＿030929同时被BRCA，GBM，KIRC共享，kazald1和Kazal－type丝氨酸肽酶抑制结构域1有关。基因SLC7A2的转录本NM＿001008539，为BRCA所特有，SLC7A2与溶质载体家族7有关。fgf1基因的转录本NM＿033136，为GBM所特有，fgf1与纤维原细胞生长因子1有关。基因EPHA7的转录本NM＿004440，为KIRC特有，EPHA7与EPH受体7有关。基因SERPINE2的转录本NM＿006216，为KIRP特有，SERPINE2与serpin肽酶抑制剂有关。分别研究调控所有数据的癌症基因，以及癌症特异的基因。

图7 综合癌症数据、BRCA、GBM、KIRC、KIRP数据韦恩图Fig.7 Comprehensive cancer data，BRCA，GBM，KIRC，KIRP for venn

2.7 筛选强、弱以及无相关位点－转录本对

筛选出所有单个癌症数据及所有综合数据的符合条件的弱相关位点对、强相关位点对以及无相关位点对。在所有癌症中强相关位点数有186个，弱相关位点数16 280个。

表2 筛选出的各癌症强相关位点数与弱相关位点数Table 2 Select all the cancer related points and weak related points

2.8 筛选强相关位点、弱相关位点的甲基化与表达数据的差异位点数

筛选强相关位点、弱相关位点的甲基化与表达数据的差异位点数，在表达数据里差异表达的位点数，在甲基化数据里也是差异的。

表3 强相关位点、弱相关位点的甲基化差异位点数与表达数据差异位点数Table 3 M ethylation and expressing differences sites of strong related sites and weak related sites

2.9 各癌症GO注释描述及GO分类富集分析

对有癌症样本与正常样本的三套癌症数据进行表达数据的差异筛选。其中乳腺癌强相关表达数据差异位点数为77个，弱相关表达数据差异位点数为4 731个；卵巢癌强相关表达数据差异位点数为54个，弱相关表达数据差异位点数为1 000个；直肠癌强相关表达数据差异位点数为97个，弱相关表达数据差异位点数为3 583个。从基因层面分析各个癌症的转录本，将癌症转录本以癌症为单位放入DAVID中进行GO注释描述、GO分类富集分析。

如表4所示，BRCA癌症里有4 808个转录本，READ癌症里有3 680个转录本，OV癌症有1 026个转录本。对乳腺癌癌症在生物学过程BP＿1层面的GO注释，发现乳腺癌表达数据转录本富集在细胞过程等基本生物学过程上（P，Benjamini＜0.01），除此之外，可以看到这些基因转录本对乳腺癌癌症有发育的作用，影响癌细胞的增殖；使得癌细胞有附着力，能附着在组织或者器官上；影响着再增殖的过程，使得癌细胞无限繁殖继而继续生长；富集的基因有能使癌细胞移动的能力，使得癌细胞扩散到其他组织中，并且会导致死亡等。

如表5所示，对于卵巢癌癌症表达数据在生物学过程BP＿2层面的GO注释（P，Benjamini＜0.01），发现卵巢癌癌症表达数据集合不仅分别富集在细胞周期，细胞分裂，细胞凋亡等功能上，而且有些基因对生物学过程，细胞过程有负调控作用。一些基因注释在细胞扩散的功能上，解释了卵巢癌癌症癌症细胞在病人身上发生扩散和转移的现象。在表中，还可以看到有些基因注释为刺激细胞产生反应，也就是说当癌症发生时，这些基因的作用为刺激癌细胞，使得癌细胞产生各种机体反应。

表4 BRCA癌症在BP＿1层面的GO注释描述Table 4 GO annotation description of BRCA in BP＿1

表5 OV癌症在BP＿2层面的GO注释描述Table 5 GO annotation description of OV in BP＿2

对于直肠癌癌症表达癌症在生物学过程BP＿1层面的GO注释（P，Benjamini＜0.01），可以看到直肠癌表达数据基因集合注释在生长，增值的功能上，在富集基因的作用下，促进癌细胞不断增值，发育。注释在粘附的功能上，使得癌细胞粘附在器官或组织上，得以进一步的分裂，增值，又可以看到，癌细胞增长的同时，机体对刺激发生了反应，又促进了免疫学的过程。注释在运动的功能上，这些富集基因的功能促进了癌细胞的转移和扩散到其他器官和组织上。

3 结论与讨论

3.1 结论

本研究通过整合TCGA的泛癌症基因组的高通量甲基化谱和基因表达谱，识别癌基因相关的DNA甲基化调控位点；结果表明仅有一半的CpG位点对下游基因具有调控作用；且存在癌症特异的甲基化调控位点；并揭示这些位点调控的基因确实参与了癌症发生发展相关的功能。

3.2 讨论

不同癌症Batch的PCC值分布曲线是相似的并且稳定，但是在直肠癌与卵巢癌里有3套Batch的PCC分布曲线出现异常，这可能是数据量过少或者数据不完善的原因造成的。

癌症的强相关位点数量远小于弱相关位点，并且大部分的位点是无相关的，这代表在27 K芯片测的启动子区域数据是有许多遗漏的，仍需完善。研究分析可知，有一些位点稳定的调控着所有的癌症，与所有的癌症都有这关联。有些位点是癌症特异的，只与这些癌症有关联，可以进一步分析这些启动子区域的位点是如何影响这些癌症的发生。有些位点与若干个癌症相关联，有的位点没看出对其他癌症有调控作用。

被所有癌症共享的基因很少，而癌症的发生往往不是只受到一个基因的影响，而是分别由几个基因共同作用而产生的。同时癌症也也受到特异的基因的影响。后续研究可以分别挑选这些不同类别基因进行研究。

通过将筛选出的癌症差异基因，放入DAVID中，做以癌症为单位的三套癌症表达数据的GO注释，得出结论，这些筛选出来的差异基因在生物学过程上，确实是与各癌症有着密切的关联。本研究的结果是对当前甲基化调控位点集的重要补充，也是识别癌症新型分子标记特征的重要资源。

（

）

［1］ JONESP A.Functions of DNAmethylation：islands，start sites，gene bodies and beyond［J］.Nature Reviews Genet⁃ics，2012，13（7）：484－492.

［2］ FAN G.DNA methylation and its basic function［J］.Neu⁃ropsychopharmacology Reviews，2012，38（1）：23－38.

［3］ SHAMESD S.DNA methylation in health，disease，and cancer［J］.CurrentMolecular Medicine，2007，7（1）：85－102（18）.

［4］ DAY J J，SWEATT JD.DNA methylation and memory formation［J］.Nature Neuroscience，2010，13（11）：1319－1323.

［5］ WU H，ZHANG Y.Reversing DNA methylation：mecha⁃nisms，genomics，and biological functions［J］.Cell，2014，156：45－68.

［6］ REA M，ZHENGW，CHEN M，et al.Histone H1 affects gene imprinting and DNA methylation in Arabidopsis［J］. Plant Journal，2012，71（5）：776－786.

［7］ ZALA D，HINCKELMANN M V，YU H，et al.Vesicular glycolysis provides on⁃board energy for fast axonal trans⁃port［J］.Cell，2013，152（3）：479－491.

［8］ SUN H S，KENNEDY P J，NESTLER E J.Epigenetics of the depressed brain：role ofhistone acetylation andmethy⁃lation［J］.Neuropsychopharmacology Official Publication of the American College of Neuropsychopharmacology，2013，38（1）：124－137.

［9］ LRY T J，LID，WALTER M J，et al.DNMT3A muta⁃tions in acutemyeloid leukemia［J］.New England Journal of Medicine，2010，363（25）：2424－2433.

［10］RUSICIO A D，EBRALIDZE A K，BENOUKRAF T，et al.DNMT1⁃interacting RNAs block gene⁃specific DNA methylation［J］.Nature，2013，503（7476）：371－376.

［11］GUO X，WANG L，LIJ，et al.Structural insight into au⁃toinhibition and histone H3⁃induced activation of DNMT3A ［J］.Nature，2015，517（7536）：640－644.

［12］COPPIETERSN，DIERIKS B V，LILL C，et al.Global changes in DNA methylation and hydroxymethylation in Alzheimer′s disease human brain［J］.Neurobiology of Ag⁃ing，2014，35：1334－1344.

［13］AKHAVAN⁃NIAKIH，SAMADANI A A.DNA methyla⁃tion and cancer development：molecular mechanism［J］. Cell Biochemistry＆Biophysics，2013，67（2）：501－513. ［14］ARAND，SABATOS，HELLMAN A.DNAmethylation of distal regulatory sites characterizes dysregulation of cancer genes［J］.Genome Biology，2013，14（3）：2242－2254.

［15］RICKETTSC J，MORRIS M R，GENTLES D，et al. Methylation profiling and evaluation of demethylating ther⁃apy in renal cell carcinoma［J］.Clinical Epigenetics，2013，5（1）：16－16.

［16］BAEK S J，YANG S，KANG TW，et al.MENT：Methyl⁃ation and expression database of normal and tumor tissues ［J］.Genes，2013，518（1）：194－200.

Identification of cancer DNA methylation regulatory sites

WEIYunzhen1，LIU Xiaojuan2，WANG Fang1，SU Jianzhong1，ZHANG Yan1∗，LIU Hongbo1∗
（1.College of Bioinformatics Science and Technology，Harbin Medical University，Harbin 150081，China；
2.Department ofRehabilitation，The First Affiliated Hospital ofHarbin Medical University，Harbin 150001，China）

DNA methylation is an important epigenetic modification，which plays an important role in the regulation of gene transcription.Abnormal DNA methylation may lead to cancer and disease，and identifying oncogene⁃related DNA methylation gene regulatory sites is important for the development ofmechanisms to resolve the occurrence of cancer and identify new cancermarkers.In this study，we integrate high⁃throughput DNA methylation profiling and gene expression profiling of pan⁃cancer genome in TCGA，then identify oncogene⁃related DNA methylation regulation sites.For each cancer，we calculate the correlation betweenmethylation of CpG sites and gene expression，and filter the CpG sites，which regulate downstream genes（including strong regulatory sites，weak regulatory sites and not regulatory sites）.The results show thatonly half of the CpG sites regulate the downstream genes.Analyzing of regulatory sites that is shared between cancers show that regulatory sites are not necessarily the same in different cancer，and the presence of cancer⁃specificmethylation regulatory sites.Moreover，gene function enrichmentanalysis of differential DNA methylation and differentially expressed genes show that genes regulated by methylation are indeed involved in the development of cancer⁃related functions.The results of this study are an important supplementation to the current DNA methylation regulatory sites set，and an important resource to identify new molecularmarkers characteristics of cancer.

DNA methylation；Gene expression；Transcriptional regulation；Cancer

R73；Q7

B

1672－5565（2015）03－170－09

2015－05－27；

2015－07－20.

国家自然科学基金项目（61403112，31371334）。

韦云真，女，本科生，研究方向：计算表观遗传学；E⁃mail：weiyunzhen＠yeah.net.

∗

张岩，女，教授，研究方向：计算表观遗传学、生物信息学；E⁃mail：tyozhang＠ems.hrbmu.edu.cn；刘洪波，男，讲师，研究方向：计算表观遗传学、生物信息学；E⁃mail：hongbo919＠gmail.com.