CN102776270A - 检测dna甲基化的方法和装置 - Google Patents

Abstract

本发明公开了检测DNA甲基化的方法，用于检测DNA甲基化的装置。所述方法包括获得细胞和组织的甲基化DNA富集样本的测序数据的步骤和筛选可信甲基化位点的步骤；所述装置包括用于执行获得细胞和组织的甲基化DNA富集样本的测序数据的步骤的组件；和用于执行筛选可信甲基化位点的步骤的功能模块。

Description

检测DNA甲基化的方法和装置

技术领域

本发明涉及DNA检测领域，更具体地说，本发明涉及检测DNA甲基化的方法和装置。

背景技术

一.DNA甲基化与癌症

5-甲基胞嘧啶是在DNA转甲基酶的催化作用下，某些胞嘧啶被转移上一个甲基所得的产物(图1)，许多人称为“第五种碱基”[1]，它的形成过程被称为“DNA甲基化”。研究证实，在哺乳动物和植物基因组中，包含约1％-6％的5-甲基胞嘧啶[2]，它通过参与基因表达的调控而影响生物体的发育、遗传、细胞的新陈代谢等过程[3]，它的重要作用在近年来已经逐渐被人们重视，成为新的研究热点。

根据目前的研究结果，5-甲基胞嘧啶在转录水平抑制基因表达的机制主要有：1.阻碍转录复合物(CF、CA等)与转录因子(E2F等)结合[4]；2.与甲基化结合蛋白相结合，激发去乙酰化酶作用，改变染色质的结构，阻止转录进行[5]；3.在其他条件的诱导下，甲基化区域发生去乙酰化作用，组蛋白空间上变得紧密，使染色质结构改变从而抑制转录[6]。

DNA甲基化与癌症发生有密切的关系，早在上世纪80年代，研究人员利用高效液相色谱方法(High performance liquid chromatography，HPLC)，比较研究肿瘤细胞和正常细胞的5-甲基胞嘧啶的含量，发现肿瘤细胞，尤其是开始发生癌转移的恶性肿瘤细胞基因组中，DNA甲基化水平显著低于正常细胞[7]。LDubeau，Efiala等人检测了19组卵巢癌与正常卵巢表皮细胞的甲基化水平，有17组癌细胞的全局DNA甲基化水平低于正常细胞[8]，同样，Cheng等发现卵巢癌细胞中的5-甲基胞嘧啶含量与正常细胞相比，平均少了25％[9]。

基因组的全局低甲基化现象在许多不同种类的癌症中都有报道，包括前列腺癌[10]，白血病[11]，肝癌[12]，宫颈癌[13]等等，但同时局部的高甲基化也与疾病相关。

在哺乳动物中，DNA甲基化只发生在CG二核苷酸的胞嘧啶上[14]，在他们的基因组中，存在一些CG特别密集的区域，被称作CpG岛(CpG island)，国际上，对于CpG岛有一个统一的定义[15]：1.长度大于500bp；2.CG含量超过55％；3.CG出现数目的观察值与期望值的比大于0.65。

据统计，人类基因组中含有约3万个CpG岛。在正常组织中，位于基因启动子区域的CpG岛通常是非甲基化的，然而在肿瘤细胞中，经常会观察到这些区域的过甲基化现象[14]。Costello等对98对原发性肿瘤细胞(包括乳腺癌，结肠癌，睾丸癌，肝癌，肺癌等)与正常组织细胞中1200个随机选取的CpG岛做研究[16]，发现了不同程度的高甲基化现象，在同类型的肿瘤细胞中，甲基化水平也存在极大差异，甲基化程度最高的是结肠癌，而睾丸癌相对较低，只有12个CpG岛过量甲基化。发生高甲基化的CpG岛通常都位于基因的5’端附近，或是转录起始位点，并且，不局限于那些已知的癌症相关基因，一些与肿瘤似乎并不相关的基因也有过甲基化，如MYOD等[17]，另外，乳腺癌细胞中，核糖体基因5’端CpG岛也被报道存在高水平甲基化[16、18]。

二.DNA甲基化与非侵入性产前诊断

非侵入性产前检验对于包括妊娠期间并发症和胎儿遗传缺陷在内的妊娠相关状况的早期检测至关重要，因为这样能对母亲和胎儿进行早期医学干预。然而，常规检测方法是侵入性的，对母亲和胎儿均有风险。国家卫生局(TheNational Health Service)最近指出在侵入性羊膜穿刺和绒膜绒毛取样(CVS)检验后有1-2％的流产率。Lo等[19]和美国专利6,258,540指出可在母体血浆和血清中检测循环无细胞胎儿核酸，已经开发了这些侵入性方法的备选方法进行产前筛选，例如检测胎儿异常。PCT专利申请国际公布号WO2010/033639描述了胎儿和母亲之间显示不同CpG甲基化模式的基因组DNA，因此，可根据母体样品中核酸的甲基化状态检测和定量测定所述样品中的胎儿核酸。

三.DNA甲基化的研究方法

在过去几年，人们已经发明出一系列方法用于DNA甲基化的研究，这些方法可以大致概括为两步，第一步是对基因组DNA进行处理，甲基化位点和非甲基化位点处理后的结果会有区别，而第二步是用实验手段让这种区别体现出来。

1.基因组DNA处理方法：

目前的基因组DNA处理方法大致有以下三种：1)甲基化特异性酶切法。其依据是，一些限制性内切酶的敏感位点包含CG位点，当胞嘧啶被甲基化修饰后，会影响酶对位点的识别，因此，只有非甲基化区域会发生酶切作用，甲基化位点完全不受影响。此方法所使用的两个经典酶对是Hpa II-MspI(CCGG)和Sma I-Xma I(CCCGGG)[20]，由于第二对限制酶识别序列非常罕见，所以一般使用Hpa II-MspI(CCGG)。这种方法简便易行，但受到酶切位点序列的限制，如对于Hpa II-MspI来说，只能检测CCGG形式的序列区域，应用范围比较狭窄。2)亚硫酸盐转换法。亚硫酸盐能够将未甲基化的胞嘧啶转为尿嘧啶[1]，经过PCR扩增反应后，尿嘧啶即转为胸腺嘧啶，而对于甲基化修饰过的胞嘧啶则不会发生改变，因此从序列的差异上就能判断出哪些CG中的胞嘧啶为5-甲基胞嘧啶。此法可靠性与精确性极高，亚硫酸盐测序方法称为目前DNA甲基化研究的金标准。然而如果要检测全基因组的位点，则需耗费大量的试剂，时间与人力，成本问题是该方法的瓶颈。3)DNA甲基化结合蛋白富集法。此法类似于染色体免疫共沉淀(Chromatin ImmunoPrecipitation，ChIP)方法，先将DNA打断，用能够识别甲基化位点的DNA结合蛋白[21](如5-甲基胞嘧啶抗体，或MBD蛋白等)对DNA片段进行富集。富集区域所在染色体坐标通常就是甲基化位点。这种方法常用于全基因组甲基化谱的绘制，本发明对基因组DNA处理采用的是此方法。该方法弥补了Bisulfite-seq的缺陷，能够高通量检测全基因组的位点，然而位点识别的准确度却大大降低了。主要原因是富集度达不到预期要求，芯片质量参差不齐，同时缺少有效地芯片数据分析方法，因而只有当存在足够量的生物学重复的前提下，该方法的结果才能令人满意。

2.检测位点的方法

检测位点大致可分为低通量和高通量两种方法：

2.1低通量方法

此法包括：电泳、一代测序和PCR等，较为常用的亚硫酸盐测序(BSP)及甲基化特异性PCR(MSP)[22]即属于这类方法。

BSP方法指在DNA经过亚硫酸盐转换后，用一代测序检测序列与参考序列的差异：如果参考序列为胞嘧啶(C)，测序结果为胸腺嘧啶(T)，表示该处为非甲基化位点，若测序结果仍为胞嘧啶(C)，此处可以认定为甲基化位点[1]。

MSP方法也是由亚硫酸盐处理，然后设计两种不同的引物，分别针对甲基化和非甲基化两种情况下的序列，对DNA进行扩增，根据两种引物条件下的扩增情况即可判断是否被甲基化修饰[22]。

此类方法规模较小，但是准确性高，通常用于个别基因或CpG岛区域的甲基化水平检测。

2.2高通量方法

包括芯片检测和高通量测序两类。各大基因芯片制造商如Nimblegen、Agilent、Affymetrix等都有成熟的商业化芯片可以检测甲基化富集位点。他们主要采用的是TileArray技术，即对预先设计好的基因组区域，用平铺法设计探针[23]。杂交后，富集区域的探针信号会显著高于其余探针，用合适的统计学方法对数据进行分析，即可得到甲基化位点所在坐标。

Xiaoyu Zhang等2006年用5-甲基胞嘧啶抗体结合TileArray芯片的方法(称为MeDIP-chip)绘制了第一副拟南芥全基因组甲基化图谱[24]，此后，人们将该技术及其改进方法用于各物种全基因组DNA甲基化的研究，得到了许多重要成果。

然而，芯片检测法受到设计区域的限制，对于大型哺乳动物，目前制作工艺还不能使芯片覆盖到整个基因组，因而通常选取CpG岛和基因启动子区域进行设计。近年来高通量测序技术逐步成熟，在许多方面的研究都有取代芯片的趋势。对于DNA甲基化来说，高通量测序没有检测区域的限制，并且比芯片检测结果更为准确[25]。

目前高通量测序技术较为成熟的有Roche公司的454，Illumina公司的产品GenomeAnalyzer，和Applied Biosystem公司的SOLiD，后两种产品一次测序数据都可以达到30G以上，相当于覆盖了基因组10倍，同时错误率也可以控制在0.01％以内，这已经完全满足了我们检测需求。

Cokus等用亚硫酸盐处理后结合高通量测序方法(称为BS-seq或MethylC-seq)绘制了拟南芥的单碱基分辨率级别的甲基化图谱[26]，Meissner等人用同样的方法绘制了小鼠基因组部分区域的甲基化谱[27]，这些研究人员的尝试让人们看到全基因组高分辨率甲基化谱绘制的可能性。Ryan Lister等首次使用该方法绘制了人类H1干细胞和IMR90胎儿肺纤维原细胞的甲基化谱，他们对每个细胞系产生了将近90G数据，包含12亿条短序列(Short Reads)，覆盖了86％的基因组范围[28]。

BS-seq作为目前通量最高，准确性最好，分辨率最大的技术，其成本问题是它主要的瓶颈，对于大型哺乳动物，需要进行多轮测序才能获取一组样本的甲基化谱，如果有大量样本需要检测，该方法所耗的人力及财力是无法令人接受的。此时人们想到了一些较为折中的方法，如Capture-BS-seq，MeDIP-seq，MBD-seq等，他们牺牲了诸如覆盖度，分辨率等一些特性，来换取成本的降低。Emily等人用自行设计的微阵列，随机抓取了324个CpG岛的DNA片段，用BS-seq分析其甲基化水平，这种方法即是Capture-BS-seq。当研究者只关心基因组上一部分区域，如启动子区，重复区或是CpG岛的甲基化位点时，没有必要对整个基因组进行深度测序，用此法既能获取到单碱基分辨率的甲基化位点，又大大节省了成本[29]。MeDIP-seq是用5-甲基胞嘧啶对DNA富集，再结合高通量测序的方法，他所产生的全基因组甲基化谱分辨率较低，通常用若干可能包含甲基化位点的区域来表示，无法得到具体的5-甲基胞嘧啶的坐标，因此如果研究人员只关心DNA甲基化水平存在差异的区域，MeDIP-seq将是比较合适的选择。MBD-seq与MeDIP-seq非常类似，差别在于用甲基化结合蛋白MBD代替5-甲基胞嘧啶抗体对DNA片段进行富集。

四.甲基化结合蛋白家族(Methyl-CpG-Binding Domain，MBD)

20世纪90年代初，人们发现了两种对甲基化DNA有亲和性的蛋白质，称为MeCP1和MeCP2[30]。MeCP1是由许多复杂亚结构组成的复合蛋白，他的具体结构至今仍存在争议[31]。MeCP2的质量约为55kDa，它能特异性的结合于甲基化修饰的CpG位点[32]。通过对MeCP2结构的仔细研究，发现在其N端存在一个大约70个残基构成的功能域，MeCP2与DNA的结合主要依赖于该结构，于是将它称作甲基化结合蛋白(MBD)[33]。

研究人员用MBD的氨基酸序列在蛋白数据库中搜索，发现了一种蛋白，其N端具有一个与MBD非常相似的结构，开始人们称其为PCM1(含MBD的蛋白质)，后来简称为MBD1，体外实验表明，该蛋白能够与甲基化DNA结合，并抑制其转录[33]。通过进一步对EST数据库的搜索，人们发现了另外三种基因，他们编码出的蛋白也具有类似MBD的结构，分别被称作Mbd2、Mbd3和Mbd4[34]。序列的对位排列结果(图2)显示Mbd4和MeCP2的MBD结构接近，而Mbd1，Mbd2，Mbd3的MBD更为接近[34]。

Hendrich等人用凝胶滞后分析(gel retardation assay)方法，研究五种蛋白对甲基化DNA的亲和能力，发现在体外，只有Mbd3无法与甲基化DNA形成聚合物[34]，这个结果非常出乎人们的意料，因为Mbd3与另一种有甲基化DNA亲和能力的MBD2b蛋白有极高的相似度，尤其是他们的MBD功能域，其原因目前尚未研究清楚，不过有一些假说认为Mbd3可能是专门识别一些特征序列。

五.DNA甲基化富集测序分析方法

DNA甲基化富集测序根据所使用的蛋白不同可分为，MeDIP-seq及MBD-seq，前者使用的是5-甲基胞嘧啶抗体，后者使用的是甲基化结合蛋白(MBD)。由于其实验步骤与免疫共沉淀测序(ChIP-seq)大致相同，数据分析也大多采用ChIP-seq的方法。

1.DNA甲基化富集测序分析基本步骤

Shirley Pepke等在其2009年发表在Nature上的一篇综述中，将ChIP-seq数据分析算法总结为五个基本步骤[35]：

1.1平滑信号谱

测序所得短序列与参考序列比对后，得到每条序列的基因组坐标，扫描整个基因组可获取每个碱基的覆盖度，该值与ChIP-chip中的探针信号值相类似，于是通常我们也习惯地将覆盖度称为信号。但是，覆盖度由整数表示，是离散的，为了将其连续化，我们会采用一些平滑算法，比较常用的是移动平均值法，即计算每个位点周围一个固定长度窗口内的信号平均值，作为平滑后的信号。

1.2估计背景值

背景值是由于富集不完全造成，并且也是无法避免的。如果实验设计中有对照样本，则可以直接把对照的信号作为背景，否则，需要对信号的分布建立模型，估计出每个位点的噪声范围。

1.3设定阈值，进行初步筛选

得到背景之后，将信号值减去背景，高于预先设定的阈值，那么将此位点作为候选进入下一步分析

1.4依据一定原则过滤初筛结果

双链对称性是比较常用的原则，它的依据是富集片段属于正义链和负义链的概率相等，于是比对结果中两条链上序列数目应该接近。一些双链数目明显不等的区域很可能是由于比对错误，或扩增时引入的序列偏性导致的，可以从候选位点中剔除。

1.5计算每个位点的可信度(P-value)

可信度的评估是为了进一步筛选富集位点，通常根据信号的分布，信号与背景的差异程度来计算得出。设定合适的可信度阈值，可以减少最后结果中的假阳性率。

多数ChIP-seq分析工具及算法都遵循这五个基本步骤，只是在每一步中采用了不同的算法，以期获得更高的检出率，同时降低假发现率。MACS(Model-based Analysis of ChIP-seq)是使用较为广泛的一个算法[36]。它首先根据序列的分布计算出位点的中心，然后将序列平移至该中心附近，这一做法符合ChIP-seq的数据特点，因为通常测序长度小于富集片段的长度，因此需要将正义链序列向3’平移，负义链序列向5’端平移。平移之后，计算每个位点的信号值，用泊松分布对信号建模，计算每个位点属于富集区域的可信度。

除MACS外还有诸如Cisgenome，ERANGE，GLITR，PeakSeq，QuEST，SiSSRs，spp，USeq等其他算法，各自具有不同特点[37、38]。

ChIP-seq分析方法统计整个基因组被序列覆盖的情况，按照一定的算法找出覆盖度特别高的区域，称为峰(Peak)，而这些区域就被认为是可能的甲基化位点。在DNA富集程度较高，测序质量较好的情况下，所找到的位点确实有比较高的准确度。然而，实际应用过程中，这种方法存在以下不足：第一.没有考虑序列构成。DNA甲基化发生在胞嘧啶，对于哺乳动物，则只发生在CG双核苷酸的胞嘧啶上。但是即使某一区域不含CG，用ChIP-seq方法也有可能被认定为甲基化位点；第二.在覆盖度不高的区域，检出率相对较低。多数算法采用的是通过比较局部与总体覆盖度的差异，来判断是否是富集区域，当这种差异不是非常显著时，通常算法拒绝认为这是一个甲基化位点；第三.当需要比较两组样品的甲基化情况时，由于找到的位点位置并不对应，会带来一定的困难。总结以上三点不足，我们认为最主要的原因在于，ChIP-seq的分析方法是一种相对通用的方法，对于甲基化数据的研究固然适用，但是这样就没有利用好甲基化数据本身的特征。

2.其他算法

并非所有的DNA甲基化富集测序数据分析都采用与ChIP-seq相同的方法。Yoshinao等用MeDIP-seq研究人类乳腺癌细胞的DNA甲基化水平变化情况[39]，他们测得了将近1亿条序列，测序深度达到12倍。序列比对之后，他们将基因组划分成若干个100kb的区域，统计每个区域内的序列数目，通过与对照组对比，并设定2倍作为阈值，得到片段富集区域。这种方法有利于多个样本之间的比较，但只适用于数据量充足，基因组覆盖较完全的情况下使用。

Tomas A Down等人通过观察MeDIP-seq数据和BSP验证得到的结果，发现MeDIP-chip和MeDIP-seq数据的信号值在一定程度上与CG的密度相关，他们用分别用线性模型和二次模型对CG密度因子和信号值进行拟合，发现在多数区域都能得到较好的拟合结果。于是他们用贝叶斯概率模型计算在试验观察到的数据下，各位点的甲基化程度(用后验概率表示)，通过与BSP数据对比，发现结果吻合，他们的算法被称作BATMAN(Bayesian tool for methylationanalysis)[40]。

综上所述，本领域迫切需要一种能够快速、高通量且可靠地检测并分析全基因DNA甲基化位点的方法。

发明内容

本发明的一个目的是提供一种检测全基因组DNA甲基化位点的方法。

本发明的另一目的是提供一种检测肝癌细胞或组织的全基因组DNA甲基化位点的方法。

本发明还有另一目的是提供一种用于检测全基因组DNA甲基化位点的装置。

本发明还有另一目的是提供一种用于检测肝癌细胞或组织的全基因组DNA甲基化位点的装置。

因此，在一方面，本发明提供一种检测全基因组DNA甲基化位点的方法，该方法包括：

1)获得细胞或组织的甲基化DNA片段的富集样本以及获得片段信息的步骤；和

2)筛选可信的甲基化位点的步骤；

其中步骤2)包括：获得甲基化DNA片段在基因组上的坐标，即，染色体编号，起始位点，终止位点和方向(正义或反义链)；再进行位点识别与可信度打分，根据预先设定好的阈值，筛选出可信的甲基化位点；

其中位点识别与可信度打分包括：

a-1.划分CpG簇

将基因组上CpG位点划分为N个子集(P_i，i＝1...N)，要求每个子集内部的位点在染色体上连续排列，且相邻两个位点之间距离小于等于L，而任意两个子集中的位点之间的距离大于L，按照上述方法划分出的一个子集定义为CpG簇，CpG簇中相距最远的两个位点在染色体上所跨越的范围定义为CpG簇区域；

a-2.获得各CpG簇区域的覆盖度

对于一个CpG簇，其CpG簇区域在染色体上的起点为S_i，终点为E_i，获取所测序列中比对起点位点落于[S_i-21，E_i+21]范围内的序列子集R，1为超声破碎后片段的片段长度的期望值，获得CpG簇区域的覆盖度；

a-3划分子区域

将基因组划分为N’个区域(Q_j，j＝1...N’)，每个区域内，相邻两个CpG位点间距离小于等于L’，而任意两区域内的CpG位点间距离大于L’，并且满足L’＞L；

a-4获得各子区域的覆盖度

具体步骤同a-2；

a-5获得各CpG簇包含甲基化位点的概率

对于一个子区域Q_j，其覆盖度为C，它包含了k个CpG簇区域P_i(i＝1...k)，各CpG簇区域覆盖度为c_i，获得P_i包含甲基化位点的概率。

在优选的实施方式中，所述方法用于检测肝癌细胞或组织的全基因DNA甲基化位点。

在优选的实施方式中，所述方法的步骤1)包括：

1-1.提取和纯化待测动物细胞及组织的DNA；

1-2.DNA的片段化；

1-3.亲和层析柱的准备；

1-4.甲基化DNA的富集；和

1-5.富集样本的测序。

在另一优选的实施方式中，步骤1-3包括利用MBD2b蛋白。

在另一方面，本发明提供一种用于检测全基因DNA甲基化位点的装置，所述装置包括1)用于执行获得细胞或组织的甲基化DNA富集样本的测序数据的步骤的组件；和2)用于执行筛选可信甲基化位点的步骤的功能模块。

在一优选实施方式中，所述装置用于检测肝癌细胞或组织的全基因DNA甲基化位点。

在另一优选实施方式中，所述装置中用于执行所述细胞或组织的甲基化DNA富集样本的测序数据的步骤的组件包含MBD2b蛋白。

本发明的其它方面由于本文的公开内容，对本领域的技术人员而言是显而易见的。

附图说明

图1显示胞嘧啶在DNA转甲基酶的催化作用下，被转移上一个甲基，成为5-甲基胞嘧啶。

图2显示五种甲基化DNA结合蛋白(MeCP2、MBD1、MBD2、MBD3和MBD4)的MBD区域的序列对位排列。

图3显示覆盖度5X的情况。

图4是CpG簇一致性系数的经验累积分布图，曲线从上到下分别是L＝300、275、250、225、200、175、150、125、100、75、50和25的情况。

图5是在不同L值的条件下高一致性的CpG簇所占比例。

图6是韦恩图。

图7显示序列比对情况随错配数的累积变化。

图8是各染色体上的覆盖比例柱状图，黑色为肿瘤组织，灰色为癌旁组织(图中“chr”＝染色体)。

图9是在不同P值的阈值条件下，所预测出的甲基化CpG位点所占比例。

图10显示差异位点在各区域的比例。

图11是过甲基化修饰的外显子根据编号分类的分布直方图，编号为1的柱形对应于第一外显子的数目，之后的柱形以此类推。

图12是低甲基化修饰的外显子根据编号分类的分布直方图，编号为1的柱形对应于第一外显子的数目，之后的柱形以此类推。

具体实施方式

发明人发现利用易于获取且价格低廉的MBD2b蛋白，结合二代测序，能够高通量，准确地将富集结果与基因组位置对应，再利用我们开发的TAMD方法，所得的甲基化位点的假阳性率可以降低到0.1％。在此基础上，完成了本发明。

除非另行定义，文中所使用的所有专业与科学用语与本领域熟练人员所熟悉的意义相同。此外，任何与所记载内容相似或均等的方法及材料皆可应用于本发明中。文中所述的较佳实施方法与材料仅作示范之用。

一.定义

本文所用的术语“TAMD方法”表示MEDIP-seq和MBD-seq数据分析工具(Tools for Analyzing MeDIP-seq and MBD-seq Data)。

本文所用的术语“二代测序”又称高通量测序，其主要特点是，能够同时对输入的序列进行大规模并行测序，并且所得结果为大量的(一般为2千万左右)长度不超过200bp的短序列(454测序仪例外，他的读长大于1kbp，但序列较少，在200万左右)。

本文所用的术语“覆盖度”指将二代测序所得序列与基因组序列进行比对，得到每条短序列的基因组坐标。于是，我们可以统计每个基因组位点有多少条序列比对上，该数目便称为覆盖度，例如图3所示的是覆盖度5X的情况。

本文所用的术语“CpG簇”指基因组上连续分布且两两之间的距离不超过L的一群CpG位点。参数L控制CpG簇内位点的数目，其取值范围从1到无限大，L等于1时，一个CpG位点即为一个CpG簇，随着L增大，一些分布较密的位点逐渐合并。根据文献报道，一定区域内的CpG甲基化状态高度相关[41]，因此，当L取合适的值时，CpG簇内的绝大多数位点的甲基化状态一致，在理想情况下，我们便可用整个CpG簇包含甲基化位点的概率来代替其中某个CpG被甲基化修饰的概率。

经过MBD富集后，有一些并不含有甲基化位点的片段会残留在结果中，这些片段我们在文中称作噪声片段。而本文所用的术语“富集度”指的是除去噪声片段后，余下的片段数与总的片段数之比。一般来说，在相同的实验条件下，用同一批次蛋白得到的富集度比较稳定，它的数值可以通过实验测得，我们的实验中富集度是90％。

二.TAMD方法介绍

1.贝叶斯概率模型：

对于某个CpG簇，我们想知道，如果观察到n条片段(即，下文所述超声破碎断裂后得到的片段)落入该簇，那么有多大的概率可以认为该CpG簇被甲基化修饰。我们将这个概率记做P(M|n)，M表示甲基化，根据贝叶斯公式

$P\left(M\right|n)=\frac{P\left(n\right|M)P\left(M\right)}{P\left(n\right|\stackrel{\‾}{M})(1-P\left(M\right))+P\left(n\right|M)P\left(M\right)}$

其中

表示非甲基化位点。全局的甲基化概率P(M)，作为本模型的一个输入参数，可以从已报道的文献，或同源生物的情况推知。

对于非甲基化片段，我们假设其在整个基因组上是随机均匀分布的，因而对于一条片段来说，它落入某一区域的概率等于

d为区域的长度，D表示整个基因组的长度，因此，某一区域内落入的非甲基化片段数目符合二项分布，即

对于甲基化片段，我们假设其落在每个CpG簇内的概率与位点内甲基化胞嘧啶数目有关，根据我们对H1细胞的BS-seq和MeDIP-seq的观察，这一假设也确实成立。因此，某区域内甲基化片段数目也同样符合二项分布，即其中，nc为位点内甲基化胞嘧啶数目，NC为甲基化胞嘧啶总数。

在此，nc与NC都是未知的，我们观察到，距离较近的两个胞嘧啶甲基化状态是高度相关的，这与之前文献报道的结果也相符，因而对于一个甲基化位点来说，如果其中的胞嘧啶或CG的间隔较小，也就是L值较低的情况下，我们可以认为这个位点的甲基化胞嘧啶数目就是胞嘧啶总数，这样便大大简化了模型。接下来，对于NC，应该等于基因组内的甲基化位点数目，这在事先我们并无法获知，为此，我们采用类似于遗传算法中的淘汰机制，先假定所有位点都属于甲基化位点，计算所有位点的概率，根据预先设定的阈值，去除那些属于甲基化位点可能性较低的，剩下的数目作为新的NC，参与计算，如此循环，直到没有更多的位点被淘汰为止。

2.TAMD方法

TAMD方法的基本思想是，对全基因组上所有的候选甲基化位点一一排查，用贝叶斯推断的方法计算每个位点甲基化的可能性。这里所指的候选位点是指可以被甲基化修饰的胞嘧啶位置。例如，在哺乳动物中，指的即是所有的CG位点。用这种方法来分析甲基化富集测序数据就可以弥补之前提到的ChIP-seq分析方法的不足：首先，由于我们只专注于胞嘧啶位置，就完全避免了将不可能发生甲基化的位置认为是甲基化位点的错误；第二，通过调节模型参数，可以让方法识别出覆盖度较低区域的甲基化位点，同时也能将假阳性率控制在合理的范围。第三，对于同一个基因组，我们找出的位点属于一个固定的集合，从不同样本找出的位点之间有明确的对应关系，可以十分方便地找出差异位点。

本发明的TAMD方法包括以下步骤：

1)划分CpG簇

将基因组上CpG位点划分为N个子集(P_i，i＝1...N)，要求每个子集内部的位点在染色体上连续排列，且相邻两个位点之间距离小于等于L，而任意两个子集中的位点之间的距离大于L。按照上述方法划分出的一个子集定义为CpG簇，CpG簇中相距最远的两个位点在染色体上所跨越的范围定义为CpG簇区域。

2)计算各CpG簇区域的覆盖度

对于一个CpG簇P_i，其CpG簇区域起点在染色体上的位置为S_i，终点为E_i，获取所测序列中比对起点位点落于[S_i-21，E_i+21]范围内的序列子集R，1为超声破碎后片段的片段长度的期望值

CpG簇区域的覆盖度C_i按如下公式计算

C_i＝Floor(∑_r∈R Γ)，其中

s等于序列起始位点到CpG簇区域边缘的距离，z为超声破碎后片段长度的标准偏差，Ф为标准正态分布的分布函数，Floor(x)函数表示取小于x的最大整数

3)划分子区域

将基因组划分为N’个区域(Q_j，j＝1...N’)，每个区域内，相邻两个CpG位点间距离小于等于L’，而任意两区域内的CpG位点间距离大于L’，并且满足L’＞L。

根据定义，可推知，每个CpG簇区域必定包含于唯一的一个子区域

4)计算各子区域的覆盖度

计算方法同步骤二

5)计算各CpG簇包含甲基化位点的概率

对于一个子区域Q_j，其覆盖度为C，它包含了k个CpG簇区域P_i(i＝1...k)，各CpG簇区域覆盖度为c_i，则P_i包含甲基化位点概率计算如下：

$P_i\left(M\right|\mathrm{Obs}.)=\frac{B(c_i,C,\frac{n_i}{T})P\left(M\right)}{(1-B(c_i,(1-\mathrm{\α})S,\frac{d_i}{D}))(1-P\left(M\right))+B(c_i,C,\frac{n_i}{T})P\left(M\right)}$

$T={\mathrm{\&Sigma;}}_{i-1}^k{t}_i$ $t_i=\left\{\begin{array}{cc}{n}_i& B(c_i,(1-\mathrm{\&alpha;})S,\frac{d_i}{D})<\mathrm{\&epsiv;}\\ 0& B(c_i,(1-\mathrm{\&alpha;})S,\frac{d_i}{D})\&GreaterEqual;\mathrm{\&epsiv;}\end{array}\right.$

P(M)为基因组全局甲基化概率，n_i为CpG簇P_i内的位点数，m为Q_j的内包含的位点数，S为所测序列总数，α为富集度，ε为非甲基化情况下后验概率的阈值；

B(x，n，p)为二项分布的分布函数；

$B(x,n,p)=\underset{i=1}{\overset{i\&le;x}{\mathrm{\&Sigma;}}}{c}_n^i{p}^i{(1-p)}^{n-i}.$

下面结合具体实施例，进一步阐述本发明。应理解，这些实施例仅用于说明本发明而不用于限制本发明的范围。下列实施例中未注明具体条件的实验方法，通常按照生物化学或分子生物学常规实验的常规条件(例如参考本领域可用的教科书，比如《最新分子生物学方法》、《分子克隆实验室手册》等等)，或按照制造厂商所建议的条件，利用生物化学或分子生物学常规实验的常规试剂。除非另外说明，否则百分比和份数按重量计算。

实施例

实施例1.评估基因组上连续分布的一群CpG位点两两之间的距离L对CpG簇的甲基化状态一致性系数的影响

我们使用人类表观遗传学计划(Human Epigenome Project，HEP)所测得的人类第6，20和22号染色体的DNA甲基化谱数据，评估随L变化所产生的影响。

HEP检测了43个不同细胞样本(12种不同组织)中190万个CpG位点的甲基化状态，我们从中选取一组肝脏细胞中的数据，包含226548个CpG位点，以50作为位点甲基化状态的分数阈值(该数据中分数取值从0至100)，即超过50则认为该CpG位点被甲基化修饰。根据前述的CpG簇定义，我们分别用25，50，75，100，125，150，175，200，225，250，275，300这12个数作为L的值，将所有CpG位点划分为若干CpG簇，即，从DNA的3’端向5’端逐个计算每个CpG与相邻的前一个CpG之间的距离，如果小于L，则划分为同一个CpG簇，如果大于L，则认为是一个新的CpG簇。计算每个簇的甲基化状态一致性系数。该系数定义为，从簇内任意挑选的两个CpG位点的甲基化状态相同的概率：

$P=\frac{m(m-1)+(N-m)(N-m-1)}{N(N-1)}$

N为CpG簇内位点总数，m为被甲基化的CpG数目。

从一致性系数的经验累积分布图(图4)来看，在L小于等于200的情况下，超过一半的CpG簇内的甲基化状态完全一致(一致性系数等于1)，并且一致性较低(一致性系数小于0.5)的CpG簇数目所占比例并不高，小于10％，即，图4中虚线与曲线交点所对应的纵坐标值。

假如以一致性系数大于0.8为标准，随L值的增加，一致性较高的CpG簇的比例由80％降为70％，并且趋势逐渐变缓(图5)，如果提高一致性的标准，虽然符合要求的CpG簇比例会随之减少，但变化趋势都非常相似，L大于100之后，这一比例基本不变。

本实施例证明L取100是可行的。

实施例2.利用TAMD方法分析H1细胞的MeDIP-seq结果

我们采用TAMD方法对一组公共数据库中的H1细胞的MeDIP-seq结果进行分析，并与相应的MethylC-seq结果进行比较，评估该方法的灵敏性(用检出率表示)和准确性(用假阳性率表示)。

数据来源：http://www.ncbi.nlm.nih.gov/geo/roadmap/epigenomics/

MeDIP-seq：GSM456941

MethylC-seq：GSM429321

具体结果：

根据MethylC-seq的结果，有4,116,051个CpG簇被甲基化。我们用TAMD方法(具体步骤如下所述)从MeDIP-seq中筛选出甲基化状态的CpG簇2,391,192个：

步骤一：划分CpG簇

取L等于100，将人类Hg19基因组序列划分为8,220,526个CpG簇

步骤二：计算各CpG簇区域的覆盖度

根据DNA片段化，例如超声破碎设定的参数，设置片段平均长度1为200，标准偏差z为50，按如下公式计算CpG簇区域的覆盖度C_i：

C_i＝Floor(∑_r∈R Γ)，其中

步骤三：划分子区域

取L’为1000，将基因组划分为882,227个子区域

步骤四：计算各子区域的覆盖度

计算方法同步骤二

步骤五：计算各CpG簇包含甲基化位点的概率

设定全局甲基化覆盖度P(M)为0.1，富集度α为90％，阈值ε为0.1，按以下公式计算包含甲基化位点概率：

$P_i\left(M\right|\mathrm{Obs}.)=\frac{B(c_i,C,\frac{n_i}{T})P\left(M\right)}{(1-B(c_i,(1-\mathrm{\&alpha;})S,\frac{d_i}{D}))(1-P\left(M\right))+B(c_i,C,\frac{n_i}{T})P\left(M\right)}$

$T={\mathrm{\&Sigma;}}_{i=1}^k{t}_i$ $t_i=\left\{\begin{array}{cc}{n}_i& B(c_i,(1-\mathrm{\&alpha;})S,\frac{d_i}{D})<\mathrm{\&epsiv;}\\ 0& B(c_i,(1-\mathrm{\&alpha;})S,\frac{d_i}{D})\&GreaterEqual;\mathrm{\&epsiv;}\end{array}\right.$

B(x，n，p)为二项分布的分布函数；

$B(x,n,p)=\underset{i=1}{\overset{i\&le;x}{\mathrm{\&Sigma;}}}{c}_n^i{p}^i{(1-p)}^{n-i},$

其中有2,076,533个与MethylC-seq的结果一致，因此可以认为TAMD在对这组数据的应用中，检出率为49.55％，假阳性率为13.16％。考虑到低覆盖度区域的检测的误差较大，我们观察了随着覆盖度的增大，方法性能的变化情况。发现MeDIP-seq所测序列对应到MethylC-seq检出的4,116,051个甲基化CpG簇区域中，覆盖度超过5X的有1,465,673(N)个CpG簇区域；而MeDIP-seq检出2,391,192个甲基化CpG簇中，覆盖度超过5X的有1,393,170(M)个。根据检出率＝(M∩N)/N和假阳性率＝(M-M∩N)/M(如韦恩图(图6)所示)，可以看出，随着覆盖度超过5X，检出率提高到95％，假阳性率也降低到0.1％。

本实施例证明本发明方法具有优秀的灵敏性和准确性。

实施例3.利用TAMD方法获得人类肝癌组织全基因组甲基化谱

一.待测样本信息：

肿瘤组织和癌旁组织由上海生物芯片有限公司提供，取自同一肝癌患者，病理信息如下：

表1肝癌患者病理信息

二.实验过程

1.提取和纯化待测动物细胞及组织的DNA

A.对于细胞：取约5×10^6-7细胞(已用TBS洗过)，转移到一个1.5mL离心管中，加入约700μl的裂解缓冲液；

对于组织：取约50mg组织，放入5mL离心管中，加入700μl裂解缓冲液，匀浆

B.加入蛋白酶-K至终浓度100μg/ml，50℃下，在约3-5小时期间，不时搅拌，直至溶液清亮为止；

C.加入RNA酶至终浓度20μg/ml，37℃温育30分钟；

D.加入700μl Tris饱和酚(pH 8.0)，轻缓倒转摇匀10分钟。室温下，以13000rpm离心15分钟；

E.吸上清层入另一1.5ml离心管；

F.重复D-E步骤一次(任选的)；

G.上清液加650μl酚∶氯仿∶异戊醇(25∶24∶1)，轻缓倒转混匀10分钟。室温下，以13000rpm离心10分钟，吸上层入另一1.5ml离心管；

H.上清液加650μl氯仿∶异戊醇(24∶1)，轻缓倒转混匀10分钟，以13000rpm离心10分钟，取上清入另一管；

I.加1/10体积3M醋酸钠(注：样品量充足的情况下可不加3MNaAC，以得到更纯的样品)，混匀后再加1.1倍体积的异丙醇，轻柔振摇；出现白色絮状物(DNA)；于-20℃放置20～30分钟，12000rpm低温(4℃)离心10分钟，沉淀DNA，弃去上清；

J.加1ml 70％乙醇洗涤，12000rpm下低温(4℃)离心10分钟，弃去上清，重复1次；

K.自然干燥后，用pH 8.01/10TE溶解，保存在-20℃备用。

2.DNA的片段化；

取约200ng DNA用0.8％Agrose电流检测，电流Marker为“1kb”。基因组DNA片段大小应大于15K。将DNA定量到200ng/μl，取200μl，超声破碎到300-1000bp片段大小。

3.亲和层析柱的制备

A.用剪口枪头吸取“MBD2b蛋白-琼脂糖凝胶-4B”50μl到200μlPCR管中；

B.吸取制备缓冲液A(preperation buffer A)200μl到管中，轻微振动以悬浮MBD2b蛋白-琼脂糖凝胶-4B；

C.放入亲和层析仪，4℃下旋转2分钟；

D.500g离心30秒，吸去上清；

E.吸取制备缓冲液B(preperation buffer B)200μl到管中，轻微振动以悬浮MBD2b蛋白-琼脂糖凝胶-4B；

F.重复C-D步骤操作；

G.重复B-E步骤操作；

H.吸取结合缓冲液(binding buffer A)200μl到管中，轻微振动以悬浮MBD2b蛋白-琼脂糖凝胶-4B；

I.重复C-D步骤操作；

J.吸取结合缓冲液(binding buffer A)180μl到管中，轻微振动以悬浮MBD2b蛋白-琼脂糖凝胶-4B。

4.甲基化DNA的富集

A.将超声好的DNA样品取500ng约20μl加入干净离心管中，再加入内参DNA-M及U各2μl；

B.混匀后取后2μl备用；

C.其余样本加入放有蛋白的PCR管中；

D.轻微振动以悬浮MBD2b蛋白-琼脂糖凝胶-4B；

E.4℃下在亲和层析仪上结合2小时；

F.500g离心1分钟，将上清缓冲液转移到一新的1.5ml离心管中(备用)；

G.吸取洗涤缓冲液(wash buffer)200μl到管中，轻微振动以悬浮MBD2b蛋白-琼脂糖凝胶-4B；

H.500g离心30秒，吸去上清；

I.重复F-G步骤一次。

J.吸取洗脱缓冲液(elution buffer)200μl到管中，轻微振动以悬浮MBD2b蛋白-琼脂糖凝胶-4B；

K.500g离心30秒；

L.收集离心后的洗脱缓冲液，为富集得到的甲基化DNA；

M.得到的洗脱缓冲液用QIAquick PCR purification kit回收(具体步骤参见QIAquick PCR purification kit protocol)，60μl水洗脱。

5.富集样本的测序

遵照SOLiD3的标准流程对富集样本测序，其中用ABI的SOLiD3测序仪代替Illumina的GAII。

6.筛选可信的甲基化位点的步骤

A.测序数据的初级分析

用Illumina的OLB、CASAVA软件(http://www.illumina.com/software/)对所得数据进行图像分析、碱基识别与序列比对，最终获得片段在基因组上的坐标，包括染色体编号，起始位点，终止位点和方向(正义或反义链)。

B.位点识别与可信度打分

采用我们的TAMD方法，对整个染色体上的候选甲基化位点进行扫描，并给出甲基化水平的打分，根据预先设定好的阈值，筛选出可信的甲基化位点。

三.实验结果

(一)测序与比对结果

1.基本数据信息

我们从肝癌组织中共测得序列49,536,458条，与人类基因组hg18比对，有36,824,482(74.34％)条至少找到一个比对结果，同时，从癌旁组织中，测得序列52,153,693条，有36,921,304(70.79％)条能够比对上hg 18。

表2.测序数据基本统计信息

比对采用BioScope默认参数，允许的最大错配数为测序长度的20％，即10个位点，图7展示了随错配数增长，比对结果的累积变化，在两个样本中，错配数小于3的序列都超过超过了80％，可见测序质量正常。

为提高甲基化位点识别的准确度，我们对序列做初步筛选，保留比对位置唯一，错配数小于3的序列，同时过滤掉PCR扩增产生的重复片段。如此，我们最终得到可用的肿瘤组织序列18,061,423条，癌旁组织18,702,032条。

2.覆盖度统计

人类基因组hg18总长为3,080,436,051bp，在肿瘤组织测序结果中，有580,788,017bp(18.85％)的区域至少被一条序列覆盖，而在癌旁组织中，这一数字为651,473,523bp(21.15％)。

为检验测序是否存在偏性，我们计算了每条染色体上两组数据的覆盖比例(图8)，结果发现，除第Y号染色体和线粒体染色体(ChrM)外，其余23条染色体大致相似，两种组织之间也并无明显的差别。

从而证明这些序列的质量正常，且在基因组上无偏性。

(二)甲基化位点检测

1.用TAMD方法识别甲基化位点

对前文所述最终得到的可用的肿瘤组织序列18,061,423条，癌旁组织18,702,032条进行以下操作：

1)划分CpG簇

取L等于100，将人类Hg18基因组序列划分为8,211,882个CpG簇

2)计算各CpG簇区域的覆盖度

根据DNA片段化，例如超声破碎设定的参数，设置片段平均长度1为200，标准偏差z为50，按如下公式计算CpG簇区域的覆盖度C_i：

C_i＝Floor(∑_r∈R Γ)，其中

3)划分子区域

取L’为1000，将基因组划分为877,158个子区域

4)计算各子区域的覆盖度

方法同步骤二

5)计算各CpG簇包含甲基化位点的概率

设定全局甲基化覆盖度P(M)为0.1，富集度α为90％，阈值ε为0.1，按以下公式计算包含甲基化位点概率：

$P_i\left(M\right|\mathrm{Obs}.)=\frac{B(c_i,C,\frac{n_i}{T})P\left(M\right)}{(1-B(c_i,(1-\mathrm{\&alpha;})S,\frac{d_i}{D}))(1-P\left(M\right))+B(c_i,C,\frac{n_i}{T})P\left(M\right)}$

$T={\mathrm{\&Sigma;}}_{i=1}^k{t}_i$ $t_i=\left\{\begin{array}{cc}{n}_i& B(c_i,(1-\mathrm{\&alpha;})S,\frac{d_i}{D})<\mathrm{\&epsiv;}\\ 0& B(c_i,(1-\mathrm{\&alpha;})S,\frac{d_i}{D})\&GreaterEqual;\mathrm{\&epsiv;}\end{array}\right.$

B(x，n，p)为二项分布的分布函数；

$B(x,n,p)=\underset{i=1}{\overset{i\&le;x}{\mathrm{\&Sigma;}}}{c}_n^i{p}^i{(1-p)}^{n-i}.$

计算出各位点的甲基化概率P后，我们统计了不同概率阈值的条件下，预测出的甲基化位点所占比例(图9)，可以发现，P值呈两极分布，接近三分之一的位点的P值等于0(肿瘤组织中为38.6％，癌旁组织中为34.99％)，同时也有近三分之一位点P值大于0.95(肿瘤组织中34.86％，癌旁组织中37.08％)，余下的P值在(0，0.95)的区间内呈均匀分布。在之后的分析中，我们设定P值大于0.95的CpG位点为可信的甲基化位点，于是我们得到肿瘤组织中甲基化CpG9,816,730个，属于1,499,129个CpG簇区域，癌旁组织中甲基化CpG 10,443,541个，属于1,803,417个CpG簇区域。可见，肿瘤组织的总体甲基化水平低于癌旁组织，这与过去文献的报道相符。

(三)肿瘤细胞中的DNA甲基化异常

我们已经从两种组织细胞的DNA甲基化位点分布情况观察到了一些肿瘤组织的甲基化修饰异常现象，为进一步挖掘其中的特征和规律，我们根据TAMD算出的可信度和覆盖度，定义了差异甲基化位点，并分析他们在基因组上的覆盖规律。

1.差异甲基化位点

我们定义肿瘤细胞中甲基化可信度大于0.95，覆盖度大于3，癌旁细胞中可信度小于0.1，覆盖度小于0.5的CpG为过甲基化位点(Hypermethylation)，相应的，癌旁细胞可信度大于0.95，覆盖度大于3，肿瘤细胞中可信度小于0.1，覆盖度小于0.5的CpG为低甲基化位点(Hypomethylation)。

根据上述定义，我们共得到了107,374个过甲基化位点，属于19,624个CpG簇，77,745个低甲基化位点，属于22,840个CpG簇。从簇内位点数分布可以看出，在CpG密集区域，过甲基化现象发生的比例高于低甲基化。

2.差异位点的分布

根据差异位点的坐标，我们分别对外显子，内含子，启动子，基因间，CpG岛，及各类重复片段区域的过甲基化和低甲基化位点数目进行统计，并计算在两者位点总数内所占比例(图10)，可以看到分布最多的仍然是内含子和基因间区域，并且，除了基因间和部分种类的重复序列区，其他注释区域的过甲基化的比例都大于低甲基化比例。

在UCSC Hg18注释数据库所定义的35280个基因的编码区中，有2388(6.77％)个基因含有过甲基化修饰位点，2005(5.68％)个含有低甲基化位点。我们将包含差异位点的外显子根据其编号进行分类，统计频数并绘制直方图(图11与和图12)，可以发现，第一外显子被过甲基化和低甲基化修饰的频率都远远超过其他外显子，这可能与基因表达的调控相关。

在CpG岛区域，过甲基化的位点数目为15490，低甲基化数目为1913，两者相差近8倍，与之前由簇内位点数分布情况所得结论相符。许多文献报道过，癌症相关的过甲基化现象多发生在启动子区域的CpG岛，在我们的数据中，这部分位点数为4,433，占所有CpG岛区域过甲基化位点的28.62％，有37个基因在启动子区域存在CpG岛且含有过甲基化位点(表4)，其中BHLHB9，EPCAM，GNAS，GRIK2，KCNK2，KRT7，PDX1，PPARG，RASSF5，TF，TP73在现有文献中已被报道与肿瘤相关，且发现存在有启动子区的CpG岛过甲基化现象[42-51]。

ACP1

EPCAM

HOXD9

LOC100130148

MIR196B

PEAR1

RHBG

ADAMTSL3

GMIP

HPDL

LOC144571

MPHOSPH10

PPARG

TEX264

ADCY2

GNAS

HSPA1B

MAP1LC3B

MYO15B

PRRC1

TF

BHLHB9

GRIK2

KCNK2

MAPK4

NEUROD6

RAB36

TP73

ENTPD3

HIST2H2BA

KRT7

MAPT

PDX1

RASSF5

VSX1

ZMYM2

ZNF326

表4.启动子区存在CpG岛过甲基化的基因，其中深色标注的基因已经在文献中被报道与肿瘤发生相关

综上所述，从这幅图谱中，许多之前已经发现的规律得到了验证：如肝癌组织的整体甲基化水平低于正常组织，但在CpG岛区域，则有明显的过甲基化现象；CpG位点稀疏区域的甲基化水平高于密集区域；一些与癌症相关基因的启动子区CpG岛存在过甲基化现象等等，从而证明本发明方法所得图谱的准确性。

本发明优点：

第一.本发明所采用的方法解决了检测准确度和检测成本的矛盾，具有通量高，准确性好，成本低廉的特点。

第二.本发明采用的MBD2b蛋白易于获取，价格低于甲基化抗体，结合二代测序，能够高通量，准确地将富集结果与基因组位置对应，再利用我们开发的TAMD方法，所得的甲基化位点的假阳性率可以降低到0.1％。

第三.由于本发明能够用较低的成本，在较短时间内获取一个样品完整的DNA甲基化谱。对于大量样品的甲基化研究已经不再是难题，如果能够对获得的数据进行比较分析，得到某些疾病的特征谱，那么此方法将来可以作为一种疾病诊断的辅助手段。

以上所述仅为本发明的较佳实施例而已，并非用以限定本发明的实质技术内容范围，本发明的实质技术内容是广义地定义于申请的权利要求范围中，任何他人完成的技术实体或方法，若是与申请的权利要求范围所定义的完全相同，也或是一种等效的变更，均将被视为涵盖于该权利要求范围之中。

参考文献

1.Lister，R.和J.R.Ecker (2009).″Finding the fifth base：genome-wide sequencingof cytosine methylation.″Genome Res 19(6)：959-66.

2.Montero，L.M.，J.Filipski等，(1992).″The distribution of 5-methylcytosine in thenuclear genome ofplants.″Nucleic Acids Res 20(12)：3207-10.

3.Bestor，T.H.(2000).″The DNA methyltransferases of mammals.″Hum MolGenet 9(16)：2395-402.

4.Tate，P.H.和A.P.Bird(1993).″Effects of DNA methylation on DNA-bindingproteins and gene expression.″Curr Opin Genet Dev 3(2)：226-31.

5.Keshet，I.，J.Lieman-Hurwitz等，(1986).″DNA methylation affects the formationof active chromatin.″Cell 44(4)：535-43.

6.Newell-Price，J.，A.J.Clark等，(2000).″DNA methylation and silencing of geneexpression.″Trends Endocrinol Metab 11(4)：142-8.

7.Gama-Sosa，M.A.，V.A.Slagel等，(1983).″The 5-methylcytosine content ofDNA from human tumors.″Nucleic Acids Res 11(19)：6883-94.

8.Grunau，C.，C.Sanchez等，(2005).″Frequent DNA hypomethylation of humanjuxtacentromeric BAGE loci in cancer.″Genes Chromosomes Cancer 43(1)：11-24.

9.Cheng，P.，C.Schmutte等，(1997).″Alterations in DNA methylation are early，butnot initial，events in ovarian tumorigenesis.″Br J Cancer 75(3)：396-402.

10.Bedford，M.T.和P.D.van Helden(1987).″Hypomethyl ation of DNA inpathological conditions of the human prostate.″Cancer Res 47(20)：5274-6.

11.Wahlfors，J.，H.Hiltunen等，(1992).″Genomic hypomethylation in humanchronic lymphocytic leukemia.″Blood 80(8)：2074-80.

12.Lin，C.H.，S.Y.Hsieh等，(2001).″Genome-wide hypomethylation inhepatocellular carcinogenesis.″Cancer Res 61(10)：4238-43.

13.Kim，Y.I.，A.Giuliano等，(1994).″Global DNA hypomethylation increasesprogressively in cervical dysplasia and carcinoma.″Cancer 74(3)：893-9.

14.Lethe，B.，S.Lucas等，(1998).″LAGE-1，a new gene with tumor specificity.″IntJ Cancer 76(6)：903-8.

15.Takai，D.和P.A.Jones (2002).″Comprehensive analysis of CpG islands inhuman chromosomes 21 and 22.″Proc Natl Acad Sci U S A 99(6)：3740-5.

16.Costello，J.F.，M.C.Fruhwald等，(2000).″Aberrant CpG-island methylation hasnon-random and tumour-type-specific pattems.″Nat Genet 24(2)：132-8.

17.Samuelsson，J.K.，S.Alonso等，″DNA fingerprinting techniques for the analysisof genetic and epigenetic alterations in colorectal cancer.″Mutat Res 693(1-2)：61-76.

18.Esteller，M.，P.G.Corn等，(2001).″A gene hypermethylation profile of humancancer.″Cancer Res 61(8)：3225-9.

19.Lo等，Lancet 350：485-487，1997

20.武立鹏，朱.(2004).″DNA甲基化的生物学应用及检测方法进展.″中国检验医学杂志27：7.

21.Weber，M.，J.J.Davies等，(2005).″Chromosome-wide and promoter-specificanalyses identify sites of differential DNA methylation in normal and transformed humancells.″Nat Genet 37(8)：853-62.

22.Herman，J.G.，J.R.Graff等，(1996).″Methylation-specific PCR：a novel PCRassay for methylation status of CpG islands.″Proc Natl Acad Sci U S A 93(18)：9821-6.

23.Heisler，L.E.，D.Torti等，(2005).″CpG Island microarray probe sequencesderived from a physical library are representative of CpG Islands annotated on the humangenome.″Nucleic Acids Res 33(9)：2952-61.

24.Zhang，X.，J.Yazaki等，(2006).″Genome-wide high-resolution mapping andfunctional analysis ofDNA methyl ation in arabidopsis.″Cell 126(6)：1189-201.

25.Shendure，J.和H.Ji(2008).″Next-generation DNA sequencing.″Nat Biotechnol26(10)：1135-45.

26.Cokus，S.J.，S.Feng等，(2008).″Shotgun bisulphite sequencing of theArabidopsis genome reveals DNA methylation patterning.″Nature 452(7184)：215-9.

27.Meissner，A.，A.Gnirke等，(2005).″Reduced representation bisulfite sequencingfor comparative high-resolution DNA methylation analysis.″Nucleic Acids Res 33(18)：5868-77.

28.Lister，R.，M.Pelizzola等，(2009).″Human DNA methylomes at base resolutionshow widespread epigenomic differences.″Nature 462(7271)：315-22.

29.Hodges，E，A.D.Smith等，(2009).″High definition profiling of mammalianDNA methylation by array capture and single molecule bisulfite sequencing.″Genome Res19(9)：1593-605.

30.Ballestar，E.和A.P.Wolffe(2001).″Methyl-CpG-binding proteins.Targetingspecific gene repression.″Eur J Biochem 268(1)：1-6.

31.Meehan，R.R.，J.D.Lewis等，(1989).″Identification of a mammalian proteinthat binds specifically to DNA containing methylated CpGs.″Cell 58(3)：499-507.

32.Nan，X.，R.R.Meehan等，(1993).″Dissection of the methyl-CpG bindingdomain from the chromosomal protein MeCP2.″Nucleic Acids Res 21(21)：4886-92.

33.Cross，S.H.，R.R.Meehan等，(1997).″A component of the transcriptionalrepressor MeCP 1shares a motif with DNA methyltransferase and HRX proteins.″NatGenet 16(3)：256-9.

34.Hendrich，B.和A.Bird (1998).″Identification and characterization of a family ofmammalian methyl-CpG binding proteins.″Mol Cell Biol 18(11)：6538-47.

35.Pepke，S.，B.Wold等，(2009).″Computation for ChIP-seq and RNA-seqstudies.″Nat Methods 6(11 Suppl)：S22-32.

36.Zhang，Y.，T.Liu等，(2008).″Model-based analysis of ChIP-Seq(MACS).″Genome Biol 9(9)：R137.

37.Jiang，H.，F.Wang等，″CisGenome Browser：a flexible tool for genomic datavisualization.″Bioinformatics 26(14)：1781-2.

38.Rozowsky，J.，G.Euskirchen等，(2009).″PeakSeq enables systematic scoring ofChIP-seq experiments relative to controls.″Nat Biotechnol 27(1)：66-75.

39.Ruike，Y.，Y.Imanaka等，″Genome-wide analysis of aberrant methylation inhuman breast cancer cells using methyl-DNA immunoprecipitation combined withhigh-throughput sequencing.″BMC Genomics 11：137.

40.Down，T.A.，V.K.Rakyan等，(2008).″A Bayesian deconvolution strategy forimmunoprecipitation-based DNA methylome analysis.″Nat Biotechnol 26(7)：779-85.

41.Eckhardt，F.，J.Lewin等，(2006).″DNA methylation profiling of humanchromosomes 6，20and 22.″Nat Genet 38(12)：1378-85.

42.Gebhard，C.，C.Benner等，″General transcription factor binding at CpG islandsin normal cells correlates with resistance to de novo DNA methylation in cancer cells.″Cancer Res 70(4)：1398-407.

43.Huang，K.T.，A.Dobrovic等，″DNA methylation profiling of phyllodes andfibroadenoma tumours ofthe breast.″Breast Cancer Res Treat 124(2)：555-65.

44.Ibragimova，I.，I.Ibanez de Caceres等，″Global reactivation of epigeneticallysilenced genes in prostate cancer.″Cancer Prev Res(Phila)3(9)：1084-92.

45.Mantovani，G.，A.G.Lania等，″GNAS imprinting and pituitary tumors.″MolCell Endocrinol 326(1-2)：15-8.

46.Pancione，M.，L.Sabatino等，″Epigenetic silencing of peroxisomeproliferator-activated receptor gamma is a biomarker for colorectal cancer progression andadverse patients″outcome.″PLoS One 5(12)：e14229.

47.Sproul，D.，C.Nestor等，″Transcriptionally repressed genes become aberrantlymethylated and distinguish tumors of different lineages in breast cancer.″Proc Natl AcadSci U S A 108(11)：4364-9.

48.Wu，C.S.，Y.J.Lu等，″Glutamate receptor，ionotropic，kainate 2silencing byDNA hypermethyl ation possesses tumor suppressor function in gastric cancer.″Int JCancer 126(11)：2542-52.

49.Jacinto，F.V.，E.Ballestar等，(2007).″Discovery of epigenetically silencedgenes by methylated DNA immunoprecipitation in colon cancer cells.″Cancer Res 67(24)：11481-6.

50.Taylor，K.H.，K.E.Pena-Hernandez等，(2007).″Large-scale CpG methylationanalysis identifies novel candidate genes and reveals methylation hotspots in acutelymphoblastic leukemia.″Cancer Res 67(6)：2617-25.

51.Steinmann，K.，A.Sandner等，(2009).″Frequent promoter hypermethylation oftumor-related genes in head and neck squamous cell carcinoma.″Oncol Rep 22(6)：1519-26.

Claims (7)

1.一种检测全基因组DNA甲基化位点的方法，该方法包括：

1)获得细胞或组织的甲基化DNA片段的富集样本以及获得片段信息的步骤；和

2)筛选可信甲基化位点的步骤；

其中步骤2)包括：获得甲基化DNA片段在基因组上的坐标，即，染色体编号，起始位点，终止位点和方向(正义或反义链)；再进行位点识别与可信度打分，根据预先设定好的阈值，筛选出可信的甲基化位点；

其中位点识别与可信度打分包括：

a-1.划分CpG簇

将基因组上CpG位点划分为N个子集(P_i，i＝1...N)，要求每个子集内部的位点在染色体上连续排列，且相邻两个位点之间距离小于等于L，而任意两个子集中的位点之间的距离大于L，按照上述方法划分出的一个子集定义为CpG簇，CpG簇中相距最远的两个位点在染色体上所跨越的范围定义为CpG簇区域；

a-2.获得各CpG簇区域的覆盖度

对于一个CpG簇，其CpG簇区域在染色体上的起点为S_i，终点为E_i，获取所测序列中比对起点位点落于[S_i-21，E_i+21]范围内的序列子集R，1为超声破碎后片段的片段长度的期望值，获得CpG簇区域的覆盖度；

a-3划分子区域

将基因组划分为N’个区域(Q_j，j＝1...N’)，每个区域内，相邻两个CpG位点间距离小于等于L’，而任意两区域内的CpG位点间距离大于L’，并且满足L’＞L；

a-4获得各子区域的覆盖度

具体步骤同a-2；

a-5获得各CpG簇包含甲基化位点的概率

对于一个子区域Q_j，其覆盖度为C，它包含了k个CpG簇区域P_i(i＝1...k)，各CpG簇区域覆盖度为c_i，获得P_i包含甲基化位点的概率。

2.如权利要求1所述的方法，其特征在于，所述细胞或组织是肝癌细胞或组织。

3.如权利要求1或2所述方法，其特征在于，步骤1)包括：

1-1.提取和纯化待测动物细胞或组织的DNA；

1-2.DNA的片段化；

1-3.亲和层析柱的准备；

1-4.甲基化DNA的富集；和

1-5.富集样本的测序。

4.如权利要求3所述的方法，其特征在于，所述步骤1-3包括利用MBD2b蛋白。

5.一种用于检测全基因组DNA甲基化位点的装置，所述装置包括：

1)用于执行获得细胞或组织的甲基化DNA富集样本的测序数据的步骤的组件；和

2)用于执行筛选可信甲基化位点的步骤的功能模块。

6.如权利要求5所述的装置，其特征在于，所述细胞或组织是肝癌细胞或组织。

7.如权利要求5或6所述的装置，其特征在于，所述组件1)包含MBD2b蛋白。

Images