论文阅读(十四)：贝叶斯网络在全基因组DNA甲基化研究中的应用
1.论文链接：Bayesian Networks in the Study of Genome-wide DNA Methylation
摘要：
本章探讨了贝叶斯网络在基因组规模的脱氧核糖核酸（DNA）甲基化研究中的应用。它首先描述了DNA甲基化的基因组规模注释的不同实验方法。详细介绍了Methyl-seq方案，并描述了该技术引起的偏倚，这对进一步分析构成了许多挑战。这些挑战通过引入贝叶斯网络框架来分析Methyl-seq数据来解决。对之前的模型进行了扩展，以纳入来自基因组序列的更多信息。基因组结构被用作甲基化状态的先验。一个反复出现的主题是用于从技术中收集信息的模型与驱动模型规范的甲基化观点之间的相互作用。最后，描述了一项研究，其中使用了这样的模型，导致有趣的生物学结论和见解的性质甲基化。
关键词：DNA甲基化，甲基序列，贝叶斯网络
本章探讨了贝叶斯网络在基因组规模DNA甲基化研究中的应用。它首先描述了DNA甲基化的基因组规模注释的不同实验方法。详细介绍了甲基-Seq方法，并描述了该技术引起的偏差，这些偏差构成了进一步分析的许多挑战。这些挑战的解决引入贝叶斯网络框架的甲基-Seq数据的分析。这个先前的模型被扩展以包含来自基因组序列的更多信息。基因组结构被用作甲基化状态的先验：未甲基化位点倾向于聚集，并且未甲基化位点比其他位点更保守。一个反复出现的主题是用于从技术中收集信息的模型与驱动模型规范的甲基化观点之间的相互作用。最后，一项研究中描述了这样的模型被使用，导致有趣的生物学结论和见解的性质甲基化。
14.1表观遗传学导论
表观遗传机制通过可遗传但可能可逆的基因表达调控来影响表型。这些机制在许多水平上起作用，包括DNA甲基化，组蛋白修饰，核小体定位和复制定时[9]。所有这些都已被证明是可遗传的跨细胞分裂和跨不同的组织和细胞类型[14，18，19，21，43，48]。有一个正在进行的研究对这些不同的机制，在本节中，我们总结了两个最研究的表观遗传特征：DNA甲基化和组蛋白修饰。
DNA甲基化DNA甲基化涉及甲基与胞嘧啶（C）核苷酸的共价连接，取代嘧啶环上5位上的氢原子，或甲基与腺嘌呤（A）的加成。腺嘌呤甲基化仅在细菌基因组中发现[18]，术语“DNA甲基化”通常指胞嘧啶甲基化。细胞中的主要事件，如细胞分化，X染色体失活和逆转录转座子沉默仅举几例，其特征在于DNA甲基化，在小鼠中，已经表明甲基化调节酶的抑制导致胚胎死亡[33，55]。广义地说，DNA甲基化与异染色质状态和转录沉默有关。在脊椎动物中，DNA甲基化主要限于CpG位点（C，后面是鸟嘌呤（Gs））。大多数脊椎动物基因组是甲基化的，未甲基化的位点倾向于沿着基因组沿着聚集在一起[4，55]。许多这样的未甲基化簇发生在启动子区域，并且它们的甲基化与附近基因的沉默相关[60]。
组蛋白修饰组蛋白的N端尾部受到各种类型的翻译后共价修饰，包括赖氨酸和精氨酸甲基化，赖氨酸乙酰化，泛素化和丝氨酸磷酸化[30]。一个组蛋白可能在给定的时间内具有许多不同的修饰，从而产生许多可能的修饰构象，有时与组蛋白密码有关[52，57]。组蛋白修饰可以通过调节组蛋白的DNA结合能力和招募进一步的染色质重塑复合物来影响DNA的包装组装[30]。通过这样做，这些修饰可以影响基因表达;一些修饰与转录抑制有关，另一些则与转录激活有关。这两种类型的修饰可以存在于启动子或基因内区域[8]。
表观遗传现象可以影响基因调控，并通过一种与它们启动时不同的机制在细胞分裂之间遗传。这导致了在细胞分裂中维持调节程序的能力，从而实现了一种形式的细胞“记忆”。因此，表观遗传学处于细胞分化研究的最前沿[10，27，39]，许多表观遗传突变与各种疾病有关[46]。特别是，许多研究报告了甲基化状态改变与各种癌症之间的关联[3，13，28]。
另一个活跃的研究领域涉及破译环境因素对表观遗传状态的影响，沿着表观遗传状态变化的随机程度。例如，已经表明，单卵双胞胎在一生中积累了DNA甲基化的差异，并且积累的差异影响了他们的基因表达谱[15]。另一项有趣的研究表明，产前烟草烟雾暴露会影响胎儿的DNA甲基化[7]。
在小鼠的几个基因座中观察到DNA甲基化的跨代遗传[6，41，44]。在植物中，跨代遗传的DNA甲基化的变化已被证明是相当频繁的，似乎是使基因调控适应不断变化的环境的常见方式[20，23，47]。这些现象对遗传和进化理论的影响目前正在探索中[51]。
高通量测序技术（也称为“下一代”测序）的出现使得在全基因组范围内研究表观遗传特征成为可能，从而产生了新的研究领域。全基因组数据集允许表征所有已知基因的表观遗传现象，并允许基因，组织，个体和物种之间的广泛比较研究[26，27，35，36，62]。基因组规模的表观遗传数据的存在使关联研究能够在寻找疾病的因果变体时同时观察基因组和表观基因组特征，并且正在开发表观基因组范围关联研究的方法[45]。
在本章中，我们将探讨贝叶斯网络在基因组规模DNA甲基化研究中的应用。我们在第14.2节中开始描述了DNA甲基化基因组规模注释的不同实验方法，并重点讨论了在一种称为甲基测序的特殊方法的分析中存在的挑战。在第14.3节中，我们介绍了一个贝叶斯网络来分析甲基-序列数据，在第14.4节中，我们扩展了这个模型，以纳入更多来自基因组序列的信息。一个反复出现的主题是用于从技术中收集信息的模型和驱动模型规范的甲基化观点之间的相互作用。我们在第14.5节中描述了一项使用这些模型的研究，并在第14.6节中展示了这些模型是如何得出有趣的生物学结论和对甲基化本质的深入了解的。
14.2下一代测序和DNA甲基化
高通量测序对许多生物学领域的变革性影响不能被夸大。最初出于降低产生个体或物种DNA序列的成本的目标，以每碱基对（bp）的低成本对许多短DNA片段进行测序的能力在分子生物学中发现了意想不到的应用。在过去的几年里，测序技术引发了一场革命，其中一个领域是表观遗传学，特别是DNA甲基化。在本节中，我们描述了高通量测序用于研究DNA甲基化的不同方法，并解释了为什么今天使用的方法相对较多。在本节中，我们将重点讨论通过使用高通量测序研究DNA甲基化的可能性和挑战。我们将首先概述能够在全基因组范围内测量DNA甲基化的不同技术，然后描述使用这些技术的主要方法，详细解释甲基测序协议，这是本章的重点。
有几种方法可以检测DNA甲基化，不同技术的详细描述可以在[31]中找到。目前可用的主要技术如图14.1所示，包括：
酶消化：这种技术使用的限制性内切酶除了序列特异性外，还对甲基化敏感（只有当它们的识别位点未甲基化时才被切割）。一个值得注意的例子是HpaII，它位于第二个胞嘧啶未甲基化的CCGG位点（在本章中，我们将一个C后面跟着一个鸟嘌呤（G）表示为CpG，为了便于标记，将一对C后面跟着一对G表示为CCGG，而不是CpCpGpG）。这种酶的消化模式可以提供DNA甲基化的读数。这种技术的一个缺点是可能不完全消化，对于某些酶，缺乏位点特异性分辨率。
亚硫酸氢钠转化：用亚硫酸氢钠处理DNA分子将未甲基化的Cs转化为尿嘧啶[24，59]。在PCR步骤之后，尿嘧啶被测序为胸腺嘧啶（Ts）。这个过程能够将表观遗传标记转化为基因组序列中的修饰。该技术的缺点包括序列复杂性的降低（参见下文的细节）、需要化学处理DNA（这可能导致降解）以及不完全转化的可能性。
A组富集：在这个过程中，基因组被消化，并且消化物富集甲基化区域。这是通过抗体[40]或甲基-CpG结合结构域蛋白[63]识别甲基化区域并将结合的DNA区域与未结合的DNA区域分离来实现的。这种技术的缺点包括酶的结合偏差和缺乏位点特异性分辨率。
14.2.1测定全基因组DNA甲基化
用于DNA甲基化的全基因组测量的方法通常使用与高通量测序或阵列杂交结合的上述技术之一。潜在地，方法可以结合多于一种的技术，但是当前的方法不这样做。这些方法可以广泛地分类为甲基分型与甲基组测序，类似于DNA的基因分型与基因组测序。在基因分型中，仅测定个体核苷酸的一小部分（SNP位置），而在基因组测序中，对个体的整个基因组进行测序。基因分型的优势在于与全基因组测序相比，其成本显著降低，允许在固定成本的实验中包含更多的个体。类似地，甲基分型技术允许通过对CpG位点的子集进行采样来调查基因组规模的甲基化模式，并且以牺牲高分辨率为代价强调低成本。在本节中，我们讨论甲基化组测序和甲基化分型的不同方法，沿着与每种方法相关的优点和缺点。
用于高通量基因组学的不同方法，包括用于测量DNA甲基化的方法，减少了对使用阵列或测序获得的消化物中的DNA片段进行计数的能力。阵列被认为不如用于定量目的的测序准确，主要是由于杂交的可变性在阵列中引入的偏差。所描述的每一种技术都可以进行阵列杂交，但到目前为止，丰度富集是最适合阵列的技术[31]。使用阵列的优点是成本低，但缺点包括在杂交步骤中引入的偏差，以及基因组规模的基于阵列的方法不是位点特异性的（有基于阵列的方法，使用酶消化或亚硫酸氢盐转化并且是位点特异性的，但由于它们采样的CpG数量少得多，我们不认为它们是基因组规模的甲基化方法）。
在过去的几年中，高通量测序的进展带来了大量增加的可用方法在全基因组范围内定位DNA甲基化。全基因组亚硫酸氢盐测序（WGBS）涉及随机消化基因组，然后进行亚硫酸氢盐处理，扩增和测序，从而能够以单核苷酸分辨率测量DNA甲基化的绝对水平。虽然这种方法已经在不同的研究中使用[11，34，35，62]，但它的使用是有限的，因为它是最昂贵的DNA甲基化注释方法，因为它需要对整个基因组进行测序。此外，测序整个甲基化组比测序整个基因组要昂贵得多，因为甲基化现象的连续性（我们想确定消化物中胞嘧啶被甲基化的细胞比例），需要比基因组测序更高的覆盖率。因此，在可预见的将来，这种方法不能用于大规模的比较研究，如人口研究和表观遗传关联研究。该方法在全基因组测序中不存在的其他缺点包括在聚合酶链反应（PCR）扩增步骤中引入的偏差和序列复杂性的降低。PCR扩增步骤引入了偏差，这是由于未甲基化的情况下引入AT丰富的序列。序列复杂性降低，因为读段中测序的每个T都可以映射到参考基因组中的T或C，从而减少了产生唯一映射读段的位置的数量。注释整个基因组中甲基化的第二种方法是高通量测序后的高通量富集。这种方法已经在[37]中使用，但需要大量的测序，不是位点特异性的，并且易于结合偏倚。
高通量测序的使用使得多种甲基化分型方法成为可能，本文讨论了其中两种：甲基Seq和简化代表性双硫酸盐测序（RRBS）。在甲基化-Seq中，基因组用甲基化敏感性限制酶HpaII消化，该限制酶HpaII消化CCGG位点内的未甲基化CpG。随后是片段的大小选择、PCR扩增步骤和片段末端测序。将测序读数映射回参考基因组揭示了未甲基化的CpG。methyl-Seq是一种方便的甲基分型策略，因为它是成本有效的，因为它只对未甲基化的位点（脊椎动物中的少数CpG位点[55]）进行测序，只需要少量的材料，并且避免了亚硫酸氢盐转化。然而，尽管该实验相对简单，但测序数据的解释并不简单，因为该方案导致基因组的非随机分割，随后是大小选择步骤。我们将在下一节中更详细地描述甲基-Seq方案以及与之相关的偏倚。
第二种甲基分型方法RRBS是基于用甲基化不敏感的酶消化，然后进行亚硫酸氢盐测序[38]。该程序是酶消化，然后是片段的大小选择，亚硫酸氢盐处理，PCR扩增和片段末端测序。通常用MspI进行消化[22，39]，其使CCGG位点富集，以富集富含CpG位点的基因组区域。RRBS是有利的，因为它比全基因组双亚硫酸盐测序便宜得多，但也有与使用双亚硫酸盐转换的方法相同的缺点：PCR扩增步骤引入的偏倚和降低的序列复杂性。
当设计一个在全基因组范围内测量DNA甲基化的实验时，必须考虑方法成本和它产生的覆盖类型之间的权衡。除此之外，不同方法的应用和分析还因许多其他问题而变得复杂。上面讨论了基于他们使用的技术的不同方法的主要并发症。最重要的是，不同检测的分析要求在难度上各不相同，测序领域的快速发展要求对方法之间的比较进行频繁的重新评估。出于这些原因，随着测序技术的变化，方法的利弊不断变化。在[31]中，有人建议甲基-Seq是最有利的利弊分析方法，就选择用于比较的措施而言（见[31]的表2）。鉴于该领域的当前阶段，我们在这里比较不同的方法，但读者应记住，随着该领域的发展和新技术的引入，比较标准将发生变化。
表14.1总结了方法比较的主要特征。它分为甲基分型方法和全基因组方法。由于领域的性质不断变化（测序成本正在迅速下降），我们省略了一列成本比较，但似乎可以肯定的是，在可预见的未来，全基因组方法将需要比甲基测序和RRBS更多的测序，RRBS将需要比甲基测序更多的测序。Methyl-seq检索的信息比RRBS跨越更多的基因组，因为相对于MspI，HpaII产生的片段大小更有利[49]。
14.2.2甲基测序法
在本节中，我们将描述甲基-Seq实验及其引入的实验偏倚。在甲基测序实验中，基因组DNA被甲基化敏感的限制性内切酶消化;在我们的例子中，是HpaII。HpaII在第二个C未甲基化的CCGG位点上酶切。该步骤获得DNA片段的消化物，使得在所有片段末端存在未甲基化的HpaII位点，并且片段内的所有HpaII位点（不包括末端）被甲基化（假设完全消化）。然后，使用在琼脂糖凝胶上的电泳对片段进行大小选择，其中可接受的常见长度在50和300 bp之间。当使用Illumina机器时，这种大小选择对于实现良好的测序通量是重要的。从通过大小选择步骤的片段构建配对末端文库，并且在PCR扩增步骤之后，对片段的末端进行测序。然后使用支持双端读段的比对器（例如Bowtie）将配对测序读段映射到参考基因组[32]。请注意，通过构建配对末端文库，测序实验返回成对的序列，其中每对是对一个片段的末端进行测序的产物。使用该方案，可以推断哪些片段存在于消化物中。
消化步骤完成后，所有未甲基化的HpaII位点都存在于消化片段的末端（见图14.1），但测序方案所需的大小选择步骤将测序限制在窄大小范围的片段。这导致在许多情况下，不能仅基于其读段计数确定位点甲基化的程度。图14.2（第370页）显示了三种不同类型的表观等位基因和由它们产生的片段，这些片段通过了大小选择步骤并被测序。当确定突出显示位点的甲基化状态时，如果仅考虑源自该位点的片段，则情况2和3是不可区分的：在这两种情况下，没有测序的片段在任一末端具有突出显示的位点。然而，在情况2中，我们看到突出显示的位点两侧的位点都是未甲基化的，因此如果突出显示的位点已经甲基化，我们将测序长度为60 bp的片段，对应于中间的片段（如情况3中测序的）。由于这样的片段没有被测序，我们可以得出结论该位点未甲基化。在情况3中，由于突出显示的位点存在于片段的内部，我们可以得出结论，它是甲基化的。虽然本图中的例子假设了二进制甲基化状态，但当涉及细胞群体时，甲基化状态被假设为连续变量，因为即使在单个细胞类型的群体中，甲基化状态也可能是异质的[1，64]。可以容易地看出，上述实例可以推广到连续变量的情况，其中不能确定单独给定该位点处的读段计数的位点被甲基化的程度，但是当利用从该位点的邻域产生的片段时可以这样做。
14.3贝叶斯网络在甲基序列分析中的应用
我们在前一节中已经看到，虽然甲基-Seq是一种有利的甲基分型方法，但该方法中存在的偏倚需要校正。当进行甲基测序实验时，我们希望了解每个HpaII位点甲基化的程度。更确切地说，我们想知道每个HpaII位点在该位点被甲基化的消化物中细胞的比例。然而，我们观察到的是实验中每个片段被测序的次数。使用我们的知识的实验程序，我们采取生成的方法来模拟的程序，通过该程序的HpaII位点的甲基化状态的影响的次数，每个片段被测序。在此生成模型的基础上，我们可以推断出预期的甲基化状态的HpaII网站，观察到的片段计数。
在下一节中，我们将讨论如何将甲基测序实验建模为贝叶斯网络，以及如何使用贝叶斯网络来推断给定片段计数的所有HpaII位点的甲基化程度。我们开始通过引入基于贝叶斯网络的数据生成模型。然后，我们描述了如何使用期望最大化（EM）算法来学习模型参数，以及我们如何将其与观察到的片段计数一起使用来推断HpaII位点的甲基化状态。在下面的部分中，我们将讨论如何通过添加一组额外的变量来扩展模型，以实现更好的位点特异性推断，以及推断未甲基化簇的位置。
14.3.1符号
14.3.2生成模型
14.3.3参数学习和后验概率的推断
EM方法
14.4基因组结构作为甲基化状态的先验
EM算法的应用
14.5应用：人类神经元甲基化分型
在前面的章节中，我们已经描述了用于校正高通量测序实验中引入的偏倚的贝叶斯网络，其中基因组以非随机方式消化，然后进行大小选择步骤，并且重点关注甲基-Seq实验。在本章中，我们详细介绍了一项研究的结果，该研究应用这种模型来表征来自四个人类个体的中性粒细胞的甲基型[49]。
在[49]中，作者表明，甲基测序应用的甲基分型方法结合贝叶斯网络校正程序（称为MetMap）足以调查整个基因组的甲基化状态并提供了显着的洞察甲基化，内部和外部的CpG岛，在网站特异性分辨率。在手头的研究中，使用标准甲基-Seq方案确定了四个男性个体的甲基型。位点特异性甲基化值可以分配给人类参考基因组中4.8%的CpG，其中20%在CpG岛内。
作者表明，通过使用MetMap的校正所实现的精度比使用原始读段计数的精度大大提高。这是通过使用直接亚硫酸氢测序确定46个HpaII位点的真实甲基化状态，并观察到MetMap推断得分与亚硫酸氢验证的Pearson相关性（0.90）显著高于原始片段计数估计的甲基化值与亚硫酸氢验证的Pearson相关性（0.67）。
14.5.1未甲基化的簇
正如我们在本章中所讨论的，脊椎动物基因组中的未甲基化位点往往聚集在一起。虽然这使得能够从基因组序列预测这样的簇的位置，但如果实验甲基化数据可用，则期望直接注释未甲基化的岛。这是因为DNA甲基化在不同的组织和条件下是不同的，实验数据将揭示特定实验的甲基化状态的更精确的图像。
在[49]中，作者确定了实验特定的未甲基化岛屿，并将其称为SUMI（强烈未甲基化岛屿）。艾德。这些未甲基化的簇对手头的实验是特异性的，因此对人类中性粒细胞也是特异性的。一个有趣的测试是研究在中性粒细胞实验中发现的SUMI中有多少也存在于基于序列的方法推断为未甲基化的区域中（见图14.5）。作者将他们发现的20986个SUMI与通过两种基于序列的方法推断为未甲基化的区域进行了比较：使用加州大学圣克鲁斯（UCSC）浏览器发现的CpG岛（称为UCSC岛）和“bona de”CpG岛（称为BF岛）[5]。UCSC岛是最常用的基于序列的岛集合，而BF岛使用更复杂且可能更精确的方法进行注释。在20986个SUMI中，有4652个不与UCSC岛重叠，7055个不与BF岛重叠。
显示从人嗜中性粒细胞样品之一推断的位点特异性甲基化评分（上图）和未甲基化簇（SUMI，第二图）的基因组部分。对于位点特异性评分，评分0确定位点为完全甲基化。第三和第四幅图分别显示了[5]注释的BF岛屿和UCSC岛屿。虽然SUMI和通过基于序列的方法推断的岛之间存在大量重叠，但在该图中可以看到一些新的SUMI，其中一个在转录起始位点。RefSeq表示国家生物技术信息中心参考序列数据库中注释的基因。
与基于序列的方法相比，作者可以通过使用直接亚硫酸盐测序来验证SUMI集的灵敏度和特异性。使用这种方法，可以在几百个bp的长度上注释所有CpG位点的甲基化状态。作为对特定区域的验证，这些区域被注释为UCSC岛和BF岛的一部分，并且不与SUMI重叠，进行双测序。在中性粒细胞样品中发现所有区域都被甲基化。作为灵敏度的验证，对与SUMI和UCSC岛重叠但与BF岛不重叠的四个区域进行了双硫代测序。在中性粒细胞样品中发现所有区域均未甲基化。这些结果表明，特定组织中真正的甲基化状态可能不同于单独由基因组序列推断的甲基化状态，也可能不同于其他细胞类型中存在的甲基化状态。
有趣的是，3797个SUMI与BF岛或UCSC岛不重叠，揭示了中性粒细胞中未甲基化的新区域。这些区域与其他功能特征相关的区域，如开放染色质，高序列保守性，基因区域和基因的5 ′端，表现出强烈的关联。表14.2显示了不同岛屿与这些特征的关系。作者指出，与所有SUMI相比，新型SUMI在距离转录起始位点一定距离处映射的SUMI比例更大，但新型SUMI与开放染色质的关联程度与所有SUMI相似;这表明新型SUMI可能经常代表远端调控序列。
14.6结论
在本章中，我们已经展示了贝叶斯网络如何用于DNA甲基化的基因组规模研究。通过设计一个生成模型，该模型将生物状态的隐藏变量和生成的实验数据的观察变量结合起来，我们能够预测生物状态，即使这种推断并不简单。
对贝叶斯网络的需求来自感兴趣的实验（甲基-Seq），该实验包括基因组的非随机消化，然后是大小选择步骤。贝叶斯网络适用于这种情况，因为它们可以用来模拟相邻站点之间的依赖关系。该模型可以进行调整，以解释其他几种甲基化方法中的类似偏倚，包括使用几种不同的甲基化敏感性限制性内切酶。它也有望在其他容易出现此类偏差的高通量测序实验中证明是有用的，例如最近开发的一些用于确定RNA二级结构的方案[58]。
参考文献
略