纳米孔测序揭示意大利地中海河流水牛乳腺炎的DNA甲基化表观遗传特征

本研究首次应用纳米孔测序技术探究意大利地中海河流水牛乳腺炎的表观遗传调控机制，鉴定出22个差异甲基化胞嘧啶位点（DMCs），发现健康组与乳腺炎组存在显著不同的甲基化模式，为开发基于DNA甲基化的新型诊断工具和抗病育种策略提供了重要见解。

　　

摘要

乳腺炎是水牛养殖业面临的主要挑战之一，对动物健康、奶品质和农场经济效益具有重要影响。其中亚临床乳腺炎（SCM）由于缺乏标准化诊断工具而难以检测。表观遗传机制如DNA甲基化在牛中已被证明与免疫反应和疾病易感性相关，但在水牛中尚未得到充分探索。本研究利用纳米孔测序技术（ONT）对意大利地中海河流水牛进行全基因组DNA甲基化分析，旨在识别与乳腺炎抗性相关的潜在表观遗传特征。结果显示CpG位点呈现单峰分布，大多数位点表现出高甲基化水平。共鉴定出22个差异甲基化胞嘧啶（DMCs），其中68%在对照组表现为低甲基化，32%为高甲基化。基因组注释显示，高甲基化DMCs主要位于内含子区域，而低甲基化DMCs则大量富集在远端基因间区域。本研究首次利用纳米孔测序探索水牛乳腺炎相关的DNA甲基化变化，揭示了健康与患病动物之间不同的表观遗传模式，为后续更大样本的研究奠定了基础。

引言

选择性育种在提高家畜生产力、可持续性和动物福利方面发挥着关键作用。在水牛养殖中，遗传选择日益关注抗病性、繁殖效率以及奶肉生产性状。特别是奶品质和产量作为意大利奶水牛产业的核心选育性状，直接影响奶制品的工艺特性以及如莫扎里拉水牛奶酪（MBC）这类具有原产地保护命名（PDO）的高附加值产品生产。

然而，乳腺炎——一种由微生物主要引发的乳腺内感染——会严重损害这些性状。它引发强烈的炎症反应，导致重大经济损失，危及动物福利，若管理不当还可能带来公共卫生风险。乳腺炎主要表现为临床型（有可见症状）和亚临床型（SCM），后者因缺乏明显症状而常被漏检，但会增加体细胞计数（SCC）。两种形式均影响产奶量和奶品质，需要兽医干预，有时甚至导致淘汰，增加养殖成本。

亚临床乳腺炎尤其隐蔽，是乳腺内感染的储库，导致产奶量显著下降。牛奶中体细胞计数（SCC）升高是乳腺炎的标志，会降解乳蛋白、缩短巴氏杀菌奶的保质期，并最终影响莫扎里拉奶酪的感官特性。目前，诊断指南包括微生物培养和SCC检测以量化乳房炎症。最近，差异体细胞计数（DSCC）被提出作为评估乳房健康状态的新指标。然而，由于水牛奶缺乏标准化的SCC阈值以及监测实践不足，使得水牛乳腺炎的管理更为复杂。

在牛中，表观遗传调控特别是DNA甲基化在乳腺炎感染中的作用受到越来越多的关注，但其在水牛中的相关研究仍十分有限。DNA甲基化是重要的表观遗传标记，参与基因调控。它通过在DNA甲基转移酶催化下将甲基基团转移到DNA序列的胞嘧啶上，形成5-甲基胞嘧啶（5mC），主要出现在CpG二核苷酸中。相邻CpG位点的甲基化改变被称为差异甲基化区域（DMR）。基因表达和转录激活通常与缺乏甲基化相关，而启动子区域的DNA甲基化则指示转录抑制。

监测这些变化为早期检测和诊断乳腺炎提供了宝贵机会。此外，特定位点的甲基化也可作为畜牧育种计划中的分子标记。纳米孔测序能够详细分析这些表观遗传修饰，并检测整个基因组中的CpG位点，包括基因间和调控区域，这对表观遗传影响健康与抗病性尚未充分探索的水牛育种具有重要意义。

方法

样本采集

从意大利南部一家奶牛场采集了100头雌性水牛的外周血和奶样。所有动物均在秋季、产后150天同一泌乳阶段采样。为减少胎次变异，仅选择第3或第4胎次的动物。记录的功能控制数据包括动物信息（ID号、产犊日期、泌乳阶段）、采样日期、日产奶量（kg/d）、乳成分（脂肪和蛋白质百分比）、SCC（cells/mL）和DSCC（%）。

牛奶分析

采样前，使用含氯己定的一次性毛巾清洁和消毒乳头及乳头末端。弃去最初的前奶，从每个水牛乳房收集约70 mL奶样于无菌塑料管中并立即冷藏（4°C）。奶样在采集后12小时内由拉提纳的Istituto Zooprofilattico Sperimentale del Lazio e della Toscana ‘M. Aleandri’实验室进行SCC、DSCC和致乳腺炎微生物的培养检测。

使用Fossomatic 7 DC仪器（Foss Electric, Hiller?d, Denmark）根据ISO UNI EN ISO 13366-2:2007测定个体奶样中的SCC和DSCC。操作按制造商建议进行，排除DSCC结果在<50×10³ cells/mL或>1,500×10³ cells/mL范围内以及SCC结果中呈现良好分离（GOSE）等于0的样本。

细菌学检查遵循标准程序进行。将每个奶样10微升涂布于血琼脂平板，37°C有氧培养48小时。使用API生化鉴定板（bioMerieux, Marcy, France）对分离菌落进行鉴定。分离出至少10个同一类型菌落视为样本阳性，而分离出三种或以上类型菌落则视为环境污染，该样本被排除研究。因此，根据这些标准对16头动物进行分类：8头阴性动物，特征为细菌培养阴性、SCC≤200,000 cells/mL、DSCC≤60%；8头阳性动物，定义为细菌培养阳性、SCC>400,000 cells/mL、DSCC>60%。

DNA提取

使用Wizard DNA提取试剂盒（Promega, Madison, WI, USA）从全血中提取基因组DNA，操作按制造商说明进行。使用Eppendorf BioPhotometer（Eppendorf AG, Hamburg, Germany）测量浓度和质量，确保产量>100 ng/μL，260/280比值>1.7，260/230比值>1.8。

文库制备与纳米孔测序

使用2 μg DNA起始输入量制备纳米孔测序文库。使用Native Barcoding Kit 24 V14 for gDNA（SQK-NBD114.24）生成带条形码的文库，遵循制造商协议。阳性和阴性组样本在文库制备过程中混合，确保每个测序运行中两种条件的平衡代表。将文库加载到R10.4.1流动槽上，每个流动槽 multiplex 3个样本。使用PromethION 2（P2）Solo设备（Oxford Nanopore Technologies, Oxford, UK）进行测序。

使用Dorado basecaller（https://github.com/nanoporetech/dorado）v0.7.3.1对原始测序数据（pod5文件）进行碱基识别。

甲基化识别

分析在西班牙加利西亚的CESGA超级计算中心进行。将Pod5文件索引并比对到水牛（Bubalus bubalis）参考基因组（NDDB_SH_1）。使用Dorado中的dna_r10.4.1_e8.2_400bps_sup@v5.0.0超精准模型进行甲基化识别，参数为–modified_bases 5mCG-5hmCG。使用Dorado中的demux命令将碱基识别数据集按条形码分类和解复用，生成存储在输出目录中的多个BAM文件。

使用modbam2bed（https://github.com/epi2me-labs/modbam2bed）v0.10.0提取甲基化信息，参数–extended以包括规范、修饰和过滤碱基的计数。

差异甲基化分析

在R（版本4.3.2）中使用Bioconductor包DSS v2.50.1处理数据。输入文件格式包含染色体名称、基因组坐标、总读段计数和甲基化读段计数的列。CpG位点按覆盖度阈值>7x（甲基化和未甲基化读段计数之和）过滤，仅保留所有样本中均存在的常见位点进行分析。最终共纳入110,450个CpG位点。

使用makeBSseqData函数创建BSseq对象，并分配两个组标签（阴性/阳性）以定义实验条件。

实施DMLtest函数并设置平滑选项=‘True’以估计平均甲基化水平。获得的结果进一步处理以识别差异甲基化胞嘧啶/区域（DMC/DMRs）。DMC基于后验概率阈值（delta设为0.1）视为显著，随后以FDR 0.05过滤剩余位点。相同标准应用于使用callDMR函数识别DMRs，该函数识别统计显著的CpG区域，并在至少3个CpGs位于50个碱基对内时合并位置。以阴性动物作为对照组，阳性动物作为乳腺炎组，将DMCs和DMRs分类为低甲基化或高甲基化特征。在此背景下，当对照组基因组中的甲基化频率较低时，我们称之为低甲基化DMC/DMR；而当该位置在对照组中的甲基化高于乳腺炎组时，则为高甲基化。

注释

使用R包ChIPseeker v1.40.0中的annotatePeak函数在默认模式下将差异甲基化胞嘧啶或区域链接到最近基因。注释考虑转录起始位点（TSS）周围±3 kb的区域。根据相对于最近基因的位置，每个注释位点或区域被分类为‘启动子’、‘外显子’、‘内含子’或‘远端基因间’。这些分析使用的参考基因组为Bubalus bubalis NDDB_SH_1，补充基因组注释由相关GFF3文件（GenBank组装登录号：GCF_019923935.1）提供，该文件获取自NCBI基因组数据库。

结果

测序结果

测序结果如表1所示。使用NanoPlot（v1.36.2）进行样本质量评估。样本的平均读长约为6 kb，平均质量（Q分数）均>13，符合牛津纳米孔测序指标。保留测序深度大于4X的样本（n=15），排除1个（乳腺炎）未通过质量检查的样本。表1还报告了其他测序指标，包括总读段数、生成的总碱基数、测序深度以及两组中在≥7×覆盖度下检测到的修饰碱基数量。

所有分析样本的CpG位点甲基化百分比分布呈单峰模式，大多数值集中在较高甲基化水平（80–100%）。CpG位点分布表明大多数CpG位点表现出高甲基化水平，与乳腺炎状态无关。

阳性和阴性组内的甲基化谱呈现Pearson相关系数，范围在0.91至0.95之间。每个样本的甲基化频率分布一致且呈左偏态。相比之下，组间（阳性和阴性样本之间）相关性显示略低的相关系数值，范围在0.86至0.89。

差异甲基化分析

使用delta 0.1（代表两组估计平均甲基化水平之间的最小差异）共识别出22个差异甲基化胞嘧啶（DMCs）。负delta值对应对照组相比患乳腺炎动物甲基化较低，而正delta值表示对照组甲基化较高。在识别的DMCs中，68%在对照组表现出低甲基化模式，而32%显示高甲基化。就DMRs而言，仅检测到一个低甲基化区域，位于8号染色体上。

这些发现通过小提琴图（图1）展示，该图显示了低甲基化和高甲基化位点之间甲基化差异的分布。每个小提琴图内部的中央箱线图显示甲基化差异的四分位距和中位数，而周围的形状勾勒出分布模式，突出显示甲基化差异更集中的区域。小提琴图的较宽部分指示具有相似甲基化差异的DMCs浓度较高。

基因组注释揭示了高甲基化和低甲基化胞嘧啶的独特模式。在高甲基化DMCs中，15%位于启动子区域，57%位于内含子区域，其余部分位于远端基因间区域。相比之下，低甲基化DMCs绝大多数富集在远端基因间区域（93%），仅7%位于第一个内含子中（图2a）。

还评估了差异甲基化胞嘧啶（DMCs）与最近转录起始位点（TSS）的距离。高甲基化DMCs主要富集在距离TSS 10–100 kb处，较小部分位于TSS 3–10 kb范围内，表明其对远离启动子近端区域的调控元件有偏好。低甲基化DMCs也显示主要富集在距离TSS 10–100 kb处。然而，显著一部分（约20%）出现在距离大于100 kb处，可能反映了远程增强子活性或其他远端调控元件（图2b）。

讨论

乳腺炎是影响全球乳业最重要的疾病之一，具有重大的经济和动物福利影响。据我们所知，这是首个使用纳米孔测序技术研究水牛乳腺炎中DNA甲基化作用的研究。

我们在水牛对照组和乳腺炎组之间进行了比较分析以探索差异甲基化模式。先前报告表明，水牛中SCC >200,000 cells/mL会对产奶量和乳糖百分比产生负面影响。此外，欧盟指令为用于生乳制品的牛奶设定了400,000 cells/mL的SCC阈值。基于这些指南，我们的研究将呈现SCC≥400,000 cells/mL、DSCC≥60%且细菌学检测阳性的水牛分类为乳腺炎。相反，对照组具有SCC<200,000 cells/mL、DSCC<60%且细菌学结果阴性。

检测甲基化标记的传统方法通常涉及亚硫酸氢盐测序，因其能够以单碱基分辨率精确绘制DNA甲基化图。然而，这些方法存在某些局限性，包括需要PCR扩增、链合成和化学处理，这可能引入偏差、降解DNA并导致甲基化检测不完全。相比之下，纳米孔测序具有若干显著优势，使其成为研究DNA甲基化的强大工具，例如即使在低测序深度下也能在同一测序过程中同时进行基因分型和表观基因分型，且无需PCR扩增或化学转化，从而保留天然DNA结构并实现甲基化标记的直接检测。实时和单分子分辨率检测此类修饰的能力提供了更全面的表观基因组视图。鉴于这些优势，我们为本研究选择了纳米孔测序。

我们的结果与哺乳动物系统中更广泛的研究一致，其中70–80%的CpG位点通常被甲基化。甲基化百分比的左偏（负偏）分布表明分析的组织整体高甲基化，这是分化或静止组织的典型特征，其中大部分基因组转录沉默。有趣的是，这种模式不同于通常报道的DNA甲基化谱中的典型双峰分布。

我们观察到DNA甲基化模式主要拦截基因间和内含子区域，与先前研究一致。已经证明，基因间区域的甲基化通过替代起始位点或microRNA表达参与转录调控，有助于基因组稳定性。此外，第一个内含子——一个通常富含增强子的基因区域——的甲基化与基因表达调控密切相关。因此，这些区域的表观遗传修饰可能作为动物生理状态的有价值生物标志物，因为它们影响超出常见启动子相关研究的基因活性。

虽然整体模式在样本间似乎一致，但这种全球趋势可能表明存在位点特异性变异，这些变异仍可能促进乳腺炎的发展。特别是，我们的发现与先前研究一致，表明平均甲基化频率在对照组中高于乳腺炎组。这可能反映了在亚临床乳腺炎发展过程中潜在参与宿主反应的特定基因的上调。对照组表现出更高比例的低甲基化DMCs，主要位于距离转录起始位点（TSS）更远处，这种模式可能表明对照组与患乳腺炎动物相比倾向于基因组稳定性。乳腺炎与TSS附近DNA甲基化模式的改变相关，可能影响基因表达。我们识别出八个差异甲基化基因（DMGs），大多数富集在远端区域（4）、内含子区域（3），仅一个在启动子区域。预测的基因（LOC102394483、LOC123331926、LOC123331980、LOC112582810、LOC123331776、LOC112581131、LOC123333069和LOC123331969）显示出不同的功能作用。前四个涉及在乳腺炎期间平衡免疫反应、管理细胞应激、促进组织修复和消退炎症。其余基因仍未表征，目前尚无可靠信息。

LOC102394483编码一种泛素羧基末端水解酶17样蛋白，是泛素-蛋白酶体系统（UPS）中的关键酶。UPS对于降解氧化或错误折叠的蛋白质至关重要，在细胞应激反应和免疫监视中发挥关键作用。该基因的活性支持更有效的蛋白质周转，在非疾病条件下减少细胞应激并促进稳态。值得注意的是，该基因在对照组中被发现低甲基化。

LOC123331926是一种内源性逆转录病毒衍生的假基因，与病毒包膜多蛋白（特别是PABLB组）具有序列同源性。内源性逆转录病毒（ERV）是古代病毒感染残留物，在进化过程中已整合到宿主基因组中。尽管通常无功能，但假设某些ERV序列在调节基因表达和调节免疫反应中保留作用。LOC123331926在对照组的内含子区域被发现高甲基化。这种增加的甲基化可能反映LOC123331926的转录抑制，作为一种保护机制以限制ERV激活并降低免疫失调风险。相比之下，在乳腺炎组中观察到的低甲基化可能与ERV表达增加有关，假设其有助于慢性炎症或免疫过度激活。然而，需要进一步研究以阐明这些观察的功能相关性。

人类中已描述了类似机制，其中ERV被认为通过涉及先天和适应性免疫的多种途径触发炎症过程。

LOC123331980编码一种黑色素瘤相关抗原B4样蛋白。MAGE家族已知调节关键细胞过程，包括自噬、凋亡、DNA修复和mRNA稳定性，同时在应激反应中发挥关键作用。尽管该基因的激活可能最初增强细胞防御，但 prolonged 或 excessive 活性可能有助于组织损伤和慢性炎症。这一假设与MAGE蛋白的其他研究结果一致，这些蛋白在应激条件下经常上调并在病理状态期间发挥细胞重塑作用。

LOC112582810编码具有磷酸肌醇磷脂酶C X（PI-PLC X）结构域的蛋白质，通常存在于参与细胞生理学和信号转导的酶中。PI-PLC同工酶通过调节钙信号和其他对细胞增殖和分化重要的细胞内途径在细胞代谢中发挥重要作用。在我们的研究中，与患乳腺炎个体相比，LOC112582810在对照组中被发现高甲基化。

本研究使用了来自单个农场的16个样本。虽然这控制了环境影响，但小样本量可能限制统计功效和结果的稳健性。为验证这些结果并获得更准确的甲基化模式理解，增加分析动物数量至关重要。

此外，虽然识别了差异甲基化基因，但未评估它们对基因表达的直接影响。由于大多数这些基因位于基因间区域，通过RNA-seq或类似方法研究它们的转录活性至关重要。这将有助于确定这些表观遗传修饰的功能相关性及其对近端调控区域的潜在影响。

此外，Bubalus bubalis基因组尚未完全表征和注释，限制了关于基因组标记和特定生物学通路的详细数据的可用性。在此背景下，本研究提供了有价值的新见解，有助于扩展该物种的基因组知识并为未来研究奠定基础。

结论

本研究首次利用纳米孔测序探索意大利地中海河流水牛乳腺炎相关的DNA甲基化变化。这项强大技术能够同时进行基因分型和表观基因分型，使其对精准育种和理解表观遗传对经济重要性状的贡献具有价值。

我们识别出对照组和患乳腺炎动物之间不同的DNA甲基化模式。我们的发现提供了对乳腺炎表观遗传调控的宝贵见解，并可能突出与该疾病相关的潜在遗传机制。这些知识可以支持水牛育种者实施疾病管理计划，最终增强早期检测和预防策略。

本研究基于少量动物，降低了其统计功效。因此，这些发现是探索性的，未来研究将旨在增加样本量以确认这些结果并确保其可靠性。

未来研究还应调查血细胞和牛奶体细胞中DNA甲基化模式之间的相关性，以评估血液分析作为乳腺炎相关表观遗传变化的间接标志物的潜力。这种方法可以促进开发用于早期乳腺炎检测的新诊断工具。

此外，将转录组学和代谢组学数据与DNA甲基化谱整合可能提供对乳腺炎过程背后分子机制更全面的理解。将表观遗传见解融入育种和管理实践有潜力改善水牛健康并提高乳业生产效率。

伦理批准

本研究经比萨大学机构动物护理和使用委员会批准（编号：11/2021，日期：2021年3月19日）。所有样本收集均遵守欧洲规则（理事会条例[EC] No. 1/2005和理事会条例[EC] No. 1099/2009）。作者确认他们遵循了欧盟用于科学目的动物保护标准。

作者贡献

SA和OGR：概念化；OGR和ALC：方法论；TG和ED：采样；IC和TG：调查；IC和ALC数据管理；IC：可视化；SA和OGR：验证；ED和FC：监督；GC、AS、TG、VP：项目管理和资金获取；IC和ALC：撰写初稿；SA、FC、MB、OGR：审阅和编辑。所有作者均已阅读并批准最终稿件。

数据可用性声明

支持本研究结果的数据可根据相应作者要求提供。