气候因素驱动旧大陆家牛的适应性演化
iMeta主页：http://www.imeta.science
研究论文
● 原文: iMeta(IF 33.2, 中科院双一区Top)
● 英文题目: Climate factors drive the local adaptation of Old World cattle
● 中文题目: 气候因素驱动旧大陆家牛的适应性演化
● DOI: https://doi.org/10.1002/imt2.70108
● 2026年2月9日，西北农林科技大学雷初朝、陈宁博、夏小婷和宁夏大学马云等在iMeta在线发表了题为“Climate factors drive the local adaptation of Old World cattle”的文章。
● 本研究通过整合全基因组数据与多维环境信息，系统揭示了气候因素驱动旧大陆家牛局部适应的遗传机制。环境全基因组关联分析鉴定出大量与关键气候变量相关的遗传位点，这些位点注释的基因涉及能量代谢、抗逆反应等通路，部分位点在不同牛种间呈现显著的遗传分化。结合表观调控与eQTL分析，阐明了非编码变异通过调控基因表达从而影响家牛气候适应性的分子机制。本研究为培育气候适应性牛品种提供了关键靶点与理论依据。
● 第一作者：孙露洋、金良梁、孙婷
● 通讯作者：陈宁博（ningbochen@nwafu.edu.cn）、雷初朝（leichuzhao1118@nwafu.edu.cn）、夏小婷（xiaxiaoting@nwafu.edu.cn）、马云（mayun@nxu.edu.cn）
● 合作作者：边晨启、汪富文、李双、李昊、刘建勇、张继才、Johannes A. Lenstra、Artem Nedoluzhko、郑竹清、向瑞东、Wai Yee Low、巴桑旺堆、巴桑珠扎、李斌、索朗曲吉、赵丽、斯朗旺姆、次仁罗布、普布占堆、西热强马、宋颜亮、高元鹏、魏大为、宋仁德、彭旻晟、黄必志
● 主要单位：西北农林科技大学动物科技学院、宁夏大学动物科技学院
● 对来自亚欧非旧大陆的85个家牛群体进行了系统的全基因组分析，揭示了与家牛气候适应性相关的关键遗传因子；
● 环境全基因组关联分析鉴定出3,165个与气候因子显著相关的单核苷酸多态性位点；
● 在牛第18号染色体上发现可能通过调控ANKRD11、CDK10等基因表达来影响温度适应性的启动子区变异位点；
● 研究结果为应对气候变化下肉牛业的可持续发展及通过分子标记辅助育种培育气候适应性牛品种提供了重要的遗传资源。
气候变化对动物种群构成多重挑战。作为关键农业物种，家牛应对气候变化的适应性遗传机制尚未完全阐明。本研究收集并整合了来自欧亚大陆和非洲85个家牛品种/群体的全基因组数据，以及其分布区域的28项气候因子，系统揭示了与家牛基因组变异与气候因子相关的遗传基础。通过环境全基因组关联分析，共鉴定出3,165个单核苷酸多态性位点，其等位基因频率与特定气候因子显著关联。这些SNP所关联的候选基因主要涉及脂肪沉积、能量代谢与肾脏发育等生物学过程，并在雅库特牛、藏牛和恩达马牛等适应极端环境的品种中呈现显著遗传分化。此外，我们在牛18号染色体一个受强选择的区域内，发现了与温度因子显著相关的调控位点。这些位点位于ANKRD11、RPL13、CDK10和VPS9D1等基因的启动子和其它调控区域，可能通过调控基因表达影响温度适应性。进一步整合表达数量性状位点数据集发现，位点rs47611688在普通牛群体中是一个顺式eQTL，参与调控FARP2基因的表达，其基因型与帕尔默干旱严重指数显著相关。这些发现为家牛环境适应性的分子育种提供了直接资源，可用于增强牛群对气候变化的适应力，为未来环境挑战下家牛的可持续选育提供理论与实践依据。
Bilibili：https://www.bilibili.com/video/BV1MHFDz3Ew7/
Youtube：https://youtu.be/7hO69MYYrRM
中文翻译、PPT、中/英文视频解读等扩展资料下载
请访问期刊官网：http://www.imeta.science/
引 言
气候变化构成强大的选择压力，驱动着物种的局部适应性进化，主要体现在空间遗传变异与基因型-环境的相互作用。全基因组范围的气候变量关联扫描，已成为解析非模式物种适应性遗传基础的有效方法。家牛（普通牛与瘤牛）展现出卓越的气候适应可塑性：普通牛具有明显的耐寒特性（体型较大、被毛厚密），而瘤牛则表现出显著的耐热与抗寄生虫特性。历史上普通牛与瘤牛的杂交，则衍生出众多适应不同地域环境的家牛品种/群体。尽管先前基于基因组学的研究已通过关联信号识别出若干候选基因，但由于研究对象多集中于高产品种或特定地理区域，且缺乏整合的多组学数据，我们对牛种全球尺度气候适应的整体理解仍显不足，尤其涉及非编码区单核苷酸多态性（SNP）的调控机制尚未明晰。本研究聚焦于旧大陆（欧亚大陆与非洲）长期适应多样化环境条件的家牛群体，为解析家牛气候适应的遗传基础提供了理想的研究对象。
结果与讨论
旧大陆家牛的全基因组遗传特征
本研究共收集了来自旧大陆85个牛品种/群体的744个个体的基因组数据（表S1、S2），并将其划分为10个亚群：非洲瘤牛（AFI, n= 107）、普遍非洲普通牛（AFT, n= 17）、欧洲普通牛（EUT, n =72）、东亚普通牛（EAT, n= 47）、藏牛（TB, n= 64）、欧亚普通牛（EUAT, n=39）、南亚瘤牛（SAI, n= 124）、东南亚瘤牛（SEAI, n= 41）、华中地区牛（CC, n= 99）及东亚瘤牛（EAI, n= 134）（图1A）。通过变异检测，共识别出74,041,479个常染色体SNPs。全基因组核苷酸多样性分析显示，瘤牛群体的遗传多样性高于普通牛群体（图S1），这可能与普通牛系长期经历高强度商业化选育有关。群体间遗传分化分析显示，普通牛与瘤牛品种间的成对固定指数（FST）最高可达0.502（图S2及表S3）。群体结构、主成分分析与系统发育网络的分析结果均表明，除TB与CC存在混合遗传背景外，其余8个亚群呈现出清晰的遗传分化（图1B–E、图S3及图S4A）。主成分分析进一步揭示了群体遗传结构的地理分化模式：第一主成分（PC1）区分普通牛与瘤牛；第二主成分（PC2）区分南亚瘤牛与东亚瘤牛（图1C）；第三主成分（PC3）则区分欧洲普通牛与亚洲普通牛（包括东亚普通牛与藏牛）（图1D）。群体分析同时发现，东亚瘤牛与南亚瘤牛群体中存在与普通牛的杂交信号。剔除杂交个体后，群体核苷酸多样性出现小幅的下降（图S4B）。为进一步解析欧亚及非洲牛品种间的亲缘关系，我们基于FST值构建了邻接网络图（NeighborNet），直观展示了不同牛种间的遗传联系与分化格局（图1E、图S5及表S4）。
图1. 牛群体基因组学图谱、环境变量与候选SNP分布
(A) 全球牛品种/群体地理分布示意图。图中展示了85个牛品种/群体的采样位置，包括：非洲瘤牛（AFI）、非洲普通牛（AFT）、欧洲普通牛（EUT）、东亚普通牛（EAT）、藏牛（TB）、欧亚普通牛（EUAT）、南亚瘤牛（SAI）、东南亚瘤牛（SEAI）、华中地区牛（CC）及东亚瘤牛（EAI）。本图叠加了 Köppen–Geiger 气候分类地图（数据来源于国家地球系统科学数据中心， http://www.geodata.cn）。气候类型图例：A，热带气候；f，热带雨林气候；m，季风气候；w，热带草原气候。B，干旱气候；w，沙漠气候；s，草原气候。C，温带气候；s，夏干型；w，冬干型；f，无旱季型。D，寒带气候；s，夏干型；w，冬干型；f，无旱季型。E，极地气候；T，苔原气候；F，冰原气候。(B) 牛群体遗传混合分析结果，展示了K = 2和K = 7（最优K值，交叉验证误差 = 0.245）时的祖先成分比例。(C) 主成分分析（PCA）显示全部744个个体的PC1与PC2分布，以及(D) PC1与PC3的分布。颜色代表样本的地理区域归属。(E) 系统发育网络图，使用 SplitsTree4 构建，展示了各品种间的遗传关系。(F) 28个气候因素年均值的斯皮尔曼等级相关系数绝对值热图。28个气候因素及其缩写如下：地表温度（STR）、向外长波辐射（OLR）、位温（POTTMP）、气温（TMN）、最高温度（TMAX）、最低温度（TMIN）、帕尔默干旱严重指数（PDSI）、最大比湿（QMAX）、比湿（SHUM）、最小比湿（QMIN）、相对湿度（RHUM）、可降水量（PRWTR）、降水率（PRATE）、潜在蒸发率（PEVPR）、风速（WSPD）、U风分量（UWND）、地表气压（SRFP）、水汽压（VAP）、潜在蒸散发（PET）、地面霜冻频率（FRS）、土壤湿度（SOILW）、V风分量（VWND）、昼夜温度范围（DTR）、湿润日频率（WET）、云量（CLD）、雪盖（SNOW）、总降水量（PRE）和植被覆盖度（VEG）。文中重点讨论的特定相关性已在图中高亮标出：强相关性（TMIN 与 SHUM）以橙色标记，非显著相关性（QMIN 与 POTTMP）以蓝色标记。(G) 候选SNP在不同环境变量间的比较分布图。针对28个变量中的每一个，四组柱状图分别表示与其年均值、最小值、最大值及变化幅度（极差，即最大值减最小值）显著关联的SNP数量。
气候变量的空间异质性与复杂关联性
环境变量之间存在复杂的交互关系（附表S5、S6）。其年平均月分布特征揭示了不同牛群体所处气候模式的显著差异（图S6–S10）。对环境变量进行主成分分析（PCA）的结果进一步证实了牛群体间环境因子的显著分化（图S11）。地理位置上邻近的群体通常处于相似的环境之中，但也存在值得注意的例外。例如，南亚瘤牛各品种虽地理邻近，却存在跨度极大的生态位。与之相反，滇中牛（中国云南）和Abergelle牛（非洲）尽管地理距离遥远，却生活在极为相似的环境生态位中。气候因子之间的关联呈现复杂格局：最低温度（TMIN）与比湿度（SHUM）之间存在强相关性，而最小比湿（QMIN）与位温（POTTMP）之间则未呈现显著关联（图1F），这一对比凸显了气候因素间相互作用的复杂性。
气候适应性基因组变异的鉴定
我们联合运用了GEMMA、LFMM和BayPass三种全基因组关联分析方法，系统评估了每个SNP位点的等位基因频率与环境气候变量之间的关联（图S12，13）。最终鉴定出3165个与28个气候变量相关联的显著SNPs（表S7），其中约50.81%位于内含子区域，外显子区域的SNP数量相对较少（占比1.39%）（图S13B）。在所有气候变量中，与雪盖（SNOW）和地面霜冻频率（FRS）相关的SNP数量最多，总计达1269个（占总数的40.1%）。相比之下，湿润日频率（WET）和云量（CLD）各自仅检测到一个显著关联的SNP（WET：p = 6.79E-10，BF = 33.9；CLD：p = 2.22E-11，BF = 32.2；见图1G）。功能富集分析结果表明，转录调控、细菌反应、轴突导向、谷氨酸能突触以及cAMP/Rap1信号通路等生物学过程与环境适应密切相关（图S13C，表S8）。
温度适应性遗传变异的群体地理分布
温度是直接决定牛群热舒适性的关键气象指标。通过对牛群全基因组SNP与8个温度相关气候变量进行全基因组关联分析，我们共鉴定出1,584个显著关联的SNP，从中筛选到81个与温度适应高度相关的SNP，属于北极圈内雅库特牛和高原地区藏牛特有的品种特异性位点。在雅库特牛中，检测到70个特异性SNP与位温、雪盖及地表温度显著相关，涉及ABLIM1、TCF7L2等30个基因（图2A，图S14）。其中，ABLIM1参与调控肌肉发育，TCF7L2则是调控脂肪细胞代谢与适应性脂肪沉积的关键因子，其在PPARγ/lncRNA通路中具有抑制产热、维持脂肪可塑性的功能。在藏牛中，我们鉴定出11个与雪盖、霜冻频率、位温及地表温度相关的品种特异性SNP，其中包括位于甲状腺激素合成关键基因DUOX1中的变异位点——该基因此前已被报道在藏牛中受到强烈自然选择。此外，我们还发现PTK2与PDE10A两个参与血管形态发生基因中的新关联SNP。其中PTK2已在藏羊和藏鸡中被证实受高海拔环境的选择压力。
在相似环境压力下，处于地理隔离的种群可能独立演化出相似的有利等位基因。我们的多群体选择扫描分析揭示了此类趋同进化的信号。在43个与霜冻频率和雪盖显著相关的SNP中（p ＜ 1.4E-10，BF ＞ 30），东亚普通牛、欧洲普通牛与欧亚普通牛在18号染色体约200 kb区域（14.4–14.6 Mb）表现出强烈选择信号，而非洲普通牛中未见此信号（图2B，S15A）。单倍型结构在冷适应与热适应牛群间呈现明显分化（图2B），其中SNP rs460983591在群体间的等位基因频率差异最为显著（图2C）。更多近期选择的证据来源于冷适应群体中连锁不平衡衰减较慢（图S15B），以及错义突变rs472613982的扩展单倍型纯合性衰减显著减缓（图2D）。冷适应牛群呈现多个以祖先单倍型为主的连锁不平衡区块，而热适应牛群则表现出丰富的衍生单倍型（图2E），这表明寒冷环境中的选择压力倾向于维持祖先单倍型，而高温异质性则有利于衍生单倍型的多样化。该区域内包含与表皮生长因子信号（CHMP1A）、热敏感性（CDK10）、脂质过氧化（DPEP1、CHMP1A）及核糖体功能（RPL13）相关的基因。为探究该区域关联SNP的潜在调控机制，我们对欧洲普通牛（海福特牛）和东亚瘤牛（涠洲牛）的关键组织进行了RNA测序与转座酶可及染色质测序分析。研究发现，在ANKRD11、RPL13、CDK10和VPS9D1基因的启动子区域，两类牛群中检测到的染色质开放峰存在显著差异（FDR ＜ 0.01；图2F；表S9，10）。转录组分析显示这些基因表达量下降（p ＜ 0.05），而染色质免疫沉淀测序则检测到H3K27me3修饰的富集。值得注意的是，高表达RPL13可促进紫外线损伤修复，CDK10则具有温度敏感性。综合上述结果可知，位于非编码区域的显著关联与受选择SNP可能通过改变调控元件来影响基因表达，进而发挥其生物学功能。
图2. 基因内与适应性相关的区域分析
(A) 雅库特牛（蓝色）与藏牛（红色）适应性基因中群体特异性替换位点的频率分布。(B) 牛18号染色体上霜冻频率的关联分析结果及候选区域的单倍型结构。(C) 等位基因rs460983591在全球牛群中的地理分布频率。地图颜色梯度代表不同区域的地表温度差异，参考等位基因以蓝色表示，衍生等位基因以红色表示。(D) 在不同环境温度下生存的牛群中，rs472613982位点周围扩展单倍型纯合性的衰减情况。低温适应组包括东亚普通牛、欧亚普通牛与欧洲普通牛；高温适应组包括南亚瘤牛、东南亚瘤牛、东亚瘤牛、非洲瘤牛及非洲普通牛。(E) 不同环境温度下生存牛群的连锁不平衡区块结构。(F) 通过集成基因组学查看器绘制的图表，展示了18号染色体候选功能变异区间的调控效应。图中显示了通过RNA测序获取的海福特牛与涠洲牛骨骼肌中基因表达水平，并经经验贝叶斯校正t检验判定显著性。(G) 恩达马牛适应性基因中谱系特异性替换位点的频率分布。(H) 牛5号染色体上降水的关联分析结果。图中绘制了基于群体θπ的选择清除信号，并标注了显著区域的单倍型结构。其他群体包括东亚普通牛、欧亚普通牛、欧洲普通牛、南亚瘤牛、东南亚瘤牛、东亚瘤牛及非洲瘤牛。该区域的基因注释标注于底部。(I) NR4A1基因中rs466258385等位基因在全球牛群中的地理分布频率。rs466258385位点在非洲普通牛、非洲瘤牛与欧洲普通牛中的扩展单倍型纯合性衰减情况。地图颜色梯度代表不同区域的降水变化，参照等位基因以橙色表示，突变等位基因以蓝色表示。(J) 通过ATAC-seq技术分析的海福特牛、婆罗门牛与鲁西牛骨骼肌和肝脏中FARP2基因rs476116768位点的染色质可及性图谱。
湿度适应性遗传变异的群体地理分布
湿度是影响牛呼吸效率、皮肤屏障功能及病原暴露风险的关键环境因素。通过与13个湿度相关气候变量进行全基因组关联分析，共鉴定出870个显著关联的SNP位点。本研究在恩达马牛中发现了115个与总降水量（PRE）、潜在蒸发率（PEVPR）及降水率（PRATE）显著相关的候选特异性SNP（图2G与图S16）。这些位点的等位基因频率具有明显的品种特异性，可能与该品种特有的适应机制相关。其中，23个SNP位于CDH17基因附近，该基因在肠上皮细胞中表达，与肾脏发育密切相关。CDH17基因中的错义突变rs111440168在恩达马牛中表现出连锁单倍型特征及较高的扩展单倍型纯合性（图S17）。在3号染色体上，SLAMF1基因内的3个候选SNP被鉴定与降水率显著相关。SLAMF1中的两个错义突变仅存在于普通牛品种中，其等位基因频率在多雨地区显著升高（图S18），而藏牛群体中未见此现象。此外，在瘤牛中发现ARFGEF3基因（亦称BIG3）中3个发生突变的SNP位点（图S19）。该基因编码的Arf-GTP交换因子通过定位于高尔基体网络负调控胰岛素颗粒的生物合成与分泌。ARFGEF3功能缺陷可加速胰岛素原加工过程，增加胰岛素分泌，并导致葡萄糖稳态失调。
通过多群体选择扫描分析，在5号染色体26.8–27.8 Mb区域发现了趋同适应的遗传特征（图2H）。在非洲普通牛中鉴定出10个与降水量显著关联的非重叠SNP，这些位点涉及肺泡发育与免疫调节相关基因。核苷酸多样性分析显示该区域在非洲普通牛中存在强烈的选择信号和固定的单倍型特征。NR4A1基因中的错义突变rs466258385在非洲普通牛中衍生等位基因频率高达0.91（图2I），携带该衍生等位基因的单倍型表现出比祖先等位基因单倍型更高的扩展单倍型纯合性。该位点与降水变量的关联性，与区域内基因在肺功能调控及免疫调节方面的已知功能相吻合，但二者之间的直接因果关系仍需进一步的功能实验验证。
与eQTL的结合
在已鉴定的SNP中，98.6%位于非编码区域（表S7）。为探究这些通过GWAS筛选出的SNP是否可能作为调控性变异发挥作用，我们将其与已报道的牛顺式表达数量性状位点（cis-eQTL）数据库进行了整合分析，以期阐明非编码SNP的功能调控逻辑。分析结果显示，共有208个（15.36%，208/3165）潜在eQTL位点与本研究发现的关联SNP存在重叠（p < 0.05；表S11）。值得注意的是，3’非翻译区（UTR）在此类重叠关系中呈现显著富集（42%，6/14；p = 0.0001，Fisher精确检验），提示顺式调控元件（特别是3'UTR区域）的突变可能具有重要的调控功能。在这些重叠位点中，我们识别出多个重要候选基因，如SLAMF1、PTK2和ARFGEF3。其中，帕尔默干旱严重指数（PDSI）显著关联的遗传区域，与FARP2基因在垂体组织中的cis-eQTL高度共定位（rcp = 0.96），SNP rs382549743和rs476116768关联性最强（图S20A）。ATAC-seq分析显示，在普通牛（海福特牛）中rs476116768附近存在明显的染色质开放峰（图2J）。此外，ChIP-seq分析表明该SNP位于绝缘子区域，并在TFAP2A转录因子结合基序中显著富集（p = 8.7E-04；图2J，图S21）。转录因子AP-2α已知与肾脏发育过程中的肾单位分化相关。rs476116768的衍生等位基因主要存在于欧洲普通牛中（图S20B），其频率与干旱波动幅度呈显著正相关（PDSI-Distance；R² = 0.265, p = 6.3E-07）。这些结果提示，rs476116768可能通过调控FARP2的表达参与和干旱适应相关的生理过程。该适应性表型背后的具体分子与生理机制，仍是未来研究的关键方向。
结 论
家牛已成功适应了从寒冷高原到炎热干旱地带等多样且复杂的气候环境。本研究整合全基因组关联分析、气候变量、表观基因组数据及顺式表达数量性状位点数据，构建高效分析框架，系统鉴定出与28个环境因子显著相关的3,165个单核苷酸多态性位点。并阐明了编码区与非编码区变异在环境适应中的差异化功能贡献，并揭示了趋同进化的遗传特征——例如PTK2基因在不同高海拔物种中的普通进化。本研究发现的候选变异位点（例如与耐寒性相关的CDK10基因、与干旱适应性相关的FARP2基因等）为培育气候适应性牛种提供了重要靶标。综上，本研究建立的整合分析框架为系统解析牛环境适应性遗传基础提供了高效的研究范式，也为应对气候变化挑战下的畜牧育种实践提供了科学依据与遗传资源。
方 法
样品制备与重测序
本研究采集了18个地方牛品种/群体的76份血样与样本并进行DNA提取。采用标准流程为每个样本构建了双端测序文库，平均插入片段为500 bp，读长为150 bp。将新测序数据与668个已公开的83个品种/群体全基因组数据进行整合，最终获得平均测序深度为11.28×的全基因组数据集。
变异检测与质控评估
基因分型流程参照“千牛基因组计划”第八阶段指南（http://www.1000bullgenomes.com/）执行。首先使用Trimmomatic去除低质量碱基及接头序列，随后通过BWA-MEM将高质量测序数据比对至牛参考基因组（ARS-UCD1.2）。利用SAMtools对BAM文件进行排序，再通过Picard的MarkDuplicates工具识别并标记PCR重复序列。使用基因组分析工具包GATK的BaseRecalibrator和PrintReads模块进行碱基质量值校正。采用GATK的HaplotypeCaller模块生成GVCF文件，继而通过GenotypeGVCFs和SelectVariants模块进行联合变异检测。进一步过滤非双等位SNP及基因型缺失率 > 0.1的位点。利用BEAGLE进行基因型填充与单倍型定相，并过滤DR2值 < 0.9的位点。最终获得的高质量SNP通过SnpEff进行功能注释，用于后续分析。
群体遗传学分析
将含34,902,931个SNP的VCF文件通过VCFtools转换为PLINK格式。采用滑动窗口法（窗口大小50 kb，步长20 kb）估计不同地理群体牛的全基因组核苷酸多样性。通过PLINK计算品种/群体间的FST遗传距离及近交系数。基于成对遗传距离矩阵采用邻接法构建系统发育树，并通过MEGA7可视化。使用PLINK的--indep-pairwise 50 10 0.1参数对基因型数据进行连锁不平衡修剪。随后基于同一SNP数据集（7,057,330个SNP）通过EIGENSOFT中的smartpca程序进行主成分分析，并利用ADMIXTURE分析群体结构。为确定最佳聚类数（K），我们对K = 2至10进行交叉验证。计算了五个筛选品种的核苷酸多样性。最终基于731个牛样本的32,506,509个SNP，通过FST估计品种/群体间的遗传距离，并采用SplitsTree4中的邻接网络算法基于距离矩阵构建无根网络。
气候因素收集
根据各品种采样点坐标获取对应的环境变量数据。共收集28个环境变量，涵盖温度（6个）、湿度（13个）、风速（3个）、气压（1个）、辐射与云量（2个）及植被（1个）等类别（附表S5）。其中7个环境变量数据集来源于CRU，21个变量来自NOAA PSL，空间分辨率为0.5°至2.5°网格。本研究使用netCDF4从.nc文件中提取近40年的逐月气候数据，剔除无效值（如缺失数据）后，基于逐月数据计算各年份平均气候值，并统计每年的月最高值、最低值和极差，以全面评估气候变化趋势与幅度。通过核密度图展示家牛生存的气候分布特征（图S5-S9）。为探究环境变量间的相关性，我们使用RStudio对数据集中环境变量的年均值进行Spearman等级相关性检验，评估变量间的长期线性关系，并借助ggcorrplot、corrplot和viridis等包可视化相关矩阵。最后通过gord包对所有数据（包括年均值、最大值、最小值及其差值）进行主成分分析。
全基因组环境关联分析
本研究使用PLINK对731个牛样本进行最终筛选。基于计算的杂合度分布，过滤高杂合度个体（高偏差可能提示样本污染或近交）。所有样本的平均杂合度均未偏离均值±3个标准差。基于IBS/IBD分析识别可能存在亲子关系的基因组并予以剔除，最终获得32,300,781个高质量SNP用于后续分析。
以各品种气候变量的年均值、最小值、最大值及极差作为分析参数（表S6）。为在控制群体遗传结构干扰的同时识别气候适应性相关标记，我们联合应用线性混合模型（GEMMA）、潜因子混合模型（LFMM）和贝叶斯方法（Baypass）。混合线性模型设定如下：y = Wα + xβ + u + ε；其中u ∼ MVNn(0, λτ⁻¹K)，ε ∼ MVNn(0, τ⁻¹In)。模型中y为表型向量，W为协变量矩阵，α为固定效应，x为基因型向量，β为标记效应大小，u为随机效应，ε为误差项，τ⁻¹为残差方差，λ为方差组分比，K为亲缘关系矩阵，In为单位矩阵。为控制群体分层与潜在混杂因素，将前20个主成分作为协变量纳入模型，该数量可解释数据集中99.91%的遗传变异。同时使用LFMM算法计算所有SNP的|z|分数。通过基因型数据的主成分分析确定最佳潜因子数；采用R语言prcomp函数进行主成分分析，根据解释方差拐点确定K = 6。为严格控制多重检验，对所有方法均采用Bonferroni校正，全基因组显著性阈值设为0.05除以总SNP数量。在BayPass分析中，将每条常染色体划分为包含不超过100,000个SNP的子数据集，共获得377个子集。在核心模型下对每个染色体数据集进行三次独立运行，随后基于辅助协变量模型评估SNP与气候适应的关联性。为验证结果一致性，采用Forstner-Moonen距离评估三次独立运行的协方差矩阵。以三次独立运行的中位数贝叶斯因子作为评估指标，设定30 dB作为显著性阈值以控制假阳性。为排除强群体分化造成的假关联，基于BayPass核心模型的XtX统计量设定严格过滤标准：剔除FDR=0.5%阈值以上的SNP位点。最终整合三种方法的结果，筛选出3165个高置信度SNP，要求至少在两种方法中显著。
家牛气候适应性因果变异的鉴定
本研究从NCBI获取了32个已发表的ATAC-seq数据集，涵盖两个成年海福特牛、三个青春期后婆罗门母牛、三个成年涠洲牛（未发表数据）及两个成年鲁西牛的八种组织（脂肪、肝脏、脾脏、肺、肌肉、大脑皮层、下丘脑与甲状腺）（表S9-10）。使用Trim Galore进行接头检查与修剪后，通过Bowtie2将ATAC-seq高质量读段比对至ARS-UCD1.2参考基因组。利用SAMtools对BAM文件排序，Picard去除重复比对序列。通过deepTools对ATAC-seq BAM文件进行RPKM标准化。使用MACS2分别对各重复样本进行峰检测，并通过bedtools合并同一组织内三个生物学重复的峰区域。从相同样本（两个成年海福特牛）中获取61个已发表的ChIP-seq数据集（表S9-10），采用相同流程处理ChIP-seq数据。通过Integrative Genomics Viewer观察染色质开放区域峰图。
为比较冷适应与热适应牛群的基因表达差异，我们整合了两个公开的海福特牛肌肉组织RNA-Seq数据集和三个新测序的涠洲牛肌肉组织数据集（表S9-10）。涠洲牛样本总RNA使用Illumina HiSeq 2500系统测序，原始读段经FASTP处理后获得高质量数据。基于洁净数据，采用HISAT2和StringTie将序列比对至参考基因组并进行转录本组装。利用基因组注释获得的基因长度将转录本计数标准化为每百万转录本数。通过limma包进行差异表达分析，采用经验贝叶斯 moderated t检验评估两组牛群间的基因表达差异，并使用Benjamini-Hochberg方法校正p值以控制错误发现率。
为鉴定基因表达与环境适应性复杂性状的共同因果变异，我们从https://cgtex.roslin.ed.ac.uk/获取顺式eQTL数据。通过bedtools intersect工具，系统比对了乳腺、脾脏、肾脏、肌肉（瘤牛/普通牛/杂交牛）、下丘脑、垂体、子宫、巨噬细胞、脂肪、空肠、输卵管、白细胞、瘤胃、卵巢、肌内脂肪、睾丸、胚胎、血液、肝脏、肺、乳细胞、单核细胞、淋巴结及乳腺L等27个组织的显著顺式eQTL位点（表S11）。进一步使用Coloc包的coloc.abf函数进行eQTL与气候关联SNP的共定位分析，获得两种性状共享单一因果变异的后验概率，保留PP.H4 > 0.8的位点作为高置信度共享因果变异。
环境适应性相关候选基因的筛选
本研究利用bedtools intersect与SnpEff工具在基因区域内鉴定出重要单核苷酸多态性位点。为验证所筛选候选基因是否与特定功能类别相关，采用DAVID平台的京都基因与基因组百科全书（KEGG）及基因本体（GO）数据库进行功能富集分析，以FDR校正后p值小于0.1作为显著富集标准。
为识别群体特异性适应变异，我们根据不同目标群体的遗传结构特征，制定了基于衍生等位基因频率的严格筛选标准：针对雅库特牛与恩达马牛，筛选标准为目标群体DAF > 0.7且其他所有牛群DAF < 0.3；鉴于藏牛群体同时包含普通牛与瘤牛血源的遗传异质性，采用改良筛选方案，将非藏牛群体DAF < 0.3且至少两个藏牛品种DAF > 0.5的位点定义为藏牛特异性变异。为评估18号染色体候选区域的选择特征，使用PopLDdecay分析连锁不平衡衰减特征：在包含选择清除区域的1 Mb窗口内计算成对r²值，并以基因组随机选取的1 Mb区域作为对照，分别绘制冷适应与热适应群体的连锁不平衡衰减曲线。此外，采用滑动窗口计算不同群体候选区域的核苷酸多样性及Tajima’s D统计量，整合多维度选择信号验证目标基因组区域在特定群体中受到选择的作用强度。通过R语言rehh软件包进行扩展单倍型纯合性分析，计算候选SNP位点的单倍型纯合性衰减特征。最后利用LDBlockShow绘制候选基因组区域的局部连锁不平衡区块图谱以可视化其遗传结构。
代码和数据可用性：
本研究涉及的全部测序数据已上传至NCBI序列读段归档数据库，登录号如下：PRJNA905718、PRJNA1085861、PRJNA658727。相关分析所用的数据及脚本均已保存于GitHub平台（https://github.com/Luyang-Sun/Climate-adaptation-in-cattle.git）。补充材料（实验方法、图表、中文翻译版本及更新资料）可通过在线DOI或在iMeta Science期刊官网（http://www.imeta.science/）获取。支撑本研究结果的其他相关数据详见本文补充材料部分。
引文格式：
Luyang Sun, Liangliang Jin, Ting Sun, Chenqi Bian, Fuwen Wang, Shuang Li, Hao Li, et al. 2026. “Climate factors drive the local adaptation of Old World cattle.” iMeta 4: e70108. https://doi.org/10.1002/imt2.70108.
孙露洋（第一作者）
● 西北农林科技大学动物科技学院在读博士研究生。
● 研究方向为牛的基因组学和进化，以第一作者在iMeta、Frontiers in Genetics、Animal Biotechnology期刊发表论文。
金良粱（第一作者）
● 西北农林科技大学动物科技学院在读博士研究生。
● 研究方向为单细胞多组学与调控基因组学研究，以第一作者（含共同）在iMeta、Frontiers in Genetics、Biology期刊发表论文。
孙婷（第一作者）
● 博士，2021年毕业于西北农林科技大学动物科技学院。
● 研究方向为动物基因组学，主持国家自然科学基金青年项目、中国博士后科学基金各一项。以第一作者/通讯作者在GigaScience、Journal of Genetics and Genomics、BMC Genomics、Animal Genetics、Journal of Animal Breeding and Genetics等期刊发表论文。
马云（通讯作者）
● 宁夏大学动物科技学院教授，博士生导师。
● 研究方向为反刍动物遗传育种，运用系统生物学理论与多组学整合策略，从进化角度系统解析牛等反刍动物重要经济性状与适应性性状的遗传调控基础。先后主持国家自然科学基金区域创新联合基金重点项目、国家自然科学基金面上项目、宁夏回族自治区重点研发重大项目及科技创新领军人才项目等多项国家与省部级重要课题。其研究成果发表于 Nature Communications、International Journal of Molecular Sciences、Frontiers in Veterinary Science、BMC Genomics等专业期刊，并参与《养牛学》《中国黄牛学》等教材与专著的编写工作。
夏小婷（通讯作者）
● 西北农林科技大学动物科技学院，青年教授，博士生导师。
● 研究方向为地方黄牛遗传资源与功能基因挖掘。主持国家自然科学基金青年项目，以第一作者或并列第一作者在Genome Biology、Advacned Science、Nature Communications、Science Bulletin、Heredity、BMC Genomics、Animal Genetics、Stress Biology等期刊发表论文10余篇。
雷初朝（通讯作者）
● 西北农林科技大学动物科技学院教授，博士生导师。
● 研究方向为动物遗传资源与分子育种。目前为国家肉牛牦牛产业技术体系岗位专家，第四届国家畜禽遗传资源委员会牛专业委员会委员。兼任中国畜牧兽医学会畜禽遗传标记学分会副理事长，中国畜牧兽医学会动物遗传育种学分会理事，中国畜牧兽医学会养牛学分会理事。主持国家“863”计划、国家自然科学基金、国家肉牛牦牛产业技术体系专项等课题30余项，其研究成果发表于Nature Communications、iMeta、Science Bulletin、Genome Biology等权威学术期刊。
陈宁博（通讯作者）
● 西北农林科技大学动物科技学院副教授，博士生导师。
● 研究方向为中国地方黄牛遗传资源与优异基因挖掘研究，重点解析家牛群体演化规律、地方牛种质特性与抗逆性状的遗传机制，开发基于基因组信息的分子育种标记与技术体系。主持国家自然科学基金青年基金（B类、C类）、生物育种-青年专项、国家重点研发计划课题等8项。以第一作者和通讯作者（含共同）在Nature Communications (2篇)、iMeta、Science Advances、Genome Biology、Science Bulletin (3篇)、PNAS、Advanced Science和GPB等期刊发表论文20余篇，授权国家发明专利9件。兼任中国畜牧兽医学会畜禽遗传标记分会常务理事、中国畜牧兽医学会养牛学分会理事、陕西省遗传学会常务理事，BMC Biology、GPB、AROH、JGG和JIA等杂志编委及青年编委。
（▼ 点击跳转）
高引文章 ▸▸▸▸
iMeta | 引用20000+，海普洛斯陈实富发布新版fastp，更快更好地处理FASTQ数据
高引文章 ▸▸▸▸
iMeta | 兰大张东组：使用PhyloSuite进行分子系统发育及系统发育树的统计分析
高引文章▸▸▸▸
iMeta | 唐海宝/张兴坦-用于比较基因组学分析的多功能分析套件JCVI
1卷1期
1卷2期
1卷3期
1卷4期
2卷1期
2卷2期
2卷3期
2卷4期
3卷1期
3卷2期
3卷3期
3卷4期
3卷5期
3卷6期
4卷1期
4卷2期
4卷3期
4卷4期
4卷5期
4卷6期
1卷1期
1卷2期
2卷1期
2卷2期
2卷3期
2卷4期
1卷1期
1卷2期
“iMeta” 是由威立、宏科学和本领域数千名华人科学家合作出版的开放获取期刊，主编由中科院微生物所刘双江研究员和荷兰格罗宁根大学傅静远教授担任。目的是发表所有领域高影响力的研究、方法和综述，重点关注微生物组、生物信息、大数据和多组学等前沿交叉学科。目标是发表前10%(IF > 20)的高影响力论文。期刊特色包括中英双语图文、双语视频、可重复分析、图片打磨、60万用户的社交媒体宣传等。2022年2月正式创刊！相继被Google Scholar、PubMed、SCIE、ESI、DOAJ、Scopus等数据库收录！2025年6月影响因子33.2，中科院分区生物学1区Top，位列全球SCI期刊前千分之三(65/22249)，微生物学科2/163，仅低于Nature Reviews，学科研究类期刊全球第一，中国大陆5/585！
“iMetaOmics” 是“iMeta” 子刊，主编由中国科学院北京生命科学研究院赵方庆研究员和香港中文大学于君教授担任，目标是成为影响因子大于10的高水平综合期刊，欢迎投稿！
"iMetaMed" 是“iMeta” 子刊，专注于医学、健康和生物技术领域，目标是成为影响因子大于15的医学综合类期刊，欢迎投稿！
iMeta主页：
http://www.imeta.science
姊妹刊iMetaOmics主页：
http://www.imeta.science/imetaomics/
出版社iMeta主页：
https://onlinelibrary.wiley.com/journal/2770596x
出版社iMetaOmics主页：
https://onlinelibrary.wiley.com/journal/29969514
出版社iMetaMed主页：
https://onlinelibrary.wiley.com/journal/3066988x
iMeta投稿：
https://wiley.atyponrex.com/journal/IMT2
iMetaOmics投稿：
https://wiley.atyponrex.com/journal/IMO2
iMetaMed投稿：
https://wiley.atyponrex.com/submission/dashboard?siteName=IMM3
邮箱：
office@imeta.science