【前人研究进展】在鹅肥肝研究领域，DNA甲基化技术被广泛应用于深层挖掘畜禽脂质代谢等重要性状分子的调控机制中。其中，任学文等[3]系统研究了脂肪酸代谢的关键特征，利用气相色谱等方法发现脂肪酸含量随日龄增长呈现一定的动态变化趋势。王宝维等[4]综述了鹅肥肝生产研究现状与发展新思路，发现现有工作主要围绕少数脂代谢关键基因（如FAS、FABP）的转录表达规律展开。孙铭阳[5]通过组蛋白乳酸化修饰靶向调控FASN干预非酒精性脂肪肝病的机制研究，发现干预肝脏脂代谢和糖代谢途径的机制，为本研究提供了思路。Wanthanee等[6]发现营养提取物可通过信号通路缓解肝脏脂肪变性。杨柏高等[7]利用单细胞全基因组甲基化（Single-cell whole genome bisulfite sequencing, scWGBS）技术探究了基因组甲基化水平的变化规律，基因组甲基化水平可能通过扰乱胚胎发育、ATP代谢、MAPK信号通路、Wnt信号通路等相关基因的甲基化水平。
【本研究切入点】上述研究虽为营养诱导的表观遗传调控提供了重要线索，但鹅肥肝形成过程中DNA甲基化的全基因组动态图谱及其与脂质代谢通路的关联仍缺乏系统解析。近年来，随着脂质组学与表观遗传学的交叉融合，研究人员发现DNA甲基化等表观遗传修饰在脂质代谢调控中发挥着关键作用[8]。肝脏作为脂质代谢的核心器官，通过多维度调控网络维持游离脂肪酸（Free fatty acids, FFA）、甘油三酯（Triglycerides, TG）、胆固醇及胆汁酸（Bile acids, BA）的动态平衡[9]。特别是在肝脏脂质代谢研究中，DNA甲基化通过调控脂质合成酶基因（如FAS、SCD1）和转录因子（如PPARγ）的表达，影响甘油三酯合成和极低密度脂蛋白分泌的平衡[10-11]。相比传统转录组学，DNA甲基化研究具有可反映长期代谢记忆、稳定性高且可逆调控等独特优势，使其成为解析鹅肥肝这种生理性脂肪沉积过程的理想工具[12-13]。
【拟解决问题】因此，本研究以全基因组甲基化测序（Whole genome bisulfite sequencing, WGBS）技术为切入点，系统分析3个不同填饲阶段（前期、中期、后期）朗德鹅肝脏DNA甲基化变化规律，筛选出与鹅肥肝形成相关的关键通路和差异甲基化基因，以期为解析鹅肥肝形成的表观遗传分子机制提供理论依据。
1 材料与方法
1.1 伦理声明
本试验动物处理严格遵循《实验动物伦理审查指南》，该研究获得了青岛农业大学动物科学技术学院动物管理伦理委员会的批准（审批编号为DKY20200517）。
1.2 材料
试验时间为2020年7月至2021年7月，选取由山东潍坊春冠食品有限公司提供的同批次、体况相近的70日龄朗德鹅20只，于相同环境条件下预饲7 d。正式填饲开始后，分别在填饲前期（7 d）、中期（16 d）和后期（25 d）各随机选取3只鹅（分别为前期-1、2、3；中期-1、2、3；后期-1、2、3）进行屠宰，在肝大叶的相同位置采集组织样本，样品经蒸馏水洗净后，液氮速冻，并于-80 ℃保存，用于DNA甲基化分析。填饲饲料为98%玉米、1%水与1%植物油组成的混合型饲料，填饲操作在春冠食品有限公司肥肝生产厂内完成。
1.3 方法
1.3.1 组织DNA的提取、检测
采用DNA提取试剂盒（康为世纪生物科技股份有限公司）从肝大叶肝尖组织提取基因组DNA，使用分光光度计和Qubit 3.0检测浓度与纯度（OD260/280≥1.9）。
1.3.2 BS文库构建及全基因组甲基化重测序
将DNA经过片段化、末端修复、加A碱基、连接甲基化测序接头后，使用亚硫酸氢盐（Bisulfite）试剂盒（上海碧云天生物技术股份有限公司）进行转化处理后，脱盐处理，DNA序列进行PCR扩增，扩增后进行测序文库片段大小选择，通过Agilent2100库检参数合格后，构建文库。采用最新的Illumina NovaSeq测序平台进行测序，每个个体测序深度不小于30×基因组覆盖度，以保证甲基化位点分型的准确性，共计产生约276.5 Gb的基因组数据。
1.3.3 数据处理
测序得到的raw reads过滤为高质量的clean reads，需过滤含N比例大于10%的reads；低质量的reads（质量值Q<20的碱基数占整条reads的40%以上）。使用Bismark软件将clean reads比对至鹅参考基因组，计算CG、CHG、CHH上下文甲基化水平，并统计基因功能区域（启动子、上游2 kb、下游2 kb、外显子、内含子等）甲基化分布。采用主成分分析评估样本间甲基化模式差异。
1.3.4 差异甲基化分析
使用methyKit软件鉴定差异甲基化区域（Differential methylated region, DMR），标准为：至少5个C位点、覆盖深度>5×、区域长度≤1000 bp、相邻区域距离≥200 bp，差异阈值P≤0.001且甲基化水平差异≥20%。使用ANNOVAR对DMR进行基因注释，基于比对结果进行C碱基有效深度和覆盖度分析、全基因组甲基化水平趋势分析，并绘制全基因组甲基化图谱。
1.3.5 差异甲基化GO、KEGG功能富集分析
对DMR相关基因进行了GO富集分析和KEGG通路富集分析。若DMR区域中点与基因或基因的Promoter区域（TSS上游2 kb至下游2 kb）有交集，则认为该DMR与此基因有关联，对DMR关联基因的GO和KEGG注释的基因按照P-value＜0.05分别进行GO和KEGG富集分析。
2 结果与分析
2.1 不同填饲时期鹅肥肝DNA甲基化测序数据统计分析
表1 鹅肥肝DNA甲基化重测序数据产出情况
Table 1
样本 Samples | 有效reads Effective reads | 比对reads Comparison reads | C/T转化率 C/T conversion rate | 比对率 Comparison rate/% | Q30/% |
|---|---|---|---|---|---|
| 前期-1 Early-1 | 29130791700 | 22663755943 | 0.994884552 | 77.80 | 90.63 |
| 前期-2 Early-2 | 29253349800 | 22583586046 | 0.994809201 | 77.20 | 90.91 |
| 前期-3 Early-3 | 29070828300 | 21512412942 | 0.994877554 | 74.00 | 86.96 |
| 中期-1 Middle-1 | 28905671100 | 21158951245 | 0.995017668 | 73.20 | 86.09 |
| 中期-2 Middle-2 | 31352284200 | 22824462898 | 0.995071942 | 72.80 | 86.99 |
| 中期-3 Middle-3 | 29842813800 | 21755411260 | 0.995101413 | 72.90 | 86.94 |
| 后期-1 Final-1 | 29220756300 | 20805178486 | 0.995159206 | 71.20 | 86.23 |
| 后期-2 Final-2 | 30281492100 | 21772392820 | 0.995517081 | 71.90 | 87.14 |
| 后期-3 Final-3 | 29216945100 | 20890115747 | 0.995057107 | 71.50 | 86.22 |
表2 不同分布类型甲基化C的组成比例及其甲基化水平
Table 2
样本 Samples | mC | mC甲基化水平mC Methylation level | mCpG/% | mCpG甲基化水平mCpG Methylation level | mCHG/% | mCHG甲基化水平mCHGMethylation level | mCHH/% | mCHH甲基化水平mCHHMethylation level |
|---|---|---|---|---|---|---|---|---|
| 前期-1 Early-1 | 4210763360 | 3.96 | 64.947 | 19.79 | 0.707 | 4.15 | 0.719 | 1.41 |
| 前期-2 Early-2 | 4266739976 | 4.08 | 65.608 | 19.36 | 0.717 | 4.14 | 0.728 | 1.41 |
| 前期-3 Early-3 | 3907388946 | 4.0 | 66.683 | 20.25 | 0.705 | 4.15 | 0.724 | 1.41 |
| 中期-1 Middle-1 | 3869687066 | 3.82 | 64.957 | 20.59 | 0.691 | 4.17 | 0.705 | 1.41 |
| 中期-2 Middle-2 | 4201846089 | 3.90 | 65.146 | 20.09 | 0.678 | 4.16 | 0.690 | 1.41 |
| 中期-3 Middle-3 | 4074214208 | 3.82 | 62.979 | 19.90 | 0.678 | 4.17 | 0.689 | 1.41 |
| 后期-1 Final-1 | 3871948203 | 3.86 | 64.917 | 20.24 | 0.675 | 4.19 | 0.685 | 1.40 |
| 后期-2 Final-2 | 4044916377 | 3.79 | 64.593 | 20.44 | 0.646 | 4.18 | 0.662 | 1.40 |
| 后期-3 Final-3 | 3873964105 | 3.87 | 65.082 | 20.30 | 0.691 | 4.18 | 0.707 | 1.41 |
2.2 甲基化位点在不同基因功能区域的分布特征
为了解各组样本之间总体甲基化差异和样本间变异度大小，对样本进行主成分分析，3个组样本组内间距小，表明组内差异不大，3个组能明显区分开，表明组间差异明显（图1a）。
通过对鉴定到的甲基化位点进行注释和统计分析，发现编码基因的不同区域往往受到不同程度的DNA甲基化修饰。因此对基因的不同区域包括外显子、内含子、编码序列，5'非翻译区，3'非翻译区，基因上游2 kb，基因下游2 kb进行C位点序列的平均甲基化水平统计。如图1b所示，发现在全基因组范围内外显子区的总体甲基化水平最高，内含子区的总体甲基化水平不高，5'非翻译区总体甲基化水平略高于3'非翻译区，而启动子区域总体甲基化水平最低。
2.3 不同填饲阶段鹅肥肝差异DNA甲基化分析
将填饲前期与填饲中期、填饲后期分别进行比较，分别检测出951、738和778个DMR，如图2所示，3个不同填饲时期的差异甲基化都有共同特征，即差异甲基化基因大部分位于染色体1到染色体4上。对其进行分析发现，在填饲前期与填饲中期比较组中，有524个DMR甲基化水平上升，有427个DMR甲基化水平下降；在填饲中期与填饲后期比较组中，有432个DMR甲基化水平上升，有236个DMR甲基化水平下降。在3个组中，甲基化水平下降的DMR数量均多于甲基化水平上升的DMR数量。
图1
图1
全基因组甲基化基本特征与样本关系分析
Fig. 1
Basic characteristics of whole-genome methylation and sample relationship analysis
图2
图2
不同填饲阶段差异甲基化曼哈顿图与差异甲基化数量
Fig. 2
Manhattan plots and the number of differential methylation for different over-feeding stages
对填饲前期与填饲中期的差异甲基化与差异基因的分析结果，共发现有8个基因，其中高甲基化且表达下调的基因有5个，分别为SLC6A15、Mc2r、MMRN1、DYNLRB2、TENM2。低甲基化且表达上调的有3个，分别为FRMPD2、DHX32、SORCS1。详细信息见表3。
表3 各时期甲基化水平与基因表达水平负相关基因
Table 3
| 时期Time | 基因名Name of gene | 基因功能Function of genes | 甲基化Methylation of DNA | 表达量Quantity of expression |
|---|---|---|---|---|
前期与中期 Early and middle | SLC6A15 | solute carrier family 6 members 15 | hyper | down |
| Mc2r | melanocortin 2 receptor | hyper | down | |
| MMRN1 | multimerin 1 | hyper | down | |
| DYNLRB2 | dynein light chain roadblock-type 2 | hyper | down | |
| TENM2 | teneurin transmembrane protein 2 | hyper | down | |
| FRMPD2 | FERM and PDZ domain containing 2 | hypo | up | |
| DHX32 | DEAH-box helicase 32 | hypo | up | |
| SORCS1 | sortilin-related VPS10 domain-containing receptor 1 | hypo | up | |
中期与后期 Middle and final | ANGPTL5 | angiopoietin like 5 | hyper | down |
| DNAH10 | dynein axonemal heavy chain 10 | hypo | up | |
| CCDC187 | coiled-coil domain containing 187 | hypo | up | |
| PCK1 | phosphoenolpyruvate carboxykinase 1 | hypo | up |
同时填饲中期与填饲后期的差异甲基化与差异基因的分析，共发现有4个基因，其中高甲基化且表达下调的基因有1个，为ANGPTL5。低甲基化且表达上调的有3个，分别为DNAH10、CCDC187、PCK1。详细信息见表3。
2.4 差异甲基化基因的GO和KEGG富集分析
图3
图3
不同填饲阶段差异甲基化GO富集图
Fig. 3
GO enrichment map of differential methylation in different over-feeding periods
KEGG分析显示，前期与中期差异基因显著富集于121条通路，其中胰岛素信号通路、FoxO信号通路等与脂质沉积密切相关，涉及MAPK10、IL6等6个基因。中期与后期富集到96条通路，MAPK信号通路、FoxO信号通路等显著富集，关联基因包括MAPK10、PCK1、PLA2G4A等13个基因（图4）。
图4
图4
不同填饲阶段差异甲基化KEGG富集图
Fig. 4
Differential methylation KEGG enrichment maps across different over-feeding stages
3 讨 论
肝脏是脂质代谢的核心器官，在能量过剩时易发生脂质过度蓄积，导致人类非酒精性脂肪性肝病（Nonalcoholic fatty liver disease, NAFLD）的病理性脂肪变性[12]。然而，朗德鹅在短期强制填饲后可形成生理性肥肝，而不伴明显肝功能损伤，这一独特适应机制提示其脂代谢调控与哺乳动物NAFLD存在本质差异。表观遗传修饰尤其是DNA甲基化在代谢性疾病中发挥重要作用[13]，但其在鹅肥肝形成过程中的动态变化规律尚不清楚。本研究通过高深度WGBS（深度≥30×），系统描绘了填饲诱导鹅肥肝形成过程中的DNA甲基化动态图谱。结果显示，与人类NAFLD中常见的全局低甲基化趋势不同，鹅肥肝形成伴随特定功能区域（如启动子）的甲基化调控。这一发现揭示了生理性脂质蓄积与病理性脂肪变性在表观遗传层面的差异，为阐明水禽独特脂代谢适应机制提供了新的思考角度。
3.1 CG位点甲基化主导及启动子低甲基化在脂代谢中的作用
本研究利用WGBS测得的数据质量分析显示，测序深度≥30×，Q30均高于87%，且比对率与CT转化效率均符合高标准要求，为后续分析提供了可靠数据基础。研究发现，在鹅肝脏中，CG位点上下文（mCpG）的甲基化占据绝对主导地位（65.08±1.50）%，显著高于CHG和CHH类型。这一模式与在哺乳动物肝脏中的研究结果一致。本研究对甲基化位点的基因组分布特征分析发现，基因启动子区域呈现显著的低甲基化状态（启动子低甲基化通常与染色质开放和基因转录激活相关）。这一发现为解释填饲过程中，鹅肝脏为何能够高效激活脂质合成相关基因的表达提供了重要的结构基础，这与人类NAFLD中观察到的全基因组广泛低甲基化趋势不同[14]。这种差异表明，水禽鹅肥肝的形成并非一种病态的、失控的脂肪变性，可能是一种受精密调控的生理性适应过程。鹅在进化中形成的这种独特表观遗传调控机制，使其能够在短期内承受巨大的能量负荷，体现了鹅肥肝的精密生理性适应[15]。
3.2 关键信号通路通过甲基化调控脂质合成与转运失衡
KEGG通路富集的分析结果显示，胰岛素-FoxO通路轴和MAPK信号通路在脂质合成中发挥主要作用，填饲中期胰岛素信号通路显著富集，MAPK10基因的低甲基化状态促进其表达上调，通过PPARγ途径直接促进脂质合成，这一发现与Pawlak等[16]在人类NAFLD研究中报道的胰岛素信号通路异常激活有相似之处，但鹅肥肝中MAPK10通过PPARγ途径促进脂质合成的机制更为直接，不同于人类中观察到的胰岛素抵抗现象。还检测到填饲中期FoXO信号通路中的PCK1基因低甲基化上调表达，增强了糖异生-脂质转化效率，这一机制与Ye等[17]在小鼠NAFLD模型中的发现形成有趣对比，在小鼠中PCK1表达上调会促进疾病进展，而在鹅肥肝形成中却表现为生理性适应。
3.3 甲基化动态变化与肝重增加的时空关联性
本研究通过不同填饲阶段的DMR分析，发现填饲中期至后期DMR数量下降（738 vs 778），且低甲基化DMR占比高于50%，这与肝重增加趋势一致（填饲后期肝体积扩张3～5倍）。PCA分析结果显示了3个填饲阶段样本的组间差异显著，与Chen等[20]在鹅肌肉组织中的研究发现形成了有趣的对比，表明不同组织对营养干预的表观遗传响应存在显著差异。另外PPAR和胰岛素信号通路在填饲后期显著激活，相关基因（如PCK1、S1PR1）的低甲基化状态可能通过提高脂质合成效率直接促进肝重增加，本研究发现的S1PR1基因低甲基化与极低密度脂蛋白（Very low-density lipoprotein, VLDL）分泌关联的机制，与Pike等[21]在新生仔猪中观察到的中链甘油三酯调控机制有部分重叠，但鹅特有的VLDL分泌限制现象为水禽脂肪代谢研究提供了独特的视角。与哺乳动物NAFLD研究相比，本研究发现Wnt信号通路在鹅肥肝形成中富集程度较低，这与肖成等[22]在绵羊脂肪细胞分化中的发现形成鲜明对比，表明不同物种间脂肪沉积的调控网络存在重要差异。这种差异可能反映了草食动物与杂食动物在能量代谢适应策略上的进化分异。
4 结论
本研究利用全基因组亚硫酸氢盐测序技术，绘制了朗德鹅在填饲前期、中期及后期肝脏DNA甲基化的动态变化图谱。经分析，鹅肥肝形成过程中的DNA甲基化修饰具有显著的阶段特异性与区域选择性。全基因组范围内，CpG位点甲基化占据主导地位，而基因启动子区域呈现特征性的低甲基化状态，这为脂质合成相关基因的转录激活提供了关键的表观遗传基础。通过对差异甲基化区域的分析及功能富集，经研究分析发现，FoxO信号通路、胰岛素信号通路、PPAR信号通路及MAPK信号通路是调控鹅肥肝生理性脂质沉积的核心通路，关键候选基因MAPK10、PCK1、PLA2G4A等甲基化水平在填饲过程中发生动态变化，通过调控其表达水平，协同促进了肝脏内的糖异生-脂质转化效率与甘油三酯合成。该发现表明DNA甲基化作为一种关键的分子开关，通过多通路、多靶点的精密调控网络，驱动鹅肝脏在短期高能量负荷下的代谢重编程，从而形成独特的生理性脂肪肝。这些通路和候选基因从表观遗传学层面为解析水禽抗脂肪肝损伤的独特适应机制提供了新的视角与理论依据。
参考文献
四川麻兔背最长肌DNA甲基化和转录组联合分析
[J].
Combined analysis of DNA methylation and transcriptome in longissimus dorsi muscle of Sichuan Ma rabbits
[J].
Transcriptome analyses reveal genes and pathways associated with fatty acid composition traits in pigs
[J].
不同品系五龙鹅胸肌脂肪酸含量动态变化规律
[J].
Dynamic changes of fatty acid composition in breast muscles of Wulong geese across different strains
[J].
鹅肥肝生产研究现状与发展新思路
[J].
Research status and new development ideas of goose fatty liver production
[J].
鸢尾苷通过组蛋白乳酸化修饰靶向调控FASN干预非酒精性脂肪肝病的机制研究
[D].
Study on Tectoridin Targeting FASN Through Histone Lactylation Modification to Intervene in Non-Alcoholic Fatty Liver Disease
[D].
Sangyod rice extract mitigates insulin resistance in HepG2 cells and hepatic steatosis in diabetic rats via AMPK/mTOR/MAPK signaling pathways
[J].
基于scWGBS探究玻璃化冷冻对猪孤雌激活囊胚基因组甲基化水平的影响
[J].
Exploring the effect of vitrification on genome methylation level of porcine parthenogenetic activation blastocysts by scWGBS
[J].
电针“丰隆”“足三里”对非酒精性脂肪肝大鼠SIRT1/FOXO1信号通路的影响
[J].
Effects of electroacupuncture at "Fenglong"(ST40) and "Zusanli"(ST36) on the SIRT1/FOXO1 signaling pathway in non-alcoholic fatty liver disease model rats
[J].
丝裂原活化蛋白激酶信号通路的研究进展
[J].
Research progress of mitogen-activated protein kinase signaling pathway
[J].
Molecular cloning, characterisation, and expression analysis of adipocyte fatty acid binding protein gene in Xupu goose (Anser cygnoides domesticus)
[J].
Fatty acids modulate the expression of pyruvate kinase and arachidonate-lipoxygenase through PPARγ/CYP2C45 pathway: a link to goose fatty liver
[J].
Uncovering the embryonic development-related proteome and metabolome signatures in breast muscle and intramuscular fat of fast-and slow-growing chickens
[J].
Omics integrating physical techniques: aged Piedmontese meat analysis
[J].
表观遗传学在非酒精性脂肪性肝病中的研究进展
[J].
Advances in epigenetics research in non-alcoholic fatty liver disease
[J].
宫内发育迟缓相关疾病的表观遗传修饰调控机制研究进展
[J].
Research progress on epigenetic modification and regulation mechanism of intrauterine growth retardation-related diseases
[J].
Molecular mechanism of PPARα action and its impact on lipid metabolism, inflam-mation and fibrosis in non-alcoholic fatty liver disease
[J].
Deficiency of gluconeogenic enzyme PCK1 promotes metabolic-associated fatty liver disease through PI3K/AKT/PDGF axis activation in male mice
[J].
Fads1 and 2 are promoted to meet instant need for long-chain polyunsaturated fatty acids in goose fatty liver
[J].
High fat diet-triggered non-alcoholic fatty liver disease: a review of proposed mechanisms
[J].
皖西白鹅育肥期肌肉脂肪酸组成及肝PPARα、FADS2和ME1基因表达规律的研究
[J].
The fatty acid profile in muscles and expression of PPARα, FADS2 and ME1 genes in liver of Chinese wanxi white geese in fattening period
[J].
Effects of medium chain triglycerides on hepatic fatty acid oxidation in clofibrate-fed newborn piglets
[J].
