iMeta | 南科大杨亮/刘洋-开发研究细菌DNA甲基化组及其转录调控工具Bacmethy

Bacmethy：一种新颖、便捷的研究细菌DNA甲基化模式及其转录调控作用工具

iMeta主页：http://www.imeta.science

引  言

在原核生物中，DNA甲基化修饰是表观遗传调控的主要形式，该过程由DNA甲基转移酶（DNA MTase）介导。DNA甲基转移酶将甲基基团从供体S-腺苷甲硫氨酸转移至目标碱基的位置，从而形成不同的修饰形式，主要包括N-6-甲基腺嘌呤（细菌中最常见的类型）、N-4-甲基胞嘧啶和5-甲基胞嘧啶。

研究发现，DNA甲基化可通过不同的分子机制影响基因的转录过程。比如，甲基化修饰后的碱基会进入到DNA双螺旋的主沟，从而影响与DNA结合蛋白的相互作用，这种对基因转录的影响既可以是抑制的也可以是促进的。在某些情况下，蛋白质与DNA的结合可以防止其发生甲基化，从而产生低甲基化或非甲基化现象，其特点是甲基化水平较低。此外，甲基化修饰还可能会影响染色体的拓扑结构。因此，通过详细解读甲基化组的特征，可以预测DNA甲基化对转录调控以及生理活动的影响。

近十几年来，通过SMRT-seq技术，对细菌DNA甲基组功能的研究显著增加。SMRT-seq是一种在全基因组水平实现单核苷酸分辨率来检测细菌DNA甲基化修饰的工具。到目前为止，这项技术已被用来生成了大量细菌和古菌的基因组和甲基化组的数据。SMRT-seq原始测序数据采用二进制格式存储，其包括原始序列信息和脉冲间持续时间信号。这些信息包含了模板DNA内部发生的共价修饰事件。使用如SMRT-Tool和SMRTLink等生物信息学工具可对SMRT-seq原始数据进行分析，实现检测DNA甲基化修饰位置并通过机器学习获取特定基序（motif）序列。然而目前仍缺乏便捷的生物信息学工具来系统地表征DNA甲基化组特征以及基于SMRT-seq数据来预测转录调控效应。

本研究开发了分析细菌甲基化组的工具—Bacmethy。Bacmethy使非生信专业的用户能够进行全面的细菌甲基化组分析，并促进DNA甲基化修饰影响转录调控及重要表型方面的理解。

结果和讨论

Bacmethy的使用方式

用户可以在本地使用bash命令或在网站上运行Bacmethy（图1和图S2），都需要提供样本名称、甲基化修饰检测文件（motifs.csv和motifs.gff）以及完整的基因组序列（fasta文件）。在其他可选参数中，用户可以调整甲基化水平的阈值，计算得到在不同阈值的条件下的低甲基化和非甲基化分类结果。另一个可选分析是转录因子结合预测，用户需准备转录因子序列矩阵。在Bacmethy的运行时间方面，运行时间与每个基因组中的甲基化位点数量相关。比如，对于一个有4000个甲基化位点的菌株（例如P. aeruginosa TBCF），网页服务器的标准Bacmethy分析的运行时间将少于10分钟。

图 1. Bacmethy分析流程示意图

该示意图展示了包括试验前数据准备、基本分析和可选分析。结果包括“甲基化修饰图”、“甲基化相关基因表”、“转录因子结合预测整合表”、“甲基化组的Circos图”、“功能富集”等。图中TF指转录因子。

图 S2. Bacmethy网站提交页面截图

用户可以访问https://bacmethy.med.sustech.edu.cn网址依次上传甲基化修饰检测文件（motifs.csv和motifs.gff）以及完整的基因组序列（fasta文件），并添加要分析的甲基化类型和物种名称，提交运行标准的bacmethy分析。此外，可以点击Advanced options实现多种参数以及转录因子结合预测分析选项设置，依次设置后，再点击运行即可实现分析。结果将以邮件的形式返回给用户提交任务的邮箱。

不同甲基化修饰水平的分布特征

为深入探究细菌基因组中精确的甲基化组特征，，我们首先确定了细菌基因组内各DNA MTase基序在不同基因区域（调控区或编码区）的3种修饰水平（甲基化、低甲基化或非甲基化）的分布。我们的数据显示，基因组中不同MTase 基序修饰位点的分布存在显著的多样性；我们进一步绘制了不同菌株中各MTase基序的甲基化分数分布（图2A和图S3A）。对于检测的菌株中的大多数MTase基序，仅检测到一小部分低甲基化位点。然而，在4个基序中(大肠杆菌K12中的1个基序，芒吉诺特梭菌LM2中的1个基序，金黄色葡萄球菌USA300中的2个基序)观察到大量的低片段甲基化位点。

接下来，我们还生成了不同菌株中各DNA MTase的甲基化和非(低)甲基化基序在菌株间的质量分布散点图（图2B和图S3B）。其中每个位点的甲基化特征由3个变量（identificationQV、甲基化fraction和reads覆盖度）表示。各菌株中不同motif的甲基化fraction分布是相当多样化的。除了非甲基化motif位点外，各菌株还存在数量不等的低甲基化motif位点。例如，S. aureus USA300中检测到了极高水平的低甲基化motif位点。

图 2. 甲基化修饰水平和测序质量的分布

（A）细菌菌株中每个motif的甲基化修饰水平分布图，其中每个圆圈的标题包括菌株名称和识别motif，计数指示具有不同fraction的motif数量：0 - 60%（红色），60 - 75%（灰色），75 - 90%（蓝色），90 - 100%（绿色）。（B）甲基化、低甲基化和非甲基化motif在各菌株中的质量分布散点图（在P. aeruginosa TBCF和E. coli K12共有6个甲基转移酶识别的motif）。x轴显示reads覆盖度，y轴显示identification QV。甲基化fraction由点的颜色变化表示，蓝色表示高fraction甲基化修饰，红色表示低fraction甲基化修饰。QV，质量值。

此外，为了分析调控区域内低甲基化和非甲基化位点的分布特征，我们绘制了位置富集分布图（图3）。这些图显示了转录起始位点附近有明显的峰值富集，提示甲基化修饰水平与转录起始调控之间存在联系。这些甲基化水平的动态变化可能是甲基化修饰与DNA结合蛋白的竞争所引起。

图 3.在不同菌株中每个MTase识别基序的调控区域定位的甲基化位点数量

纵坐标轴显示甲基化（蓝色）、非甲基化（红色）或低甲基化（绿色）位点的计数。横坐标轴表示修饰位点相对于ATG起始密码子的位置。

不同修饰水平的甲基化基序的功能富集

为了确定DNA甲基化所影响的基因的潜在功能，我们对两株铜绿假单胞菌（P. aeruginosa）（TBCF和LYSZa7）的DNA甲基化相关基因（指在转录调节区域具有相应DNA甲基化修饰的基因）进行了功能富集分析。如图4所示，不同菌株中甲基转移酶修饰motif的功能富集结果表现出显著差异。值得注意的是，编码转录调控因子以及含有DNA结合域的基因被发现富集，这揭示了DNA甲基化对转录调控的调节作用。

图 4. 铜绿假单胞菌菌株中甲基化相关基因的功能富集结果

该图中涉及到的甲基化相关基因指的是在其转录调控区域具有DNA甲基化修饰的基因。（A）TBCF中的TRGANNNNNNTGC motif 位点；（B）TBCF中的RCCANNNNNNNTGAR motif 位点；（C）LYSZa7中的CCCGAG motif 位点；（D）LYSZa7中的AGGNNNNNRTGT motif 位点。

预测甲基化修饰与转录因子的互作及影响

为了预测甲基化修饰是否通过转录因子影响转录调控过程，我们可以扫描甲基化位点附近的转录调节因子的结合位点或基序(transcription regulation effector binding sites，TFBS)。DNA甲基化与转录因子结合位点的密切相关的存在可能意味着基因表达受到上述两个方面的共同调控。以大肠杆菌K12的Dam甲基转移酶为例，我们从含有甲基化位点的基因的调控区域找到了潜在的转录因子结合位点。通过类似的方法，我们预测了多个共同调控事件，其中一些已在之前的研究中验证（图5）。在一项基于体内甲基化酶保护的功能实验中，Tavazoie等人检测并量化了发生在GATC甲基化修饰的基因上游区域存在一些转录因子结合序列。包括转录因子Fur结合到fepB，转录因子Crp结合到mtlA，转录因子Crp结合到flhD和转录因子Crp结合到gcd。同时，我们的Bacmethy结果也展示了相似的转录因子结合和DNA甲基化修饰位置模式（图5A, 5B）。此外，我们的分析预测了许多之前未知的甲基化与转录因子共同调控位点。例如，oppA的启动子区域包含两个GATC基序，其甲基化修饰fraction分别为0.444和0.946。Fur结合位点在这些甲基化位点的附近被发现，暗示着二者的共调控模式（图5C）。此外，我们还发现DNA甲基化发生在转录因子自身的调控区域的实例。例如，我们在LexA操作子的启动子区域检测到由Dcm介导的CCWGG修饰，已有研究证明，该甲基化修饰会影响大肠杆菌中LexA介导的SOS响应（图5D）。

图 5. 基于Bacmethy分析的大肠杆菌K12中DNA甲基化与TFBS（转录因子结合位点）位置图

（A）fepB基因中GATC motif（三角形）和Fur结合位点（绿色矩形）的位置图（灰色箭头）。（B）mtlA基因中GATC motif（三角形）和CRP结合位点（绿色矩形）的位置图（灰色箭头）。（C）oppA基因中GATC motif（三角形）和Fur结合位点（绿色矩形）的位置图（灰色箭头）。（D）lexA基因中CCWGG motif（三角形）和LexA结合位点（绿色矩形）的位置图（灰色箭头）。

利用Bacmethy和RNA-seq分析建立甲基组转录组多组学调控模型及试验验证

利用Bacmethy分析和RNA-seq分析构建甲基组-转录组整合调控模型

我们对一株临床铜绿假单胞菌TBCF10839进行SMRT-seq测序来获得其甲基组。共找到两个分别由DNA甲基转移酶AP和CP修饰的甲基化motif（RCCANNNNNNNTGAR和TRGANNNNNNTGC）。转录组分析揭示了CP酶基因的敲除突变株（Δcp）中nosR的表达显著降低，该基因编码一种一氧化氮还原酶（NOR）的转录调控因子。Δcp突变体在产生一氧化氮的RAW 264.7巨噬细胞中的生存能力较低，并在大蜡螟幼虫感染模型中的毒力减弱。然而，甲基化如何调控基因转录并影响毒力的详细机制尚未完全理解。通过使用Bacmethy流程分析，我们在nosR基因的启动子及编码区域中检测到三个甲基化motif，以及转录因子Dnr或sigma因子RpoN的结合位点（图6A）。在nosR基因上，这三个甲基化位点与转录因子的位置关系分别是，第一个邻近DNR结合位点，第三个邻近RpoN结合位点，而第二个未检测到邻近的转录因子结合位点（图6A）。

甲基组-转录组整合调控模型的实验验证

我们为nosR中的3个甲基化motif构建了定点突变株（图6B）。TBCFnosRsdm1有一个G到C的突变，位于nosR的ATG起始密码子下游60个碱基；TBCFnosRsdm2有一个G到A的突变，位于ATG起始密码子下游624个碱基；TBCFnosRsdm3有一个G到C的突变，位于ATG起始密码子下游1209个碱基。接下来，我们进行了反转录-定量PCR分析，发现在靠近转录因子结合位点的定点突变株（sdm1和sdm3）中，nosR基因的表达显著下调，但在sdm2中无明显变化（图6C）。这一发现表明，甲基化位点和转录因子可能共同控制一氧化氮还原酶基因的表达，从而增强该菌的抗巨噬细胞吞噬能力和毒力。

图 6.对铜绿假单胞菌TBCF菌株的Bacmethy分析结果进行实验验证

（A）P. aeruginosa TBCF中的TRGANNNNNNTGC motif（三角形）、DNR结合位点（绿色矩形）和RpoN结合位点（蓝色矩形）在nosRZDFYL操纵子基因（灰色箭头）上的位置图。在nosR基因中检测到三个甲基化位点，并识别出其邻近的DNR和RpoN结合位点。（B）展示了点突变菌株的设计（TBCF_nosRsdm1、TBCF_nosRsdm2和TBCF_nosRsdm3）。（C）TBCF及三个点突变菌株的RT-qPCR结果。每列呈现了每组三个生物学重复的均值和标准偏差（SD）。**，p < 0.01；***，p < 0.001。使用独立的两样本t检验进行比较。

方  法

Bacmethy流程

Bacmethy流程包括以下几个关键步骤：在预备步骤中，需要用户下载或者利用SMRTLink工具基于SMRT-seq数据生成甲基化修饰motif序列识别的文件，供Bacmethy流程使用。第一步：基因注释与定位，通过PROKKA软件进行基因组注释，并提取调控区域。第二步：DNA甲基化位点的分类，将区分不同的甲基化类型及状态（甲基化、非甲基化或低甲基化）以进行详细分析。第三步：提取可能受到甲基化调控的候选基因，即甲基化发生在其调控区域或者编码区域，然后可以对在调控区域中不同甲基化修饰水平调控的候选基因进行下游分析。第四步：转录因子结合位点的结合预测分析，采用FIMO工具扫描DNA甲基化的基因调控区序列，识别潜在的转录因子结合位置，并输出同时受DNA甲基化和转录因子调控的候选基因（RTMG表格及列表）。

可选分析一：绘制甲基组Circos图 。Bacmethy使用bash中的自脚本生成所有Circos配置文件，通过Circos可视化基因组水平的甲基化分布。可选分析二：与转录组数据的相关性分析。筛选出在野生型菌株与甲基转移酶敲除菌株中差异表达的基因（DEG）与此前分析的RTMG的重叠部分，可以作为DNA甲基化介导的基因调控的证据，可特别关注紧密间隔并发挥相似作用的基因家族。可选步骤三：基因集富集分析。用户可以使用上述步骤中生成的基因列表进行功能富集分析。可使用COG和DAVID等分类注释工具。可使用BlastP和类似的比对工具进行蛋白质功能注释。并可使用R软件包ggplot2将结果可视化。

图S1. Bacmethy的整体分析流程

Docker和Webserver构建和使用

我们使用 Docker（https://www.docker.com/）平台来构建、管理和运行Bacmethy，因此，支持不同操作系统（macOS、Linux 或 Windows）进行访问使用。要运行 Bacmethy，用户可以使用 Docker pull命令从 Docker Hub（https://hub.docker.com/repository/docker/liujihong/bacmethy）下载 Bacmethy 镜像。与此同时，我们也用VUE 2.6.14 和 KOA 2.7.0开发了网页版本的Bacmethy ，支持用户在线使用，访问链接为https://bacmethy.med.sustech.edu.cn，用户须提供从 SMRT 测序获得的 DNA 甲基化结果文件（motifs.csv 和 motifs.gff），以及完整的基因组文件，格式为 FASTA。提交任务运行完成后，用户将收到一封包含 Bacmethy 结果的电子邮件。

数据分析与实验验证

本研究使用Bacmethy对来自7种菌株的13个不同甲基转移酶的甲基化序列进行了分析，数据部分来源于NCBI数据库，部分由本团队自测。简而言之，细菌菌株在Luria-Bertani（LB）培养基中培养，提取基因组DNA，然后进行SMRT测序。PacBio reads使用SMRTLink软件v9.0的HGAP4流程以默认设置进行组装成完整基因组，完整的基因组序列被用作参考序列。使用SMRTLink的Base Modification Analysis和Motif Analysis应用程序进行DNA甲基化分析。生成的输出文件（genome.fasta, motids.gff和motifs.csv）用于Bacmethy分析，并使用REBASE数据库预测和分配已鉴定的甲基转移酶的识别基序。此外，我们利用铜绿假单胞菌TBCF和及其甲基转移酶缺失突变株的RNA-seq数据集来建立细菌转录与Bacmethy分析结果中基因甲基化位点分布的相关性。接下来，基于前期对铜绿假单胞菌TBCF中甲基转移酶影响细菌反硝化作用的研究，使用pk18modsacB质粒构建了TBCF10839中nosR上三个甲基化位点的定点突变菌株，之后提取RNA并进行RT－qPCR检测相关基因转录水平的变化，进一步验证距离转录因子结合位点较近的甲基化基序能够影响相应基因转录。

讨  论

Bacmethy基于PacBio SMRT测序数据来预测细菌甲基化对转录调控的全局影响。该分析将全面扫描基因组，识别所有潜在甲基化位点及该位点所分布的基因区域，并采用三步流程对可能具有转录调控功能的甲基化位点进行分析：一是DNA甲基化修饰水平分类，二是确定甲基化位点修饰位置，三是转录因子结合预测共调控作用。

Bacmethy具备高度的通用性，能在启动子区域以外的编码序列和基因间区域识别DNA修饰，从而提供超越传统启动子调控的广泛信息。例如，甲基化修饰可能改变DNA分子曲率，进而影响基因转录。Bacmethy不仅辨识非甲基化状态，即修饰质量值低于设定阈值，还考虑了低甲基化修饰状态，提示种群内的异质性或转录因子与甲基转移酶相互竞争的动态调控模式。Bacmethy通过命令行脚本、Docker及用户友好的网页工具三种方式方便用户使用及访问。

用户在使用bacmethy的基础上，另可以通过整合上游甲基组分析工具（例如SMRTLink、REBASE）及下游功能富集软件，能够对细菌菌株的甲基组进行全面解析，并预测与DNA甲基化相关的潜在调控作用。研究者还可通过湿实验，如不同条件下的甲基转移酶敲除菌株的基因表达量测定和表型分析等，验证Bacmethy分析得到的关键调控通路。DNA甲基化与转录因子结合的相互作用可通过甲基化位点定向突变、实时定量PCR分析、电泳迁移率移位试验等方法进行实验验证。

利用公共数据集及自己的数据，本文展示了Bacmethy在揭示DNA甲基化修饰的独特性和准确预测其对转录调控影响方面的功能应用。在大肠杆菌K12中，我们的方法也成功预测了众多共调控事件，并且某些预测已由先前研究验证过。此外，我们还发现了新的调控机制，提示了以前未知的调控关系。通过结合Bacmethy和铜绿假单胞菌DNA甲基转移酶敲除菌株及对照菌株的转录组数据，我们构建了一个进行细菌甲基组-转录组整合调控分析的模型，通过预测及实验验证DNA甲基化位点与转录因子共同调节一氧化氮还原酶相关基因的表达，为该菌株抗吞噬能力提供了解释。

（▼ 点击跳转）

“iMeta” 是由威立、肠菌分会和本领域数百千华人科学家合作出版的开放获取期刊，主编由中科院微生物所刘双江研究员和荷兰格罗宁根大学傅静远教授担任。目的是发表所有领域高影响力的研究、方法和综述，重点关注微生物组、生物信息、大数据和多组学等。目标是发表前10%(IF > 20)的高影响力论文。期刊特色包括视频投稿、可重复分析、图片打磨、青年编委、前3年免出版费、50万用户的社交媒体宣传等。2022年2月正式创刊发行！发行后相继被Google Scholar、SCIE、ESI、PubMed、DOAJ、Scopus等数据库收录！2024年6月获得首个影响因子23.7，位列全球SCI期刊前千分之五(107/21848)，微生物学科2/161，仅低于Nature Reviews，同学科研究类期刊全球第一，中国大陆11/514！

“iMetaOmics” 是“iMeta” 子刊，主编由中国科学院北京生命科学研究院赵方庆研究员和香港中文大学于君教授担任，是定位IF>10的高水平综合期刊，欢迎投稿！