利用正弦波建模方法评估气候对奶牛产奶量相关性状的影响

《animal》：Assessing climatic effects on milk yield traits using sinusoidal modeling in dairy cows

　　G. Bittante | G.J. Martínez-Marín | H. Toledo-Alvarado | N. Amalfitano | A. Cecchinato | S. Schiavon | L. Gallo
帕多瓦大学（意大利莱尼亚罗，35020）农学、食品、自然资源、动物与环境系（DAFNAE）
摘要
气候变化对奶牛的影响是当前科学研究的热门话题。本研究的目的包括：（i）测试正弦函数（SIN-F）描述气候特征及奶牛生产特征年度周期性变化的能力；（ii）比较平原地区与丘陵/山区奶牛场以及不同品种奶牛的应用效果；（iii）评估夏季期间偏离SIN-F的情况，以判断其是否可能指示热应激的影响。我们使用了来自39个气象站的数据以及117个牛群在1440个测试日内的记录，来分析年度周期性SIN-F模型及月度偏差对日产量（dMY，千克）、脂肪含量（FAT%）和蛋白质含量（PROT%）的影响。牛群本身（占变异总量的29.6%至59.8%）和泌乳阶段（18.3%至53.4%）是影响奶牛生产特征的主要因素，而胎次对日产量（dMY）的影响较为显著（27.1%）。年份因素仅占变异总量的0.9%至2.6%，长期SIN-F模型占2.4%至12.8%，而月度偏差占0.4%至2%。海拔高度的影响较小（0.1%至3.3%），品种因素占1.9%至4.3%，并且与SIN-F模型存在显著交互作用，但与海拔高度无关。奶牛的脂肪含量（FAT%）和蛋白质含量（PROT%）在12月至1月达到峰值，而温度则在6月至7月达到峰值，且这一时间点比光照周期晚一个月，其变化模式与光照周期相反。日产量（dMY）的SIN-F模型峰值出现时间比光照周期晚三个月（9月至10月），比温度峰值晚两个月。平原地区的日产量峰值与最低值差异（2.84千克/天）大于丘陵/山区（2.06千克/天），这可能反映了平原地区奶牛面临更大的气候压力。与荷斯坦奶牛相比，西门塔尔奶牛和布朗瑞士奶牛的日产量SIN-F波动幅度较小，这可能部分归因于它们的日产量较低，同时也可能因为它们对气候变化的抵抗力更强；此外，日产量SIN-F模型的峰值出现时间比荷斯坦奶牛晚一个月，这也可能与它们对气候变化的适应能力更强有关。日产量在6月至7月期间显著下降，这可能是热应激的表现。年度周期性SIN-F模型及其变化被证明是评估气候条件对不同牛群、品种和养殖系统奶牛生产特征影响的有效工具。

**意义**
热应激已成为全球奶牛养殖系统面临的主要挑战。研究结果表明，利用正弦函数对年度周期性变化进行建模可以为评估不同品种、海拔高度和养殖系统下的气候影响提供实用工具。奶牛的脂肪含量和蛋白质含量与光照周期和温度呈负相关，而日产量则表现出延迟的变化模式。平原地区的奶牛产量波动幅度较大，表明它们对气候变化的敏感度更高，而荷斯坦奶牛比西门塔尔奶牛和布朗瑞士奶牛更易受气候波动影响。因此，正弦函数建模成为评估气候变化对奶牛生产特征影响的有效工具。

**引言**
季节性气候变化和热应激对多个国家的奶牛养殖系统构成了日益严峻的挑战，尤其是在全球变暖的背景下（Thornton等人，2022年）。许多研究将动物（尤其是泌乳奶牛）的生产、繁殖、生理和病理变化与环境气候条件（光照周期、空气温度和湿度）的变化联系起来。光照周期完全由天文因素决定，具体来说是地球自转轴相对于太阳黄道平面的倾斜角度。像许多振荡模式一样（Forger，2017年），光照周期的变化遵循年度正弦函数（SIN-F），这与地球的轴倾角有关。时间（t）的正弦函数（Sundstrom和Schlicker，2023年）表示为：y = A sin(2π/B(t-C)) + D，该函数由四个参数定义：A代表函数的“振幅”——即特征y的中心线D与最大值（峰值）和最小值（最低值）之间的距离；B代表函数的“周期”——即一个完整周期的长度；C代表函数的“相位偏移”——即y等于D的第一个X值；D代表正弦函数的“中心线”y? = D。在年度周期性变化的情况下，B通常假设为1年（即B=12个月），D代表特征y的平均值。因此，需要计算的唯一两个参数是振幅和相位偏移。光照周期（或白昼长度）在夏至达到最大值（峰值），在冬至达到最小值（最低值）。光照周期的长度反映了到达地球表面的太阳辐射量，从而可以作为每日地表太阳辐照度的指标。

**方法**
本研究分析了之前项目（MeetBull项目）收集的数据，该项目涉及意大利威尼托地区（意大利东北部）的131个商业奶牛场，数据覆盖了三年时间。所有农场都参与了威尼托养牛协会（ARAV – Associazione Regionale Allevatori del Veneto，维琴察，意大利；www.arav.it）的牛奶记录计划。关于养殖系统、奶牛特征和繁殖特征的信息已在先前文章中报道（Bittante等人，2020a），而关于犊牛特征及年轻公牛的表现和特征的信息则见于其他文章（Bittante等人，2020b、2021b、2023）。本研究未使用任何人类或动物，因此无需获得动物护理和使用方面的机构批准。

**数据分析**
剔除测试日观测次数少于100次的牛群后，最终保留了117个牛群用于分析。以海拔350米为分界线，将平原地区（105个农场）与丘陵/山区（12个农场）的养殖系统区分开来。表1列出了参与研究的农场海拔高度的描述性统计信息。大多数牛群饲养在现代化的集约化奶牛场中，这些农场提供的饲料包括玉米青贮饲料、部分干草和浓缩饲料。每个奶牛场的平均奶牛数量为98±9头。共有18,658头奶牛参与了这项调查。纯种奶牛分别注册在荷斯坦（14,298头，ANAFIBJ）、布朗瑞士（2,707头，ANARB）和西门塔尔（556头，ANAPRI）品种的相关牛群登记册中，而杂交奶牛（1,097头，AIA）则记录在杂交品种登记册中。每月的牛奶产量数据来自ARAV。测试日的记录包括牛奶产量（dMY，千克/天）、牛奶脂肪含量（FAT%）和牛奶蛋白质含量（PROT%），这些数据是通过Milkoscan FT7DC（FOSS Electric）傅里叶变换红外光谱仪在ARAV牛奶实验室（意大利维琴察）测得的。所有数据收集和处理均符合ICAR（国际动物记录委员会）的指南（Orlandini，2020年）。

表1. 研究区域的地理和气候特征。

| 项目 | 所有农场 | 平原 | 山丘和山脉 |
| ---- | -------- | -------- | -------- |
| 奶牛测试日数量 | 65,883 | 360,642 | 25,241 |
| 农场海拔：平均值±标准差（米） | 148±26 | 565±55 | 872±27 |
| 海拔范围（米） | 1 - 1,183 | 1 - 2,313 | 69 - 1,183 |
| 地理坐标：纬度范围 | 44.90 – 46.02 N | 44.90 – 45.97 N | 45.59 – 46.02 N |
| 经度范围 | 10.71 – 12.52 E | 10.71 – 12.52 E | 10.95 – 12.32 E |
| 气候（K?ppen-Geiger分类） | Cfa | Cfb | Cfb |
| 描述 | 潮湿亚热带 | 海洋性 | 温和 |
| 夏季天气 | 炎热潮湿 | 寒冷至温和 |
| 气象站数量 | 39 | 45 | 39 |
| 气象站海拔（米） | 147±26 | 69±59 | 82±33 |
| 农场与气象站距离（公里） | 6.3±3.4 | 6.3±3.3 | 5.5±4.2 |
| 农场与气象站海拔差（米） | -1±50 | 4±17 | -44±14 |
| 光照周期（小时） | 12月 | 8.56±0.03 | 8.57±0.03 | 8.53±0.02 |
| 6月 | 15.44±0.03 | 15.43±0.03 | 15.47±0.02 |
| 平均气温（°C） | 1月 | 3.1±2.1 | 3.4±1.8 | 1.3±2.8 |
| 7月 | 24.6±1.9 | 25.2±1.1 | 20.9±2.6 |
| 最低气温（°C） | 1月 | -0.7±2.6 | -0.5±2.3 | -2.1±3.9 |
| 最高气温（°C） | 1月 | 7.9±1.7 | 8.2±1.5 | 5.7±2.1 |
| 备注：asl = 海拔 |

气象站和气候信息
研究区域（意大利威尼托地区）的地理坐标和气候特征总结在表1中。气候分类采用K?ppen-Geiger系统（Geiger，1954年），区分了威尼托地区的平原和山丘地带。气候数据来自威尼托地区环境预防与保护局（ARPAV，2023年）。我们保留了距离农场最近的39个气象站的数据。气候数据包括最低气温、最高气温和平均气温（分别以°C为单位）。光照周期（PhotoPeriod，以小时/天为单位）是根据Allen等人（1998年）的建议公式，根据每个气象站的纬度计算得出的：

N = 24πarccos(tanφtanδ)
其中：φ是从十进制度数转换为弧度的纬度（φ=π/180[十进制度数]）；δ是太阳赤纬角（δ=0.409 sin[(2π/365) J?1.39]），J是一年中的天数（1至365）。因此，光照周期在每年是相同的，所以只保留了一年的数据。气象站根据与研究中的奶牛农场相同的标准被分类为平原或山丘地带。

正弦函数（SIN-F）
SIN-F解释了北极圈和南极圈之间给定纬度的年周期变化，不包括赤道（光照周期恒定）。光照周期的相位移动将第一个光照周期为12小时的日子（白天和黑夜时间相等）放在3月的春分点，第二个日子放在6月的秋分点（来自拉丁语aequinoctium：aequus + nox = 相等 + 夜晚，意味着白天和黑夜时间相等）。6月的至日（来自拉丁语solstitium：sol + sisto + -ium = 太阳 + 停留，意味着光照周期长度稳定，达到最大或最小）位于3月和9月春分点之间，是北半球光照周期最长的日子（南半球则相反），而12月的至日是北半球光照周期最短的日子（位于9月和3月春分点之间）。两个至日之间的光照周期差异是SIN-F振幅的两倍，并且随着纬度的增加而增加。将SIN-F纳入统计模型以研究特定性状的年周期时间模式具有以下优点：a) 当研究的性状受光照周期影响时，其原理具有科学依据；b) 它提供了描述性状变化的参数：效应的振幅及其季节性延迟（相位移动）；c) 它在统计上较为简洁，即优化了解释的方差与模型额外自由度之间的比例（Forger，2017年）。此外，我们假设动物性状的振幅可以解释为抵抗力的指标，其相位移动可以解释为动物对年度气候变化的适应能力。

数据分析
所有缺失值或异常值（超出平均值±3标准差范围的记录）都被丢弃。

建模年周期气候和牛奶生产性状
温度数据（Tmean、Tmin和Tmax，以°C为单位）涉及每个气象站四年的月平均值（12个月×39个气象站×4年=1,872个观测值）。由于光照周期在每年是恒定的，因此只参考了一年的数据（12个月×39个气象站=468个观测值）。年内的变化被认为是长期年周期模式和其中期月度偏差的结果。由于独立观测值相对较少（测试日期通常每月一次，并且所有奶牛共有的，因此只有不同的牛群才是真正独立的），因此没有对短期每周或每日变化进行建模。长期成分假设为周期为1年的SIN-F（B = 12个月），而中期成分由包含SIN-F的模型得到的日历月度偏差表示。为了估计SIN-F的参数（A、C和D），我们使用了SAS（版本9.4，SAS Institute Inc., Cary, NC）中的PROC NLMIXED程序。对于光照周期分析，模型调整了气象站的影响；而对于温度分析（Tmean、Tmin和Tmax），模型同时调整了气象站（WS）和年份的影响。光照周期分析仅使用了一年的数据，因为这一天文性状在不同年份在同一地点没有变化。

牛奶生产性状（dMY、FAT%和PROT%）的SIN-F参数估计使用了与气候性状相同的PROC NLMIXED程序，将WS效应替换为Herd（Hn，n=1至117），并在模型中包括了由于个体奶牛引起的变异来源：产次（PAp，p=第一次、第二次、≥第三次）、泌乳阶段（LSq，q=每类30个DIM）、品种（BRr，r=荷斯坦、布朗瑞士、西门塔尔和杂交品种）及其与SIN-F的交互作用（BRr×SIN_Fj）和海拔（BRr×ALl）。

量化气候和牛奶生产性状的变异来源
根据第二个统计模型，估计了不同变异来源对每个性状的相对贡献，表示为每个来源占总平方和的百分比。由于一些变异来源，特别是SIN-F和月份，高度相关，我们使用了I型方差分析（ANOVA）来逐步划分总变异，顺序为：年份、SIN-F和月份，从而从较宽的时间尺度逐步细化到较窄的时间尺度。使用I型或III型ANOVA以及是否包含SIN-F或Mo及其在模型中的位置的含义已在补充材料S1中详细描述。分析使用SAS（版本9.4）中的PROC GLM程序进行，根据以下模型进行温度性状分析：

yijklmn = μ + YEi + SIN_Fj + MOk + ALl + ALl×SIN_Fj + WSm + YEi×WSm + MOk×WSm + e
其中：y是温度性状（Tmean、Tmin或Tmax）；μ是总体平均值；YEi是日历年份的固定效应（i=2015至2018）；SIN_Fj是对SIN-F的回归（两个参数：A是振幅，C是相位移动，B假设为12个月；D等于0）；MOk是日历月份的固定效应（k=1月至12月）；ALl是气象站的海拔类别（l = 平原、山丘和山脉）；ALl×SFj是海拔和SIN-F的交互作用；WSm是气象站的固定效应（m=1至39）；YEi×WSm和MOk×WSm分别是气象站与年份和月份的交互作用；e是残差项。

光照周期数据使用与温度性状相同的模型进行分析，但去除了YEi和YEi×WSm元素。用于分析牛奶生产性状（dMY、FAT%和PROT%）的模型基于用于气候性状的模型，将WS效应替换为Herd（Hn，n=1至117），并在模型中包括了由于个体奶牛引起的变异来源：产次（PAp，p=第一次、第二次、≥第三次）、泌乳阶段（LSq，q=每类30个DIM）、品种（BRr，r=荷斯坦、布朗瑞士、西门塔尔和杂交品种）及其与SIN-F的交互作用（BRr×SIN_Fj）和海拔（BRr×ALl）。由于PROC GLM模型在MO之前包含了SIN-F，因此MO效应的解表示每个月相对于SIN-F的偏差。

结果和讨论
通过分析使用I型或III型ANOVA以及是否包含SIN-F或Mo及其在统计模型中的位置对评估气候和牛奶产量性状不同变异来源的重要性的影响（详见补充材料S1），发现同时包含这两个高度相关的季节性因素（SIN-F和Mo）在III型ANOVA下会导致它们的重要性被严重低估。使用I型ANOVA时，所有共享的方差都被分配了，但结果取决于包含的顺序。我们假设当模型中包含的因素之间存在逻辑或因果顺序时，使用I型ANOVA是合适的。我们将SIN-F放在Mo之前的分析顺序是基于SIN-F代表了由天文因素决定的预期年周期气候模式，而Mo代表了SIN-F捕捉到的长期气候趋势的可能偏差。我们在后续章节中对结果的解释基于这一假设，这一点应该予以考虑。

正弦函数用于描述气候性状以及海拔和气象站的影响
本研究的第一个目的是评估SIN-F在解释光照周期和气温性状的年周期变化方面的可靠性，以便将这些函数作为比较牛奶生产性状结果的基准。

光照周期（白天）
光照周期似乎不受海拔或个别气象站的影响（表2），因为纬度差异很小，考虑到研究区域与45度纬度的偏差不大。

表2. 光照周期和气温：描述性统计、模型中包含的变异来源的百分比发生率和显著性，以及描述气候性状年周期变化的正弦函数（SIN-F）参数的估计。

| 项目 | 光照周期（小时/天） | 气温（°C） |
| ---- | -------------- | -------------- |
| 观测值数量 | 468 | 1,872 | 1,872 | 1,872 |
| 描述性统计 | 平均值 | 11.99 | 8.91 | 13.61 | 18.9 |
| 标准差 | 2.38 | 6.97 | 7.78 | 8.5 |
| CV% | 19.88 | 77.55 | 6.64 | 5.0 |
| P值 | 18.52 | -4.5 | 0.15 | 5.2 |
| P99 | 15.47 | 20.72 | 6.83 | 3.1 |
| 变异来源百分比： | | | |
| 年份 | -0.41 | *** | 0.13 | *** |
| SIN-F | 99.91 | *** | 93.03 | *** |
| 95.26 | *** | 94.87 | *** |
| 日历月份 | 0.06 | *** | 0.88 | *** |
| 0.83 | *** | 0.94 | *** |
| 海拔 | 0.00 | n.s. | 1.71 | *** | 2.16 | *** |
| 海拔 × SIN-F | 0.03 | n.s. | 0.10 | *** | 0.14 | *** |
| 气象站 | 0.00 | n.s. | 3.47 | *** | 1.17 | *** |
| 气象站 × 年份 | 0.00 | n.s. | 0.02 | n.s. | 0.03 | n.s. |
| 气象站 × 月份 | 0.00 | n.s. | 0.35 | n.s. | 0.29 | n.s. |
| 模型R2 | 99.99 | 97.71 | 97.92 | 97.61 |

SIN-F参数
| 振幅（A） | a | 3.36 | 9.34 | 10.58 | 11.54 |
| 相位移动（C） | b | 2.66 | 3.59 | 3.52 | 3.53 |

注：P1 = 第1百分位数；P99 = 第99百分位数。*** P < 0.001；n.s.：不显著（P > 0.05）。
a：正弦曲线最高点或最低点与其中间值之间的差异，以性状测量单位表示（光照周期为小时/天，气温为°C）。
b：相位移动以1月1日至春分点之间的月份数表示（中间值是最高点或最低点之间的日子）。

因此，正如预期的那样，只有年度内的变化是显著的。统计模型中包含的SIN-F捕捉了99.91%的年周期变化，而中期月度偏差仅解释了0.06%的变化，这是由于不同月份的天数不同。图1显示了光照周期的SIN-F模式，符合年度的预期：12月21日光照时间最短（从日出到日落8.62小时，从日落到日出15.38小时），6月21日光照时间最长（光明和黑暗的时间比例相反），分别对应冬季和夏季至日。3月20日和9月21日的春分点时，光明和黑暗的时间相等（各12.00小时）。这条曲线是北半球（北纬45度左右）光周期的经典表示，而在南半球（南纬45度左右），曲线则相反，夏至和冬至也相应反转。SIN-F的振幅（夏至和春分之间的光周期差异）为每天3.36小时，最大差异（两个夏至之间）是振幅的两倍，即每天6.72小时。相位差，即1月1日到春分（3月21日）之间的间隔，为2.66个月。已知光周期的振幅从热带（北半球的北回归线，南半球的南回归线）向极圈（北极圈和南极圈）逐渐增加。

图1. 光周期（黄色实线曲线）和气温特征（Tmean，红色实线曲线；Tmin，红色虚线曲线；Tmax，红色点线曲线）：来自意大利威尼托地区39个气象站收集的数据的年循环正弦函数。

SIN-F在表示气温特征（Tmean、Tmin和Tmax）的年循环变化方面表现非常好，分别解释了总变化的93.03%、95.26%和94.87%（表2）。SIN-F的顶点和最低点之间的差异约为20°C，并且随着特征平均值的变化而增加（Tmin为18.68°C，Tmean为21.18°C，Tmax为23.08°C）。正如预期的那样，这三个气温特征的SIN-F曲线基本平行（图1）。然而，由于相位差较长（大约一个月），顶点和最低点并不与夏至和冬至重合，导致峰值（顶点）和最低值分别出现在7月和1月。这种延迟与光周期相比（见表2中的相位差），反映了这两个变量的不同性质：光周期是一个流量特征（以每天太阳辐射小时数衡量），而气温是一个存量特征，代表大气中热量的累积。

奶牛场海拔类别对温度变化的解释占1.71%到2.48%，同一海拔类别内的各个气象站对温度变化的解释比例也类似（1.01%到3.47%，表2）。从表1的描述性统计可以看出，平原地区的气温比丘陵和山区更高，尤其是在夏季，尤其是Tmax。海拔类别与SIN-F之间的相互作用显著，而气象站与年份和月份之间的相互作用则不显著。

因此，关于本研究的首要目标，SIN-F很好地描述了光周期和气温（Tmin、Tmean和Tmax，单位为°C）的年循环变化。

本研究的第二个目标是评估SIN-F描述奶牛生产特征年循环变化的可靠性，并将这些正弦曲线与气候曲线进行比较。表3报告了奶牛生产特征的描述性统计和ANOVA I型分析结果。模型的决定系数对于dMY和PROT%接近45%，而对于FAT%则低得多（16%）。奶牛特征的高残差变化是由于分析的数据是来自单个测试日的记录。因此，模型解释了超过80%的奶牛特征变化是由牛群、品种、胎次、奶牛的泌乳阶段及其相互作用造成的并不令人意外。然而，我们主要关注的是由时间（年份、年循环SIN-F和月度偏差）以及农场海拔和奶牛品种解释的变化，因此没有展示和讨论牛群和单个奶牛因素的结果。

描述性统计显示，平原地区的气温比丘陵和山区更高，尤其是在夏季，尤其是在Tmax方面。海拔类别与SIN-F之间的相互作用显著，而气象站与年份和月份之间的相互作用则不显著。

因此，就本研究的首要目标而言，SIN-F证明了它能够很好地描述光周期和气温（Tmin、Tmean和Tmax，单位为°C）的年循环变化。

本研究的第二个目标是评估SIN-F描述奶牛生产特征年循环变化的可靠性，并将这些正弦曲线与气候曲线进行比较。表3报告了奶牛生产特征的描述性统计和ANOVA I型分析结果。模型的决定系数对于dMY和PROT%接近45%，而对于FAT%则低得多（16%）。奶牛特征的高残差变化是由于分析的数据是来自单个测试日的记录。因此，模型解释了超过80%的奶牛特征变化是由牛群、品种、胎次、奶牛的泌乳阶段及其相互作用造成的并不令人意外。然而，我们主要关注的是由时间（年份、年循环SIN-F和月度偏差）以及农场海拔和奶牛品种解释的变化，因此没有展示和讨论牛群和单个奶牛因素的结果。

图2. 奶牛生产特征的年循环正弦函数（SIN-F）（a：日产量，蓝色实线曲线；b：脂肪含量，橙色实线曲线；c：蛋白质含量，绿色实线曲线）与光周期的SIN-F（左侧图表，黄色实线曲线）和气温特征（右侧图表，红色曲线）以及SIN-F的月度偏差（黑色圆圈，SE，*P<0.05；**P<0.01；***P<0.001）的比较。a）：日产量（dMY）；b）：脂肪含量（FAT%）；c）：蛋白质含量（PROT%）。

描述年循环变化的正弦函数用于描述牛奶生产特征的年循环变化。SIN-F描述了牛奶生产特征的年循环变化，以评估其可靠性。

SIN-F描述牛奶生产特征的年循环变化的第二个目标是评估SIN-F描述奶牛生产特征年循环变化的可靠性，并将这些正弦曲线与气候曲线进行比较。表3报告了奶牛生产特征的描述性统计和ANOVA I型分析结果。模型的决定系数对于dMY和PROT%接近45%，而对于FAT%则低得多（16%）。奶牛特征的高残差变化是由于分析的数据是来自单个测试日的记录。因此，模型解释了超过80%的奶牛特征变化是由牛群、品种、胎次、奶牛的泌乳阶段及其相互作用造成的并不令人意外。然而，我们主要关注的是由时间（年份、年循环SIN-F和月度偏差）以及农场海拔和奶牛品种解释的变化，因此没有展示和讨论牛群和单个奶牛因素的结果。

表3. 奶牛特征：描述性统计、模型中包含的变化来源的百分比发生率及显著性，以及描述牛奶生产特征年循环变化的正弦函数（SIN-F）的参数估计。

奶牛生产特征的年循环变化与光周期和气温变化的关系基于所采用的假设和统计模型，SIN-F描述了奶牛生产特征的大部分季节性变化。应当指出，参与调查的大多数农场全年使用成分几乎恒定的全混合日粮（TMR），并且统计模型考虑了所有可能影响奶牛生产特征年循环变化的主要变化来源（牛群、品种、胎次、泌乳阶段和农场海拔）。在许多基于牧场的奶牛养殖系统中，这种情况并不成立，因为产犊通常是季节性的，以便使奶牛的营养需求与饲料的可用性和质量相匹配。在这些系统中，季节、泌乳阶段以及饲料的可用性和质量的影响相互混淆，难以区分它们对奶牛生产特征变化的单独贡献。在我们研究的集约化室内养殖系统中，奶牛特征年内变化的高比例由正弦曲线形状表示，这支持了牛奶与气候特征之间的强烈关系，尽管需要考虑不同的相位差。但主要未解决的问题是：是光周期还是气温特征对奶牛生产特征变化有更大影响？

光周期和气温在年循环基础上显然高度相关。在我们的案例中，光周期的月平均值与Tmean之间的皮尔逊相关系数为90.2%。许多文章报告了气温特征（Herbut等人，2018年）或其衍生的温度-湿度指数（THI）（国家研究委员会，1971年）与牛奶生产特征之间的相关性，并推断出温度效应，特别是热应激的影响。可以预期，如果用光周期代替温度，也会得到类似的相关性。因此，除非同时考虑这两种气候特征，否则无法确定这些效应是依赖于温度还是光周期。在接下来的讨论中，我们改变了奶牛特征的顺序，首先讨论蛋白质含量，因为结果更为清晰。

蛋白质含量在夏季达到最低值，此时光周期持续时间和气温都处于最高值。在温带地区的集约化室内奶牛养殖系统中通常会发现类似的模式（Bernabucci等人，2015年），但在基于牧场的系统中则不是这样，因为在这些系统中，较长的光周期与饲料（和粗蛋白）的最大可用性相关。Salfer等人（2019年）使用正弦（余弦）函数对美国11个市场的牛奶平均脂肪含量和蛋白质含量进行了建模。在蛋白质含量的情况下，该函数在所有市场区域都显著，振幅范围为0.07%至0.14%（我们的平均值为0.11%），峰值出现在12月29日至1月19日期间，最低值出现在6月底至7月初期间。在这种情况下，蛋白质含量的变化预计与光周期和气温特征的变化高度负相关。在我们的案例中，蛋白质含量的变化大约比光周期早两个月，比气温特征早三个月。这意味着蛋白质含量的最低值并不是在光周期最长和气温最高的月份出现的，而是在光周期缩短的月份出现的。同样的作者在分析宾夕法尼亚州样本牛群中单个奶牛的牛奶记录数据时获得了类似的结果。研究发现，热负荷对蛋白质百分比（PROT%）的影响非常微小，以至于光周期和热负荷的综合效应几乎与光周期单独的效应重叠（Aharoni等人，1999年）。总之，PROT%似乎受到一个与光周期变化相反的年度正弦模式的影响。除了这种长期效应外，还应考虑与特定月份相关的中期变化的影响。

在我们的研究中，脂肪百分比（FAT%）与PROT%呈高度正相关（+94%），与光周期（-96%）和平均温度（-87%）呈负相关。事实上，FAT%遵循与PROT%相同的季节性变化模式，但相位滞后3.09个月（表3），介于光周期（2.66个月，表2）和空气温度特征（3.52至3.59个月，表2）之间。在之前描述美国11个市场地区牛奶成分的研究中（Salfer等人，2019年），SIN-F仅在7个地区对FAT%有显著影响。在这些显著地区，振幅范围为0.08%至0.10%（我们的平均值为0.18%），峰值出现在12月27日至1月6日之间（在我们的案例中是在1月的第一周）。同样的作者在分析单个奶牛的牛奶记录数据时也得到了类似的结果。在以色列奶牛的研究中（Aharoni等人，1999年），光周期对FAT%的影响小于对蛋白质的影响，峰值出现在11月至12月，最低值出现在4月，即相位比本研究中的PROT%和FAT%提前约两个月。

在当前研究中，日均产奶量（dMY）从秋季到春分逐渐增加，在接下来的六个月里逐渐减少，直到下一个秋分，这种变化模式也在其他研究中有所报道（Becker等人，2020年）。这意味着，在我们的研究条件下，当光周期较短时dMY增加，当光周期延长速率最大时达到峰值，当光周期缩短速率最大时达到最低值。正如预期的那样，dMY与FAT%（-17%）和蛋白质（-38%）呈负相关，而与光周期呈正相关（+38%），与平均温度（Tmean）呈轻微负相关（-5%）。这种年度模式与其他主要采用集约化室内奶牛养殖系统的国家的大规模调查结果一致。Salfer等人（2019年）发现，多胎奶牛的dMY振幅（2.39公斤/天）大于初产和二胎奶牛（1.73-1.56公斤/天）。dMY的峰值出现在3月24日（多胎奶牛）至5月3日（初产奶牛）。Aharoni等人（1999年）也观察到dMY受光周期驱动的类似正弦模式，峰值出现在5月，最低值出现在11月，并且在7月至9月由于热负荷导致dMY下降。结合光周期和夏季热负荷的影响，他们观察到dMY的表型值在9月至10月最低，这与我们的研究结果一致。值得注意的是，他们发现峰值和最低值之间的差异为3.1公斤/天（我们的研究为3.3公斤/天），而夏季最大热负荷的影响为-1.8公斤/天。

Li等人（2023年）分析了美国不同州以及不同胎次奶牛的牛奶生产时间序列，发现大多数州的牛奶产量具有季节性模式，12月至3月产量最高，7月至9月最低。他们还发现，多胎奶牛的dMY季节性变化比初产奶牛更为明显。除了大量关于光周期对哺乳期和小型反刍动物繁殖性能影响的研究外，还有广泛文献表明光周期直接影响奶牛的激素谱以及促进乳腺组织细胞增殖和活性的机制，尤其是在干奶期（Dahl，2008年）。总之，dMY似乎始终受到一个年度正弦模式的影响，其相位比光周期的逆模式滞后约两个月，比空气温度特征的相位滞后一个月。

本研究的第三个目的是调查SIN-F是否可以应用于不同海拔或不同品种的奶牛，以检测对气候特征的响应变化，这可能表明气候敏感性存在差异。奶牛场的高海拔类别占可归因于个别农场变化的一小部分但非常显著的比例（表3）。位于丘陵和山区的农场比平原地区的农场具有较低的dMY、FAT%和PROT%（分别低11%、4%和1%）（表4）。考虑到所有因素，这些差异的幅度都很小，因为参与这项研究的威尼托地区的丘陵和山区农场主要是集约化室内奶牛场，这应在解释结果时予以考虑，同时还要考虑到丘陵和山区农场数量远少于平原农场。

表4. 海拔（平原 vs 丘陵-山区）和奶牛品种对平均牛奶特征及描述年度变化的SIN-F参数的影响（A = 振幅，以特征单位计；C = 相位滞后，以月计）。
| 项目 | 平原 | 丘陵-山区 |
|------------|------------|--------------|
| 日均产奶量（公斤/天） | 35.7 | 31.6 |
| 脂肪百分比（%） | 9.1 | 8.9 |
| 蛋白质百分比（%） | 3.7 | 3.6 |
| 描述性统计 | | |
| 平均值 | 35.7 | 31.6 |
| 标准差 | 9.1 | 8.9 |
| 海拔效应 | 0.13 | 0.71 |
| 品种效应 | 0.13 | 0.71 |
| Holstein | 36.5 | 36.5 |
| Brown Swiss | 31.4 | 31.4 |
| Simmental | 29.5 | 29.5 |
| 杂交奶牛 | 33.1 | 33.1 |

SIN-F参数：
| A | 1.42 | 1.47 |
| C | 0.40 | 1.12 |
| 海拔效应 | -0.18 | -0.17 |
| 品种效应 | -0.18 | -0.17 |
| Holstein | 1.47 | 1.47 |
| Brown Swiss | 1.12 | 1.12 |
| Simmental | 0.95 | 0.95 |
| 杂交奶牛 | 1.41 | 1.41 |

SIN-F很好地适应了两个地理区域的三个牛奶生产特征，在dMY的情况下，海拔类别与SIN-F的形状之间存在显著交互作用（表3），因为平原地区的农场振幅（+38%）大于丘陵和山区农场（表4）。农场海拔对SIN-F相位滞后的影响非常微小。威尼托平原的气候（海拔1至350米）根据K?ppen-Geiger分类法（Geiger，1954年）被归类为Cfa，属于湿润亚热带气候，夏季炎热潮湿，冬季寒冷温和；而丘陵和山区的气候（海拔350至1183米）被归类为Cfb，属于海洋性气候，夏季温和，冬季凉爽（表1）。平原地区的平均气温明显高于丘陵和山区（表2），且夏季的温差大于冬季，Tmax的温差大于Tmean或Tmin的温差（表1）。正如预期的那样，由于两个气候区位于相似的纬度，光周期没有差异（表1）。

SIN-F的参数（振幅和相位滞后）在表4中报告，并在图3中以图形方式展示。平原地区饲养的奶牛表现出更大的dMY季节性变化（图3a），峰值和最低值之间的差异分别为2.84公斤/天和2.06公斤/天。由于两个地区的光周期几乎相同，饲养和管理方式也非常相似，我们推测平原地区较高的平均温度（以及可能伴随的热应激）对dMY的振幅影响更大，或者平原地区的高产奶牛对温度变化的敏感性更强。PROT%的变化模式在平原和丘陵/山区非常相似（图3c）。FAT%的曲线在相位滞后方面有一些视觉差异，丘陵和山区的滞后接近一个月（图3b），但这种差异不显著，这证实了这一特征的较大未解释变化。

我们没有发现其他使用SIN-F比较平原和山区饲养的牛群牛奶生产特征长期年度变化的研究。然而，当比较不同纬度的地区时（Caraba?o等人，2016年；Rockett等人，2023b），结果与我们的研究有一些相似之处。Salfer等人（2020年）使用SIN-F（余弦节奏测量法）分析了美国四个重要奶牛区的单个奶牛的牛奶记录数据，包括两个北部州（宾夕法尼亚州和明尼苏达州）和两个南部州（佛罗里达州和德克萨斯州）。这两对美国州之间的最高温度差异大于我们研究中平原和山区农场之间的差异。由于这四个美国州位于不同的纬度但海拔相同，它们的光周期也不同（威尼托地区的光周期年度变化与明尼苏达州相似，且大于德克萨斯州和佛罗里达州）。dMY的SIN-F振幅在较温暖的南部地区（1.8-2.1公斤/天）大于较寒冷的北部地区（1.2公斤/天），这与我们研究中温暖气候地区（平原）的振幅大于丘陵和山区的情况一致。此外，所有美国和意大利牛群的相位滞后非常相似，峰值都在4月。

最近，Li等人（2023年）比较了美国17个州的初产、二胎及更多胎次的奶牛，将它们视为不同的群体（51个队列）。使用十年的时间序列数据来分离历史年度趋势和年度季节性及残差变化。没有对年度季节性变化做出任何假设。结果显示，大多数群体的哺乳期dMY的年度季节性变化具有明显的正弦模式，初产奶牛是主要例外。峰值通常出现在1月至5月期间，最低值出现在8月至9月期间。意大利和美国关于PROT%的结果也相似：四个州之间的振幅差异（0.06%至0.11%）与纬度无关。因此，比较不同地区的研究结果时，dMY的年度变化幅度主要由光周期驱动。然而，在像本研究这样的地区，当光周期相同时，如果哺乳期奶牛暴露在亚热带夏季和可能长时间的热应激条件下，振幅会更大。

与海拔类似，奶牛品种也显著影响所有牛奶生产特征（表4）。在牛奶产量方面，Holstein和Simmental的差异最大（后者比前者低17%）；在脂肪和蛋白质浓度方面，Brown Swiss和Holstein的差异也最大（后者分别低8%和9%）。这些差异可以归因于共同环境下的品种差异，因为大多数Brown Swiss、Simmental和杂交奶牛与Holstein一起在多品种牛群中饲养，且模型中考虑了各个牛群的影响。这些dMY的差异与之前在同一品种中观察到的差异相似（Stocco等人，2018年；Bittante等人，2021a年；Martínez-Marín等人，2023年），以及在比较荷斯坦奶牛与杂交奶牛时的差异（Malchiodi等人，2014年；Saha等人，2020年）。SIN-F似乎能够很好地描述不同品种奶牛的长期周期性变化，我们发现奶牛品种对SIN-F的参数有显著影响。如表4所示和图4a所示，两种阿尔卑斯起源的品种（西门塔尔和布朗瑞士）的dMY振幅较小，并且与荷斯坦奶牛和杂交奶牛相比存在大约一个月的相位延迟。下载：下载高分辨率图片（91KB）下载：下载全尺寸图片下载：下载高分辨率图片（84KB）下载：下载全尺寸图片下载：下载高分辨率图片（77KB）下载：下载全尺寸图片图4. 奶牛品种（荷斯坦，黑色曲线；布朗瑞士，棕色曲线；西门塔尔，红色曲线；杂交，黄色曲线）对牛奶特性平均周期性正弦函数的影响（a，日产量，蓝色实线曲线；b，脂肪含量，橙色实线曲线；c，蛋白质含量，绿色实线曲线）。a）：日产量（dMY）；b）：牛奶脂肪含量（FAT%）；c）：牛奶蛋白质含量（PROT%）。奶牛品种也影响了牛奶质量特性，这证实了之前对同一品种的研究结果（Stocco等人，2018年；Stocco等人，2019年；Piazza等人，2023年），特别是对于牛奶蛋白质（Maurmayr等人，2018年；Amalfitano等人，2020年）。品种对PROT%的影响显著但不明显（表4和图4c）。相反，布朗瑞士奶牛的FAT%振幅比荷斯坦奶牛和杂交奶牛更大（图4b），而西门塔尔奶牛的相位延迟接近一个月。Salfer等人（2020年）还比较了属于4个品种类型的奶牛的正弦函数：荷斯坦、杂交、杰西和布朗瑞士奶牛。除了杰西奶牛未包含在我们的数据库中外，我们研究中发现的dMY振幅与美国这三种品种的振幅完全相当。此外，相位延迟也受到品种的影响，极端情况是荷斯坦奶牛（4月16日达到峰值）和布朗瑞士奶牛（5月10日达到峰值）。因此，Salfer等人（2020年）的结果与我们的发现基本一致。我们不知道还有其他研究比较了不同品种奶牛的周期性正弦函数。然而，研究表明不同品种之间存在耐热性差异。在意大利，Maggiolino等人（2020年）通过dMY对THI的回归分析发现，布朗瑞士奶牛的耐热阈值较高，这意味着其耐热性更强。Landi等人（2023年）还发现布朗瑞士奶牛存在遗传与环境之间的交互作用。在半热带条件下，布朗瑞士奶牛和杂交奶牛在耐热方面似乎比荷斯坦奶牛表现更好，但dMY变化的差异不显著（Cuellar等人，2023年）。不同奶牛品种对气候变化的耐受性和恢复力水平似乎是由于遗传因素造成的。在这项研究中，不同品种的dMY存在差异（荷斯坦和杂交奶牛的dMY高于布朗瑞士和西门塔尔奶牛），因此不能排除较高的代谢负担对品种差异的部分影响。压力，无论是生理的、生物的还是心理的，可以定义为生物体对压力（如不利环境条件）的反应。这是一种自然反应，有助于生物体应对挑战和威胁（世界卫生组织，2023年）。多项研究证实了奶牛对气候条件的耐受性和恢复力具有品种内的遗传基础（Santana等人，2017年；Poppe等人，2021年；Rockett等人，2023a年）。我们认为这个问题与两个不同特性有关：奶牛与周期性气候变化的同步能力，尤其是在温暖的环境中，以及它们对短期热浪的反应能力。不幸的是，这两个特性在个体动物中难以区分或测量，因为在单个泌乳期中季节效应与泌乳阶段效应重叠，而在多胎奶牛中，这也与随后几年同一季节连续泌乳的频繁开始重叠，至少对于最健康的奶牛来说是这样。需要更多的研究来建立能够区分个体奶牛短期和长期效应的统计模型（Ji等人，2020年）。本研究最后一个目的是推断牛奶生产特性的中期（夏季）偏离SIN-F的情况，作为（热）压力影响的潜在指标。在包含月份之前的SIN-F的ANOVA Type I模型中，日历月份的影响应反映个别月份与SIN-F所代表的长期模式的偏差。根据这一假设，任何给定月份偏离零的显著偏差都可能表明周期性模式的中期改变。在dMY的情况下（图4a），5月和6月观察到显著变化（分别为+1.02 kg/d和+0.64 kg/d），7月则出现下降（-0.38 kg/d）。6月到7月dMY下降约1.0 kg/d可以被认为是由严重热应激引起的偏离。8月也观察到了类似的负偏差（-0.43 kg/d），但不显著（P=0.165），因为该月份的记录数据较少。环境温度高于（或低于）热中性区可能会成为压力源（Kingma等人，2012年），尤其是在高代谢负担的动物中（生长、怀孕、泌乳）。热应激可能由于长期周期性气候变化而在炎热潮湿的夏季地区影响泌乳奶牛，或者在较温和的气候条件下，由于中期和短期热浪而影响奶牛。值得注意的是，奶牛对高温反应延迟的一个重要原因可能是能量平衡。经常有报道指出，饲料摄入量减少是对环境温度升高的短期反应（Hou等人，2021年；Souza等人，2023年）。奶牛通过动员体脂储备来补偿这种饲料摄入量的减少（Cowley等人，2015年），但随之而来的体况恶化不能持续很长时间，导致dMY延迟下降。6月到7月牛奶产量突然下降1.73 kg/d，与正弦预期值（-0.71 kg/d）相比，可能反映了热应力的中期效应（图4a）。减少饲料摄入量，从而减少代谢产热，可能会降低dMY并延迟生育能力（Tao等人，2020年）。dMY短期和中期下降的另一个原因可能是热应激奶牛无法满足增加的饮水需求，这可能导致血细胞比容和血液尿素水平升高（Becker等人，2020年）。无论如何，这个问题很复杂，因为泌乳奶牛通过代谢适应来应对热应激。压力，无论是生理的、生物的还是心理的，都可以定义为生物体对压力的反应，例如来自不利环境条件的压力。这是帮助生物体应对挑战和威胁的自然反应（世界卫生组织，2023年）。多项研究证实了奶牛对气候条件的耐受性和恢复力具有品种内的遗传基础（Santana等人，2017年；Poppe等人，2021年；Rockett等人，2023a年）。我们认为这个复杂问题与两个不同特性有关：奶牛与周期性气候变化的同步能力，尤其是在温暖的环境中，以及它们对短期热浪的反应能力。不幸的是，这两个特性在个体动物中不易区分或测量，因为在单个泌乳期中季节效应与泌乳阶段效应重叠，在多胎奶牛中，这也与随后几年同一季节连续泌乳的频繁开始重叠。需要更多的研究来建立能够区分个体奶牛短期和长期效应的统计模型（Ji等人，2020年）。本研究的主要目的是推断牛奶生产特性的中期（夏季）偏离SIN-F的情况，作为（热）压力影响的潜在指标。在包括月份之前的SIN-F的ANOVA Type I模型中，日历月份的影响应反映个别月份与SIN-F所代表的长期模式的偏差。根据这一假设，任何给定月份的显著偏离零都可能表明周期性模式的中期改变。在dMY的情况下（图4a），5月和6月观察到显著变化（分别为+1.02 kg/d和+0.64 kg/d），7月则下降（-0.38 kg/d）。6月到7月dMY下降约1.0 kg/d可以被认为是由严重热应激引起的偏离。8月也观察到了类似的负偏差（-0.43 kg/d），但不显著（P=0.165），因为该月份的记录数据较少。环境温度高于（或低于）热中性区可能会成为压力源（Kingma等人，2012年），尤其是在高代谢负担的动物中（生长、怀孕、泌乳）。热应激可能由于长期周期性气候变化而在炎热潮湿的夏季地区影响泌乳奶牛，或者在较温和的气候条件下，由于中期和短期热浪而影响奶牛。值得注意的是，奶牛对高温反应延迟的一个重要原因可能是能量平衡。饲料摄入量减少经常被报道为对环境温度升高的短期反应（Hou等人，2021年；Souza等人，2023年）。这种饲料摄入量的减少通过动员体脂储备来补偿（Cowley等人，2015年），但随之而来的体况恶化不能持续很长时间，导致dMY延迟下降。6月到7月牛奶产量突然下降1.73 kg/d，与正弦预期值（-0.71 kg/d）相比，可能反映了热应力的中期效应（图4a）。减少饲料摄入量，从而减少代谢产热，可能会降低dMY并延迟生育能力（Tao等人，2020年）。dMY短期和中期下降的另一个原因可能是热应激奶牛无法满足增加的饮水需求，这可能导致血细胞比容和血液尿素水平升高（Becker等人，2020年）。无论如何，这个问题很复杂，因为泌乳奶牛通过代谢适应来应对热应激。在夏季月份（6月、7月和8月），牛奶蛋白质含量与周期性SIN-F预测的值有显著负偏差。这表明，除了光照周期的负面影响（6月达到最低值）外，最热月份的蛋白质含量进一步降低，表明热应力的影响可能是光照周期影响的叠加。值得注意的是，Chen等人（2024年）在其关于热应激对饲料摄入量、牛奶产量、牛奶成分和饲料效率影响的荟萃分析中发现，气候室中的短期热应激（28项研究）一致降低了蛋白质含量。相比之下，通过季节变化评估的中期热应激对蛋白质含量的影响变化很大（6项研究）。然而，解释中期偏差很复杂，因为1月和2月也观察到了其他负偏差，4月、5月和10月也有正偏差，尽管程度较小。FAT%与SIN-F模式的中期偏差更难以解释。与蛋白质含量不同，夏季月份没有显示出与长期正弦趋势的显著偏差。然而，9月观察到了明显的负偏差，因为牛奶脂肪含量没有像SIN-F预测的那样迅速恢复。这引发了是否长期抑制FAT%可以归因于热应激的残余效应的问题，可能延长了光照周期引起的下降。Aharoni等人（1999年）发现7月至9月期间FAT%的中期热应激效应比长期正弦变化更强。Chen等人（2024年）在对34项研究的荟萃分析中发现，热应激对FAT%没有显著的短期和中期变化。另一方面，热应力的中期效应能够改变奶牛的脂质代谢（Faylon等人，2015年）。周期性正弦建模作为一种工具，用于研究不同品种、农场和养殖系统的环境变化和热应激效应。关于热应激对奶牛影响的研究主要有两种类型。第一种是在气候室中进行的严格实验条件下的研究，旨在比较在大多数主要变异源（如光照、通风、饮食、管理）方面同质的奶牛群体，除了温度和湿度（Cowley等人，2015年；Gao等人，2017年；Hou等人，2021年）。这些研究允许准确测量生产、行为和生理特性，但通常持续时间较短，涉及的动物数量较少，且管理条件与实际农场相差甚远。此外，它们只能推断热应力的短期效应。第二种类型，如本研究，是比较不同环境和/或季节的真实农场（Salfer等人，2020年；Rockett等人，2023b）。气候参数通常无法直接获得，因此使用最近气象站的平均数据（Maggiolino等人，2020年）。在这种情况下，假设所有靠近同一气象站的农场具有相同的气候数据。当然，这种方法没有考虑到位于同一地区的不同农场可能在不同的微环境条件下饲养奶牛，实际的热应激负荷根据室外/室内管理和牛舍的不同特性而变化（Tresoldi等人，2018年；Becker和Stone，2020年；Drwencke等人，2020年；Reuscher等人，2023年；动物密度、粪便类型和管理方式等（Espinoza-Sandoval和Calsamiglia，2023年）。只有在比较具有不同气候条件的地区时，这种假设才可接受（纬度和/或海拔）。周期性正弦建模在牛奶生产特性方面的一个有趣应用可能是评估气候条件对不同品种、个体牛群或不同养殖系统的影响程度。这对于评估热应激耐受性和恢复力指标中的基因型×环境（G×E）相互作用非常有用。许多研究根据从气象站获得的数据对牛群风险水平进行分类，尽管众所周知，即使在相同的外部气候条件下，同一地区的农场也可能由于牛舍的特性而有非常不同的内部条件（Lambertz等人，2014年）。因此，在缺乏特定牛舍内气候数据的情况下（Gernand等人，2019年），一个有趣的方法是根据动物在SIN-F振幅和相位偏移方面的生产反应对农场进行分类。总之，SIN-F几乎描述了所有长期（周期性）变化，以及大部分观察到的空气温度特性的变化。对于牛奶特性，在本研究中采用的假设和统计建模下，周期性SIN-F也代表了主要的日历变化来源，超过了年际差异或中期每月偏离SIN-F的影响。牛奶蛋白质浓度的年周期变化与光照周期（白天）呈负相关，同时也受到中期热应激条件的影响。脂肪浓度部分受到气候条件的影响。每日产奶量受到年周期气候变化的影响，其滞后时间约为光照周期的3个月，以及空气温度最高值的2个月。位于平原地区（亚热带气候）的奶牛的年周期变化幅度比位于丘陵和山区（海洋性气候）的奶牛更大。阿尔卑斯品种（西门塔尔和布朗瑞士）的奶牛在每日产奶量方面表现出不同的年周期正弦波模式，尽管它们在同一地区和同一农场饲养。这些奶牛年周期变化幅度的较小可能部分反映了不同地区和品种之间的生产力差异，但也可能意味着它们具有更强的抵抗力，以及对年周期气候变化（特别是空气温度和可能的热应激）的更好适应性（即更长的相位滞后）。年周期SIN-F及其变化可以为涉及不同农业系统的研究提供工具，尽管需要进一步的研究来验证和加强这一假设。

**伦理批准**
不适用

**数据和模型可用性声明**
没有数据存放在官方存储库中，但可根据请求提供。

**关于写作过程中使用生成式AI和AI辅助技术的声明**
作者在写作过程中未使用任何人工智能辅助技术，仅在语言修订方面使用了相关工具。

**作者ORCID信息**
Giovanni Bittante - https://orcid.org/0000-0001-7137-7049
Gustavo Javier Martínez-Marín - https://orcid.org/0000-0002-4098-7685
Hugo Oswaldo Toledo-Alvarado - https://orcid.org/0000-0001-7854-1219
Nicolò Amalfitano - https://orcid.org/0000-0001-6030-1620
Alessio Cecchinato - https://orcid.org/0000-0003-3518-720X
Stefano Schiavon - https://orcid.org/0000-0002-5539-8947
Luigi Gallo - https://orcid.org/0000-0002-8908-5105

**利益声明**
作者声明没有已知的竞争性财务利益或个人关系可能影响本文所述的工作。

**未引用的参考文献**
Bittante等人，2023年；Bittante等人，2021b年；Martínez-Marín，2024年。

**CRedIT作者贡献声明**
G. Bittante：撰写初稿、监督、资源管理、项目协调、资金获取、概念构思。
G.J. Martínez-Marín：撰写初稿、数据可视化、软件使用、方法论设计、数据分析。
H. Toledo-Alvarado：撰写与编辑、监督、软件使用、方法论设计、研究实施。
N. Amalfitano：撰写与编辑、监督。
A. Cecchinato：撰写与编辑、监督。
S. Schiavon：撰写与编辑、监督。
L. Gallo：撰写与编辑、监督。