《Tropical Animal Health and Production》：Unsupervised clustering identifies candidate metabolic profiles, including an acidotic profile detectable before calving, in periparturient dairy cows
编辑推荐：
围产期（periparturient period）奶牛的代谢状态具有高度异质性，传统基于单一生物标志物的评估方法难以全面反映全身性代谢失调。本研究旨在应用无监督聚类（unsupervised clustering）方法，结合整合性血液生理指标识别不同的代谢表
围产期（periparturient period）奶牛的代谢状态具有高度异质性，传统基于单一生物标志物的评估方法难以全面反映全身性代谢失调。本研究旨在应用无监督聚类（unsupervised clustering）方法，结合整合性血液生理指标识别不同的代谢表型（metabolic phenotypes），并评估代谢受损表型在产前的可识别性。研究从中国七家商业牧场的316头健康荷斯坦奶牛（产犊前7天至产犊后7天，即?7至+7 d）采集静脉血样，进行血气与电解质分析。经数据清洗与变量筛选后，将六项关键指标——pH、PCO?、HCO??、葡萄糖（GLU）、钾（K?）和钠（Na?）——标准化后进行K-means聚类。最佳聚类数确定为3，得到以下表型：（1）"代谢受损"组，以酸中毒（低pH）、高碳酸血症（PCO?升高）、低血糖和电解质失衡为特征；（2）"正常酸碱-高钾"组，酸碱状态正常但钾持续升高；（3）"代谢稳态"组，表现为碱中毒倾向、较高的钠和葡萄糖水平以及稳定的电解质平衡。相当比例的奶牛在产犊前即已属于代谢受损表型，其患病率在产后中期（+4至+7 d）达到峰值。在未调整分析中，该聚类奶牛产后β-羟丁酸（BHBA）浓度轻度升高（P = 0.001），乳蛋白百分比显著降低（P < 0.001）；然而，在混合效应模型中调整牧场来源后，BHBA关联不再显著（P = 0.622），而乳蛋白关联仍然稳健。这些结果表明，该聚类中代谢适应受损更一致地反映为乳蛋白降低而非BHBA升高。基于聚类中心并经探索性ROC分析验证，提出识别代谢受损表型的产前阈值：静脉pH < 7.25（AUC = 0.928）、PCO? > 58.75 mmHg（AUC = 0.956）、葡萄糖 < 27.5 mg/dL（AUC = 0.876）。钾因AUC仅为0.577而未保留为判别阈值。这些阈值为探索性的，需在独立队列中验证。上述发现表明，围产期奶牛可划分为三种不同的代谢状态，且代谢受损的酸中毒表型可在产前被识别，为早期风险分层提供了基础；未来需开展干预性研究以确定针对该表型的干预能否改善临床结局。
围产期（periparturient period）是奶牛从产前3周至产后3周的关键生理过渡阶段，伴随深刻的代谢、内分泌和免疫适应以支持泌乳启动。此阶段乳合成对营养需求的指数增长常超过干物质采食量，导致负能量平衡（negative energy balance, NEB），触发广泛脂肪组织动员，循环非酯化脂肪酸（non-esterified fatty acid, NEFA）和β-羟丁酸（β-hydroxybutyrate, BHBA）升高，不仅易诱发酮病和脂肪肝等经典代谢病，还损害免疫功能和繁殖性能。传统代谢健康评估依赖单一生物标志物（如NEFA、BHBA、葡萄糖）或产后临床疾病诊断，虽有一定信息量，但这种还原论方法无法捕捉过渡期代谢失调的系统多维性。近年证据提示酸碱稳态、电解质平衡与能量代谢之间存在交互作用，严重NEB奶牛常出现血液pH、二氧化碳分压（partial pressure of carbon dioxide, PCO?）及K?、Na?等关键离子紊乱，与肝脂沉积和胰岛素抵抗密切相关。当前管理实践仍以被动应对为主，仅在临床症状出现后干预，而早期识别代谢受损个体可显著提高营养策略（如瘤胃保护性胆碱补充）的效力。近年来，无监督机器学习开始重塑奶牛代谢健康评估，K-means和生长混合模型已应用于血液代谢参数以识别潜在代谢表型，随机森林和支持向量机等监督算法在疾病预测中展现潜力。然而，现有研究存在三方面不足：既往聚类研究主要聚焦单一场内产后奶牛，外部效度有限；表型分型多依赖孤立生物标志物类别，未整合多维生理维度；最重要的是，多数研究使用产后数据定义代谢组，限制了产前主动干预的可行性。
针对上述空白，研究人员开展了一项前瞻性观察性横断面研究，于2023年9月至2024年3月在中国宁夏银川平原七家商业奶牛场进行。研究共纳入328头荷斯坦奶牛，最终分析数据集为316头奶牛的3,731条观测值（每头牛贡献一次血样）。纳入标准为：约产前60 d经直肠触诊确认单胎妊娠、采样窗口内无临床疾病征象、采样前14 d内无治疗干预。奶牛按胎次分为初产（n = 112）和经产（n = 216）。于产犊前7天至产后7天经颈静脉穿刺采集血样于锂肝素真空管，冰上运输至场内实验室，2 h内用便携式血气分析仪（i-STAT Alinity）完成分析。产后7 d内采集晨乳样（平均3.2 ± 1.8 d），用中红外光谱测定乳脂和乳蛋白百分比，同日记录日产奶量。产后7 d内测定血BHBA浓度。
数据预处理按以下顺序进行：排除缺失值样本（n = 12）；以四分位距（interquartile range, IQR）法在原始数据上剔除极端离群值（n = 4）；共线性评估与变量筛选；对保留变量进行Z-score标准化；K-means聚类。通过Pearson相关系数和方差膨胀因子（variance inflation factor, VIF）评估共线性，最终选取六项变量——pH、PCO?、HCO??、GLU、K?、Na?——涵盖酸碱平衡、能量代谢和电解质稳态三个核心维度。HCO??与总CO?（r = 0.933）及碱剩余（r = 0.880）高度相关，保留HCO??而排除后两者；血细胞比容与血红蛋白高度相关（r = 0.960）且主要反映血液浓度而非核心代谢过程，予以排除；离子钙因与pH和HCO??共线性（|r| > 0.4）且低钙血症发生率较低而排除。最终六项变量VIF均< 2.5，最大条件指数为8.2，无显著共线性。
K-means聚类选择基于以下考虑：六项聚类变量为连续标准化变量，代表不同生理维度；分析数据集（n = 3,731）远超K-means推荐最小样本量；K-means在奶牛代谢表型研究中已有广泛应用先例；聚类中心直接对应平均生理轮廓，便于转化为探索性阈值。最佳聚类数通过肘部法（基于组内平方和WCSS）、平均轮廓宽度和生物学可解释性综合确定。肘部法在K = 2–3处显示拐点，平均轮廓宽度在K = 3时最高，间隙统计量在K = 2时峰值。鉴于轮廓分析直接评估聚类凝聚性和分离性且K = 3产生生物学可解释的表型，最终选择K = 3。共识聚类（1,000次bootstrap重采样，每次抽取80%样本不放回）验证三聚类解决方案高度稳定（平均Jaccard = 0.975；调整Rand指数ARI = 0.951）。内部验证采用Calinski–Harabasz指数（1,087.44）、Davies–Bouldin指数（1.448）、Dunn指数（0.010）和平均轮廓宽度（0.231）。随机森林分类器（1,000棵树）以聚类标签为结局、六项标准化变量为预测因子，袋外（out-of-bag, OOB）准确率达98.1%。留一变量敏感性分析逐一移除聚类变量后重新运行K-means并计算ARI。统计分析方法包括：连续变量组间比较采用Student's t检验或Mann–Whitney U检验，分类变量采用卡方检验；各聚类在各围产期时间窗（产前?7至?1 d、产后早期0至+3 d、产后中期+4至+7 d）的占比采用Cochran–Armitage趋势检验评估时间趋势；采样日与聚类归属的关联采用多项 logistic 回归（以牧场和胎次为协变量）；BHBA与聚类的关联采用Kruskal–Wallis检验及线性混合效应模型（牧场为随机效应）；候选阈值的判别性能采用探索性ROC分析。显著性水平设为P < 0.05，所有分析在R中进行。
研究结果显示以下主要发现：
**产前与产后血液参数的系统性差异**：产后奶牛在PCO?、HCO??、总CO?、碱剩余、Na?、K?、离子钙、血细胞比容和血红蛋白水平上与产前奶牛存在显著差异（P < 0.05或P < 0.001），其中PCO?、HCO??、总CO?、碱剩余、Na?、K?和离子钙差异极显著（P < 0.001）；pH、氧分压、氧饱和度和GLU水平无显著差异（P > 0.05）。这表明产后奶牛酸碱平衡、电解质代谢、能量状态和血液学特征发生系统性变化。
**变量筛选与主成分分析**：Pearson相关分析显示HCO??与总CO?（r = 0.933）、碱剩余与HCO??（r = 0.880）高度相关，保留HCO??；血细胞比容与血红蛋白高度相关（r = 0.960），均排除。最终六项指标覆盖酸碱平衡（pH、PCO?、HCO??）、能量代谢（GLU）和电解质稳态（K?、Na?）三个维度，彼此无显著共线性。PCA显示第一、第二、第三主成分方差解释率分别为40.3%、21.7%、17.8%，前三个主成分累积解释79.8%。PC1主要由PCO?（正载荷）及pH、Na?、GLU（负载荷）驱动，代表呼吸性酸中毒与碱中毒及能量/电解质状态的对比轴；PC2主要由HCO??、PCO?、K?（正载荷）定义，反映酸碱调节的肾脏和电解质组分。PCO?和pH对整体PCA结构贡献最强（贡献 > 0.45），为核心判别因子。样本在PC1–PC2平面呈连续椭圆云分布，无明显离散亚组，少数离群样本位于PC1极端值处。pH与PCO?强负相关（r = ?0.79），Na与K中度负相关（r = ?0.48），PCO?与GLU中度负相关（r = ?0.44），所有指标间无强正相关（r > 0.8），证实变量非冗余。
**最佳聚类数确定**：肘部法在K = 2–3处显示拐点；轮廓分析在K = 3时平均轮廓宽度最大（约0.36）；间隙统计量在K = 2时峰值。鉴于轮廓分析被广泛认为是最可靠的聚类质量评估方法且明确支持K = 3，最终选择K = 3。PC1–PC2空间可视化显示三个聚类具有清晰空间分离：Cluster 1（红色）主要位于右半平面（PC1 > 0），沿PC2广泛分布；Cluster 2（绿色）紧密聚集于左下象限（PC1 < 0，PC2 < 0）；Cluster 3（蓝色）集中于左上区域（PC1 < 0，PC2 > 0）。边界附近尤其是原点周围存在轻微重叠，但整体判别良好。
**聚类稳定性评估**：共识聚类1,000次bootstrap重采样显示三聚类解决方案高度稳定，各聚类Jaccard系数分别为0.971、0.989、0.965（平均Jaccard = 0.975），原始与共识聚类归属间ARI为0.951。共识矩阵呈现清晰块对角结构，组内共归属概率高、组间低，表明三种代谢表型可重复且非随机初始化或数据集特异性噪声所致。
**采样日与聚类归属的关联**：多项 logistic 回归显示采样日与Cluster 1归属显著相关（RRR = 1.18，95% CI：1.09–1.28，P < 0.001），与Cluster 2无显著关联（RRR = 1.02，95% CI：0.96–1.09，P = 0.52）。似然比检验显著（χ2 = 16.84，df = 2，P < 0.001）。模型预测概率显示Cluster 1归属概率从?7 d约5%升至+7 d约23%，Cluster 2稳定在约27–28%，Cluster 3从约74%降至53%。效应量适中：每增加一天产后，Cluster 1相对风险增加约18%，从?7 d至+7 d累积增加约2.3倍；鉴于Cluster 1患病率从产前4.4%升至产后中期23.5%（5.3倍差异），采样日仅解释部分变异。经产奶牛归入Cluster 1的风险为初产奶牛的3.42倍（RRR = 3.42，95% CI：1.87–6.25，P < 0.001），与经产牛对围产期代谢紊乱易感性更高的已知规律一致。牧场ID亦显著贡献（χ2 = 47.23，df = 12，P < 0.001），反映管理和营养的群间差异。
**三种代谢表型的特征**：Cluster 1以显著低pH（7.18）、显著升高PCO?（74.53 mmHg）、较高K?（8.60 mmol/L）和最低GLU浓度（20.54 mg/dL）为特征；Cluster 2显示中间pH（7.297）和PCO?（51.22 mmHg），但K?持续高（8.89 mmol/L），GLU中等（40.20 mg/dL）；Cluster 3显示最高pH（7.341）、正常PCO?（51.33 mmHg）、显著升高Na?（139.76 mmol/L）和GLU（51.32 mg/dL），K?在参考范围内（6.78 mmol/L）。单因素方差分析显示六项变量在三聚类间均存在显著差异（P < 0.001），Tukey HSD事后检验证实Cluster 1的pH显著低于、PCO?和K?显著高于Cluster 2和Cluster 3（P < 0.001），Cluster 3的Na?和GLU显著高于其他两簇（P < 0.001），Cluster 2与Cluster 3在pH和PCO?上无显著差异但电解质轮廓不同。生物学上，Cluster 1代表以呼吸性酸中毒（低pH、高PCO?）合并高钾血症和低血糖为特征的代谢受损状态，可能反映通气受损、能量缺乏和围产期代谢应激；Cluster 2反映正常酸碱表型伴孤立性高钾血症，可能提示肾脏或激素调节异常而无明显酸碱失衡；Cluster 3显示碱中毒-高钠-高糖轮廓，可能代表应激诱导糖异生和钠潴留驱动的代偿状态。基于聚类中心并经探索性ROC分析验证，Cluster 1的产前阈值为：静脉pH < 7.25、PCO? > 58.75 mmHg、葡萄糖 < 27.5 mg/dL；钾因AUC仅为0.577未保留为判别阈值。
**代谢表型分布的动态变化**：Cluster 3在产前晚期占优势（57.9%），产后下降（早期48.0%，中期37.6%）；Cluster 1从产前4.4%逐步升至产后中期23.5%；Cluster 2在各阶段保持稳定（约38–40%），提示持续性钾调节异常独立于急性酸碱状态。
**代谢聚类与产后血BHBA的关联**：未调整分析中，Cluster 1奶牛BHBA浓度轻度高于Cluster 2和Cluster 3（中位数：1.4 vs. 1.0 vs. 0.9 mmol/L；P = 0.001）。以牧场为随机效应的线性混合效应模型调整后，聚类对BHBA的效应不再显著（F = 0.48，P = 0.622）；对数转换BHBA敏感性分析结果一致（F = 0.31，P = 0.730）。跨牧场聚类BHBA均值检查显示，未调整分析中的表观关联主要由单一牧场（Farm1）少数Cluster 1奶牛驱动，其余五场中Cluster 1未显示BHBA升高。这表明聚类归属与BHBA的表观关联主要归因于群间变异，BHBA不能可靠区分三种代谢表型。
**乳成分与产量**：乳成分在各代谢聚类间差异显著（P < 0.001）。Cluster 3乳脂率最高（5.8 ± 0.9%），Cluster 2次之（5.2 ± 0.8%），Cluster 1最低（4.6 ± 0.7%）；乳蛋白率呈类似梯度：Cluster 3最高（3.9 ± 0.3%），Cluster 2居中（3.7 ± 0.3%），Cluster 1最低（3.3 ± 0.4%；P < 0.001）。日产奶量各聚类间无显著差异（P = 0.091），Cluster 2产奶量趋于较低（32.1 ± 6.5 kg）vs. Cluster 1（34.2 ± 7.1 kg）和Cluster 3（34.8 ± 6.9 kg）。脂蛋比各聚类间无显著变异（P = 0.21），表明乳蛋白合成减少而非乳脂分泌改变是Cluster 1严重代谢紊乱的主要成分后果。
**探索性ROC分析**：以Cluster 1 vs. Cluster 2 + Cluster 3为二分类结局，AUC值分别为pH 0.928（95% CI：0.920–0.936）、PCO? 0.956（95% CI：0.950–0.962）、GLU 0.876（95% CI：0.866–0.886）、K? 0.577（95% CI：0.561–0.593）。Youden最优阈值为pH 7.255（灵敏度0.865，特异度0.838）、PCO? 58.75 mmHg（灵敏度0.934，特异度0.853）、GLU 27.5 mg/dL（灵敏度0.932，特异度0.759）、K? 6.85 mmol/L（灵敏度0.971，特异度0.183）。pH、PCO?和GLU显示良好至优秀判别性能，K?判别性能差（AUC = 0.577），灵敏度-特异度权衡无临床意义。四变量联合多变量 logistic 模型出现准完全分离，无法稳定最大似然估计，故仅报告单变量ROC结果。
**内部验证、监督分类与敏感性分析**：三聚类解决方案Calinski–Harabasz指数为1,087.44，Davies–Bouldin指数为1.448，Dunn指数为0.010，平均轮廓宽度为0.231。Calinski–Harabasz和Davies–Bouldin指数支持良好聚类分离；Dunn指数较低在高维数据中不意外；平均轮廓宽度0.231表明弱至中等聚类结构，与存在部分边界样本一致，但整体分离得到高度稳定共识聚类结果支持。随机森林分类器OOB准确率98.1%，各类错误率分别为Cluster 1 1.61%、Cluster 2 1.66%、Cluster 3 3.45%。变量重要性分析识别K?（MeanDecreaseAccuracy = 260.9）、GLU（176.4）和PCO?（120.7）为最具信息量预测因子，HCO??重要性最低（62.5），与其与pH和PCO?的部分冗余一致。留一变量敏感性分析显示移除Na?（ARI = 0.897）或HCO??（ARI = 0.868）对聚类归属影响极小，移除pH（ARI = 0.482）、GLU（ARI = 0.593）或PCO?（ARI = 0.619）导致中度变化，移除K?影响最大（ARI = 0.329），表明K?是维持三聚类结构最具影响力的单一变量。
讨论部分指出，本研究通过多牧场围产期奶牛整合血气与电解质轮廓的无监督聚类识别了三种不同代谢表型。Cluster 1代表代谢受损状态，以低pH、PCO?升高、HCO??降低、低血糖和电解质失衡为特征，与呼吸性酸中毒合并能量代谢受损一致；Cluster 2显示正常酸碱但高钾模式，提示孤立性电解质调节异常而无明显酸碱紊乱；Cluster 3显示稳态轮廓，以碱中毒倾向、钠和葡萄糖升高及稳定酸碱平衡为特征。值得注意的是，相当比例归入代谢受损表型（Cluster 1）的奶牛在产前采样，表明显著代谢失调可在产犊前、远早于临床疾病发生时被检测。Cluster 1定量阈值的确定为场内风险分层提供了实用框架，便携式血气分析仪可在数分钟内完成这些参数的快速评估。产前酸中毒和低血糖表型的识别与既往观察一致并有所扩展：妊娠晚期亚临床酸中毒与产后真胃变位和酮病风险增加相关；低pH和高PCO?反映的酸碱平衡紊乱与围产期奶牛肝脂沉积和胰岛素抵抗发展密切相关。Cluster 1中低血糖尤为值得注意：严重NEB下过度脂解导致NEFA流向肝脏增加，促进生酮并可能损害糖异生能力；酸中毒还抑制肝葡萄糖生成并减少外周葡萄糖利用，形成加剧能量缺乏的恶性循环。Cluster 1中观察到的髙钾血症可能反映细胞应激或胰岛素缺乏，胰岛素正常促进K?进入细胞，NEB期间其减少可能导致血清K?升高。酸中毒、低血糖和高钾血症在Cluster 1中的共存反映代谢稳态的系统性失败并具有直接病理生理后果：酸中毒损害肝糖异生并减少外周葡萄糖利用，加剧能量缺乏并促进脂解；同时严重NEB下常见的胰岛素缺乏减少细胞钾摄取，导致高钾血症；这一三联征形成恶性循环：低血糖进一步抑制胰岛素分泌，酸中毒加重电解质紊乱。临床上该状态使奶牛易患酮病、真胃变位和免疫抑制。尽管未调整分析显示Cluster 1 BHBA浓度较高，但该关联对牧场来源调整不稳健，且表观关联主要由单一牧场少数Cluster 1奶牛驱动，故BHBA不应被视为本队列三种代谢表型的可靠判别指标。关键的是，数据表明该代谢受损表型可在产犊前7天被识别，为抢先管理提供了窗口。正常酸碱-高钾轮廓（Cluster 2）可能源于不同机制，如醛固酮不足或肾小管功能障碍，不直接损害能量代谢但可能通过减少采食量或改变体液平衡影响泌乳性能。本队列中乳蛋白百分比降低成为产前代谢应激更具判别力的指标，与Cluster 1奶牛低血糖状态限制乳腺氨基酸摄取一致。约30%产前样本已属于代谢受损表型，挑战了代谢危机仅在产犊后出现的传统观点，提示围产期疾病根源在于妊娠晚期，可能由过度肥胖、日粮阴阳离子差管理不当或慢性亚急性瘤胃酸中毒等因素驱动。本研究存在若干局限性：分析基于每头牛单次采血，仅提供动态围产期生理的横断面快照，需纵向分析表征表型稳定性和转换；聚类未纳入胎次和牧场为协变量；未追踪临床结局，限制了对Cluster 1预测价值的直接验证；聚类框架限于六项可现场检测变量，未包括NEFA、氯、乳酸或镁，NEFA尤为NEB cornerstone生物标志物，其纳入可能优化代谢受损表型特征描述或揭示额外亚表型；离子钙虽测量但因共线性排除，故当前表型未明确捕捉亚临床低钙血症；横断面设计无法表征个体内表型时间变化，也无法确定个体奶牛是否在三个聚类间转换，观察到的表型分布变化描述的是独立采样奶牛的人群水平模式而非同一动物轨迹，该设计也排除了评估产前酸中毒轮廓与产后结局（如BHBA升高或乳蛋白降低）时间序列的可能性，需未来纵向研究从干奶后期至泌乳早期系列血气电解质测量确认产前代谢受损轮廓是稳定个体特征、短暂状态还是进行性恶化，并建立早期干预最佳采样窗口。关于变量选择，虽用VIF和条件指数正式评估共线性，但基于相关和VIF的标准不等同于因果或生物学变量选择，某些因统计原因排除的变量可能携带独立生理信息，在更大样本或替代建模方法中可优化表型特征描述。关于数据预处理，离群值在标准化前原始尺度上移除以防极端值扭曲Z-score计算所用均值和标准差，但这也意味着最终标准化数据集以离群值移除步骤为条件，若应用不同离群标准，标准化值和聚类解决方案可能不同；使用替代离群阈值（1.5×IQR和2×IQR）的敏感性分析产生可比三聚类解决方案，支持本方法稳健性。留一变量敏感性分析显示聚类解决方案对特征移除部分敏感，K?和pH影响最大（移除时ARI分别为0.329和0.482），提示三聚类结构主要依赖这两项变量，未来更广泛生物标志物面板研究应检验纳入额外变量后表型是否保持稳定。
结论部分翻译：围产期奶牛基于整合血气与电解质参数可划分为三种候选代谢轮廓。代谢受损轮廓（Cluster 1）以产前酸中毒、低血糖和高钾血症为特征，可在产犊前一周被识别。该轮廓与乳蛋白百分比降低相关，而其与产后β-羟丁酸升高的关联对群水平效应调整不稳健，提示乳蛋白降低是该表型代谢适应受损更一致的指标。基于聚类中心和探索性ROC分析，提出以下产前阈值用于识别代谢受损轮廓：静脉pH < 7.25（AUC = 0.928）、PCO? > 58.75 mmHg（AUC = 0.956）、葡萄糖 < 27.5 mg/dL（AUC = 0.876）。钾因AUC仅为0.577未保留为判别阈值。这些阈值为探索性的，在作为诊断或管理工具实施前需在独立队列中验证。未来还需干预性研究确定针对该轮廓的干预能否改善过渡期奶牛的临床结局。
