Different reticuloruminal pH metrics of high-yielding dairy cattle during the transition period in relation to metabolic health, activity, and feed intake

The measurement of pH in the reticulorumen in combination with a time–pH threshold has been widely applied in research to diagnose subacute ruminal acidosis. However, other pH metrics also have biological value. In this study, 44 animals were monitored during the transition period using reticuloruminal pH boluses. Traditional and more complex pH characteristics were calculated to characterize the reticuloruminal pH profile: time pH <6, slope of a logistic cumulative pH curve (β 0 ), and deviations [squared error (SqEr)] from pH predictions based on a harmonic static model. In this study, we aimed to examine the associations between those pH metrics and metabolic health parameters, feed intake, and activity. Finally, to describe the reticuloru-minal pH dynamically, we also constructed a dynamic linear model. The results of this model were studied in relation to feed intake. All pH parameters were mutually correlated (particularly β 0 and SqEr; mean Pearson correlation of −0.52). pH patterns, rather than time pH <6, were associated with metabolic health and feed intake: high variation in daily pH (β 0 parameter) was reflected in higher blood concentrations of nonesterified fatty acids. Moreover, pH deviations of the harmonic model were negatively associated with feed intake and rumination behavior. This research confirms the biological importance of pH metrics focusing on pH variation and pH deviations and provides deeper insight into its associations with metabolic health status, feed intake, and activity during early lactation.


INTRODUCTION
The early lactation period of dairy cows is often characterized by feed intake depression in combination with increasing demands of energy for milk production (Ingvartsen, 2006).Hence, the energy density of the diet, particularly the supply of glucogenic precursors (van Knegsel et al., 2007), is optimized to maximally alleviate body fat mobilization (Grummer, 2008).The disproportion between a high share of rapidly fermented carbohydrates and low physically effective NDF may cause impaired ruminal health through variation in VFA concentrations and decreased ruminal fluid pH.Accordingly, early lactating cows are, as a group, considered at risk of impaired ruminal health.Jing et al. (2018) reported large interanimal variation in rumen pH characteristics in early lactation, based on SARA-indicative milk fatty acids as indirect biomarkers.Additionally, milk (fat) depression and feed intake alterations are often observed during early lactation, which also widely varies between and within herds.In trials inducing SARA, variation in such responses have been associated with changes in rumen pH (Kleen et al., 2009;Kleen and Cannizzo, 2012;Abdela, 2016), which were also associated with differences in metabolic health status (Comino et al., 2015) and rumination behavior (Khiaosa-ard et al., 2018).However, it is unclear whether interanimal variation in performance, intake, and behavioral parameters in early lactation is associated with variation in rumen pH as detailed studies monitoring rumen pH characteristics in early lactating cows under common dairy farming practice are lacking.
Classically, variation in rumen pH is described by the time that rumen pH is below a certain threshold (Abdela, 2016), although threshold levels vary between studies and monitoring methods (Zebeli et al., 2008;Kleen and Cannizzo, 2012).Nevertheless, rumen pH is characterized by important within-day variation.This pH variation can be described by the use of a logistic curve (AlZahal et al., 2007).The slope of this curve (β 0 ) gives an indication of the stability of pH.Another approach to characterize this variation is by modeling normal pH patterns and calculating deviations from this model, as applied by Denwood et al. (2018).This variation is thought to contain biologically relevant information, which could be more informative compared with global time periods during which pH is below a certain threshold.Denwood et al. (2018) showed that abnormalities in pH patterns, rather than extreme absolute pH values, were strongly associated with feed intake and milk production 2 d after pH observation.Alternatively, Colman et al. (2012) used the β 0 parameter of a cumulative logistic curve to describe the variation in rumen pH throughout the day and discovered associations of diurnal pH variation with changes in milk fatty acid composition.High pH variation was associated with high C18:1 trans-10 concentrations in milk, which are often associated with milk fat depression (Colman et al., 2012;Dewanckele et al., 2020).
The hypothesis of this research was as follows: during early lactation, variation in rumen pH and deviations from pH patterns are linked to metabolic health status, feed intake, behavior, and milk production.We also hypothesized that these pH metrics have biologically added value during the transition period compared with pH-time metrics.Accordingly, the objective in this study is to classify cows using different pH metrics (absolute pH values, pH variation, and pH patterns) and to relate production characteristics, behavior, and metabolic health parameters to these pH metrics.

Animals and Management
The experiment was performed at the research farm of ILVO (Flanders Research Institute for Agriculture, Fisheries, and Food, Melle, Belgium) from May 2019 until October 2020.The experiment (2018/329) was approved by the Ethical Committee of ILVO.Only multiparous Holstein-Friesian cows (n = 44) were monitored in this study.The cows were randomly selected from the calvings of a larger herd during the 1.5-yr experiment.The number of animals was determined based on a SARA prevalence of 20% (Plaizier et al., 2008) under practical farming conditions, which requires more animals than SARA induction experiments.Reticuloruminal pH boluses were inserted 17 ± 4.6 d before expected calving, and pH was monitored up to 3 wk after calving.The boluses were inserted using an oral balling gun.These boluses gravimetrically end up in the reticulum, as verified by Villot et al. (2018).Due to issues with data retrieval, 6 cows did not have complete data up to 3 wk after calving.
Dry cows and lactating cows were housed separately in a naturally ventilated freestall barn with a slatted floor.Stocking density was always <1 cow/stall.From the time of imminent calving (e.g., pelvic ligament relaxation, teat filling) until d 3 after calving, cows were housed in maternity pens with straw bedding within the same building.In case of disease, cows stayed in these pens for a longer period (n = 4 cows).
The rations were formulated without SARA induction but following common practice on dairy farms.From 3 wk before calving, cows received the partial mixed ration for lactating cows supplemented with a dry cow mineral premix and, on average, 1 kg of balanced concentrate per cow per day.Belgian-Dutch energy and protein evaluation systems were used: requirements and supply of protein digestible at the level of the small intestine were assessed according to the DVE system (Van Duinkerken et al., 2010), and net energy requirements and supply were assessed according to the VEM system (Van Es, 1975).The partial mixed ration for lactating cows was calculated to fulfill the needs of an average adult cow of 650 kg of BW and producing 26 kg of fat-and protein-corrected milk, and was based on maize silage, prewilted grass silage, pressed beet pulp, and soybean meal to balance digestible protein (DVE) and net energy (VEM) requirements and supplemented with a VEM-DVE-balanced concentrate.The supply of balanced concentrate changed according to lactation stage and changed slightly during the course of the experiment (running over 1.5 yr) in relation to the quality and feed values of the silages used.Concentrate intake started at 1.7 kg of balanced concentrates, 0.2 kg of formaldehyde-treated soybean meal (CovaSoy Braz., FeedValid), and 0.3 kg of soybean meal at d 3 after calving.CovaSoy was increased to 1 kg over a period of 7 d, and the balanced concentrate was increased to 6 kg over a period of 20 d.Detailed information about the ingredients, chemical composition, and concentrate build-up is given in Supplemental Tables S1 and S2 (https: / / doi .org/ 10 .6084/ m9 .figshare .19336457;Heirbaut et al., 2022).Individual feed intake was monitored throughout the trial using roughage intake control (RIC) feeding bins (Insentec, Hokofarm Group), except during the period around calving.Cows had ad libitum access to the RIC bins.During lactation, concentrate intake was monitored at the automatic concentrate providers (Greenfeed, C-Lock Inc.; DeLaval) and in the herringbone milking parlor (DeLaval).Cows had continuous access to the automatic concentrate providers.The amount of concentrates provided depended on lactation stage (Supplemental Table S2; https: / / doi .org/  10 .6084/ m9 .figshare .19336457;Heirbaut et al., 2022).The cows had access to water ad libitum.

Blood and Milk Sampling
Blood sampling took place at d 7 before expected calving and on d 3, 6, 9, and 21 after calving using a 18-gauge needle and the Venoject system (Terumo).Samples were taken from the coccygeal vessels (d −7, 3, 6, and 9) or the jugular vein (d 21) in the morning between 0915 and 0945 h.The jugular vein sampling at d 21 was required to collect additional blood samples for another experiment.However, blood sampling location was not expected to have a confounding effect for the blood parameters analyzed in the current experiment (Mahrt et al., 2014).Blood was collected in plasma NaF (4 mL) and serum blood tubes (10 mL) (SSTTM II Advance tubes; BD Diagnostics).Serum tubes were kept at room temperature for 30 min before centrifuging and NaF tubes were kept cooled until centrifuging.The NaF tubes were centrifuged at 1,000 × g for 10 min and serum tubes at 1,500 × g for 15 min at room temperature.Serum and plasma samples were divided into aliquots and stored at −20°C [serum analyses: BHB, nonesterified fatty acids (NEFA), insulin and fructosamine; plasma analysis (NaF): glucose] or −80°C (serum analysis: IGF-1).Concentrations of BHB, NEFA, and glucose were analyzed using a Gallery Discrete Analyzer (Thermo Fisher Scientific) and Randox kits (Randox Laboratories Ltd.) by the laboratory of DGZ (Torhout, Belgium).Serum IGF-1 was analyzed by a radioimmunoassay method using a nonextraction IGF-1 IRMA DSL-2800 kit (LifeSpan Biosciences) by Poznań University of Life Sciences (Poznań, Poland).The concentration of IGF-1 was determined with isotope J-125 using an automatic gamma radiation reader (Wizard2 2-Detector Gamma Counter, Perkin Elmer).In this method, the peptide to be determined is sandwiched between 2 antibodies.The tubes were coated by the first antibody.The second antibody [anti-IGF-1 (J-125) reagent] was radiolabeled for detection.The analytic peptide (IGF-1) present in blood serum, standards, and controls was bound by both of the antibodies to form a sandwich complex.Unbound compounds were removed by decanting.Fructosamine was determined using the nitro blue tetrazolium method (Johnson et al., 1983;Jensen et al., 1993) with modifications described by Kishabongo et al. (2014) and Megahed et al. (2018).
Dairy cows were milked twice daily, at 0530 and 1630 h, in a 2 × 7 herringbone milking parlor, and their milk yield (kg/d) was recorded electronically.To determine milk performance, milk samples (27 mL) were collected from the cows in a representative way during the morning milking, on a daily basis from d 3 to 23 after calving.Samples were stored at 4°C and contained preservatives (sodium azide, maximum concentration 0.02% wt/wt, and bronopol, maximum concentration 0.005% wt/wt).The milk fat and protein contents were analyzed at Qlip (Zutphen, the Netherlands) by means of Fourier transform infrared spectrometry (Milkoscan FT6000, Foss Electric).

Data Collection
In total, data from 44 cows were processed.Two types of boluses were used.First, the eCow bolus (eCow Ltd.) was used (26 cows).However, because of cumbersome and time-consuming data retrieval and data losses, this bolus type was abandoned.Hence, the data set contained 24 cows with eCow bolus data.We switched to smaXtec boluses (smaXtec Animal Care GmbH; n = 20).The boluses monitored every 15 and 10 min, respectively.The eCow boluses were, according to the manufacturer guidelines, (1) heated in a water bath to between 35 and 40°C, and (2) calibrated against pH 4 and 7 standard buffers (VWR) before insertion.Data were downloaded periodically using a handset with dongle and antenna, transmitted to the eCow servers, and later downloaded.SmaXtec boluses were calibrated in a buffer of pH 7, provided by the manufacturer, before insertion.Boluses were inserted with an oral balling gun.Data were downloaded continuously using a smaXtec base station and 2 smaXtec repeaters in the stable.Data were transmitted to smaXtec servers and downloaded later.To avoid drift, all boluses were used for a maximum of 5 wk, which is below the 50-d accuracy warranty period of smaXtec and the 60-d period of eCow.Drift normally occurs after longer periods of data logging.
The cows' walking activity was monitored using IceTag3D motion sensors (IceRobotics).The IceTag3D motion sensor is a triaxial activity sensor (16 Hz).Sensors were attached to the right hind leg of the cow.IceTagAnalyser 2010 software (IceRobotics) was used together with an IceReader wireless download unit for activating and downloading the raw data.IceTags were attached around 4 wk before expected calving and removed after d 35 in lactation.In between, the IceTags were downloaded at least once.After downloading, data were aggregated at the 24-h level and converted to comma-separated value (CSV) files.Data from the day of attachment were removed from analysis.
RumiWatch halter sensors (Itin and Hoch GmbH) were used to monitor eating and rumination behavior.Data recording occurred at a frequency of 10 Hz, and data were saved on a memory card.Raw data were retrieved using RumiWatch Manager 2 software (version 2.2.0.0).Conversion of raw data into different behav-ior parameters was done using RumiWatch Converter software (version 0.7.3.2).Data were converted to 24-h resolution and stored in CSV files for further analysis.Data from the day of attachment were removed from analysis.
pH Data.Data from eCow and smaXtec boluses were stored in .csvand .xlsx(Excel; Microsoft Corp.) files, respectively.First, data were filtered to remove incomplete 24-h data from the first and last observations.Second, data were screened for each cow individually to filter unrealistic data (one bolus had a broken pH probe, resulting in biologically irrelevant data).The raw pH data of all boluses (10-or 15-min resolution) were visually checked for any potential drift or outliers.To facilitate extended screening of the data, the hourly aggregated data were decomposed by using additive timeseries decomposition to transform the data into 3 separate components: seasonal, trend, and random noise.Until d 21 in lactation (data collected after this day were systematically removed from analysis), no systematic trends indicating sensor drift were found.A limited number of cases showed timestamps with missing data.Those were imputed using linear interpolation (package imputeTS; v3.1; Moritz and Bartz-Beielstein, 2017).After data screening, data were converted to Central European Time (CET) and aggregated at the hour level.

Calculation of pH Metrics and Classification.
To describe the rumen pH profile in detail, 3 pH metrics were calculated from the individual daily pH series.These calculations were performed for each cow and day present in the data set (data from the dry period as well as during lactation, with an upper threshold of 21 DIM).First, daily time below pH 6 was calculated for each cow.The pH in the reticulorumen is generally considered higher than that in the rumen (Beauchemin et al., 2003).Neubauer et al. (2018) found a 0.2-unit difference between the pH of the free rumen liquid and the pH in the reticulum.Therefore, we used a threshold value of pH 6.0 instead of 5.8 (AlZahal et al., 2007).Second, a logistic curve was fitted, using the drc pack-age (v3.0-1;Ritz et al., 2015).The equation of the logistic curve is given in equation [1]: where y is the accumulated time (min/d, 0-1,440) below a certain pH value x, β 2 is the upper limit (i.e., 1,440 min), β 0 is the slope, and β 1 is the inflection point.The β 0 parameter was extracted and used as the variable of interest.High β 0 values represent a more stable pH pattern throughout the day (Colman et al., 2012).Therefore, β 0 gives an indication of the diurnal variation of the pH pattern.Finally, the third metric was calculated from a static model, built using the package lme4 (v1.1-26; Bates et al., 2015) to describe the pH pattern of each individual cow, adapted from Denwood et al. (2018).This model describes the pH of each cow in function of the time in lactation and some harmonic functions (i.e., cos and sin functions).The harmonic functions allow to describe the diurnal pH patterns.
Cow was defined as a random intercept.The model is given in equation [2]: where y i is the pH value of cow i, b i is the intercept of cow i, α x are coefficients, t is hours in lactation (set at 00:00:00 CET on the day of calving), and ω is given by (2π)/24.The model was refitted every 24 h to include present and past pH measurements.Deviations from this predicted model were calculated as squared errors (SqEr); that is, a metric indicating the deviation from an expected pH pattern.
To study the association between the various pH metrics and the blood parameters, measured at discrete time points, k-means clustering procedures were used to cluster cows into 2 groups based on longitudinal pH data obtained from calving up to 3 wk in lactation.For each of the 3 metrics (time rumen pH < 6, β 0 , and SqEr), this k-means clustering was performed to group cows with similar pH characteristics.In this analysis, only cows with complete data up to 3 wk in lactation were included (n = 38).The categorical result of this clustering was used as a fixed effect in mixed effects models studying the association between the clusters and blood parameters.
Before performing the k-means clustering, the natural logarithm was taken for each of the different metrics, followed by centering and scaling of the data.Second, 3 k-means clustering procedures were performed, 1 for Heirbaut et al.: RETICULORUMINAL pH METRICS DURING THE TRANSITION PERIOD each of the 3 metrics.Clustering was done using the algorithm of Hartigan and Wong (1979).The number of clusters was based on the elbow method, a visual inspection of the total within the sum of the squares as a function of the number of clusters, using the package factoextra (v1.0.7;Kassambra and Mundt, 2020).Afterward, clustering stability was assessed using a bootstrap resampling scheme with 100 iterations and the Jaccard similarities of the original clusters to the clusters based on the bootstrapped data (package fpc; v2.2-5; Hennig, 2020).The cluster solution contained 2 clusters, denoted clusters pH + and pH − .Data from clinically diseased cows (5 cows; defined as cases that required intervention by a veterinarian or farm staff) were not taken into account in this clustering procedure but were treated as a separate group.
Association Models.After clustering, the effect of clustering on the blood parameters was investigated using a linear mixed effect model built in the package lme4 (v1.1-26; Bates et al., 2015).Cow was defined as a random effect.The day of blood sampling, pH cluster, and parity were used as fixed effects and forced into the models.Biologically 2-way interactions were evaluated and omitted if P > 0.10.The clinically diseased cows were considered as a separate group.The inclusion of an extra factor (bolus type) was considered.However, despite similar parity distribution and random selection, the 2 bolus-type groups (eCow vs. SmaXtec bolus) showed some physiological differences [e.g., differences in NEFA and IGF-1 concentration; NEFA least squares means (LSM; ± SE) on d 3 and 21: 0.566 (± 0.0714) and 0.626 (± 0.0790) mmol/L vs. 0.388 (± 0.0501) and 0.415 (± 0.0537) mmol/L; IGF-1 LSM of 4 sampling days during lactation: 92.4 (± 8.25) ng/mL vs. 155 (± 8.46) ng/mL; P < 0.05].Obviously, these physiological differences are unrelated to the type of bolus inserted in the reticulorumen, and inclusion of bolus type in the models would have would have erroneously attributed (part of) these physiological differences to the bolus type.
Other parameters (milk production and composition, feed intake, and activity) that were measured daily were studied in mixed effects models using the calculated pH metrics as continuous predictors.This allowed us to study the associations in a more profound way, taking into account daily fluctuations of the pH metrics and using the full explanatory power of the calculated pH metrics.The (evaluated) factors for the mixed effect models are summarized in equation [3]: The target variable y (e.g., total DMI) at DIM i for cow j is described by an intercept μ, the fixed effects pH metric (metric), health status (health), parity, DIM, and the interaction DIM × metric; U j is the random intercept for cow j.The coefficients are given by a 1 through a 5 .
The variable health was used as binary variable to correct for clinical diseases.Because time pH <6 was characterized by a large number of zero values (57% of the daily values), this variable was converted to a categorical variable with 4 different left-open intervals (0, 0-100, 100-360, and 360-1,440 min/d pH <6).The choice of intervals and number of intervals aimed to minimize the possible loss of information due to conversion, taking into account the limitations of the data distribution.Days in milk and the pH metrics were always forced in the models.Parity and health status, as well as biologically relevant interactions, were offered to the model.Forward model selection was used, based on a likelihood-ratio test.After modeling, variance inflation factors were checked.Normality of residuals was visually assessed by Q-Q plots.If necessary, variables were transformed by taking the natural logarithm.Reported P-values were based on Satterthwaite's method.
Dynamic Linear Modeling.To obtain deeper insight into daily deviations of pH patterns, a fourth model was constructed using a dynamic linear model (DLM).A DLM has an adaptive short-term nature with more extreme forecast errors when a process deviates from its normal status (Jensen et al., 2018).Therefore, this model could be a better choice to study daily pH patterns and detect deviations from the individual timeseries compared with the above-described methods.Data were split into training and test data.All cows with reported clinical health problems were put into the test set, together with an equal number of randomly selected nondiseased cows, resulting in a total training: test split ratio of 71: 29 (training set: 25,512 observations, data from 32 cows; test set: 10,368 observations, data from 12 cows).The DLM was constructed using the base R language, based on the training set, using an observation equation (equation [4]) and a system equation (equation [5]): Equation [4] describes the mean pH (Y) at time t (t = 0 at CET 00:00:00 on the day of calving) in relation to the parameter vector θ; equation [5] describes the evolution of this parameter vector θ over time by using the previous parameter vector θ t−1 and a system matrix G t .′ F t is a transposed design matrix.The terms v t and w t are error terms and are assumed to be internally and mutually independent.V is the observational covariance matrix, and W is the systematic covariance matrix.The true underlying pH level is estimated by

F
The observed pH values differ from this expected value by the error term v t .Each time a new observation Y t is observed, the parameter vector θ t is updated using a Kalman filter (Kalman, 1960) using the formulas as described by West and Harrison (1997).
After updating using all the prior information and the last observation, the posterior distribution of this parameter vector θ t is given by N(m t , C t ), where m t and C t are the mean and variance-covariance matrices, respectively.
The DLM describing the pH during lactation is characterized by the following matrices ′ F t (equation [6]) and G t (equation [7]).

′ = [ ]
The design and system matrices consist of a linear growth model and 2 harmonic waves (approximated using Fourier sums).The period ω is defined as 2π/24.The variance terms V and W in equations [4] and [5] are assumed to be constant.The observational variance V was estimated by making a standalone linear model (equation [8]) for each cow in the training set: where Y i,t refers to the pH value of cow i at time t, α 0,i is the intercept of cow i, α 1,i to α 5,i are fixed effects, t is the time in lactation (expressed in hours), and ω is defined as 2π/24.The parameter vector m 0 and variancecovariance matrix C 0 , needed for initializing the DLM, were estimated using an overall linear model describing the pH of all cows in the training set using the same factors as in equation [8].The parameter vector θ t con-sisted of the following elements: a 0 , a 1 , a 2 , a 3 , a 4 , and a 5 .The system variance W t was unknown and was therefore estimated using a discount factor δ (0 < δ ≤ 1) (equation [9]): In equation [9], C t−1 is the variance of the parameter vector θ t−1 .The selection of the discount factor δ was based on minimizing the root mean square error (RMSE) of the forecast errors of the DLM.
For each time step, the DLM produced a forecast f t (equation [10]): In equation [10], a t represents the prior mean given in equation [11]: The forecast errors e t = (Y t − f t ) were standardized according to equation [12]: [12] In equation [12], Q t is the one-step forecast variance, which was estimated as part of the Kalman filtering, and u t is the standardized forecast error.
After constructing the DLM based on the training set, the DLM was evaluated on the test set.During this evaluation, the RMSE was calculated and the normality of the distribution of the forecast errors assessed.After evaluation, the DLM was applied to the whole data set (training + test).
Prediction of Daily Parameters.Different prediction models were built to predict daily total DMI using data from the different pH metrics and models.Depending on the features used, different classes of models were defined (overview in Supplemental Table S3; https: / / doi .org/ 10 .6084/m9 .figshare .(5) DLM model (features of the benchmark model, hourly cumulative sums of u t , number of times u t exceeds 1 up to 5, hourly u t values, daily minimum, maximum, mean, median, first and third quantiles of u t , daily differences of first and last cumulative u t value, autocorrelation values of the u t values, daily sum of u t and absolute u t and cumulative mean of the absolute u t over the different days); and (6) full model (combining all inputs).During the dry period, the number of days relative to calving was calculated based on the expected calving date.All models were random forest models.For modeling, the package ranger (v0.13.1;Wright and Ziegler, 2017) and the caret framework (v6.0-84;Kuhn et al., 2019) were used.The number of trees was fixed at 500.Models were fitted using a 10-fold cross-validation strategy (with 10 repeats).Data from a single cow were always exclusively included in either the training set or the test set to avoid any data leakage.
An overview of the clustering characteristics during lactation according to the 3 different grouping crite-ria (time pH <6, β 0 , and SqEr) is given in Table 1.The k-means clustering resulted in differences for the variables included in the clustering as well as for SqEr (β 0 clustering) and β 0 (SqEr clustering).The number of cows in each cluster and the overlap between the different clustering procedures are reported in Figure 1.The pH − group, clustered based on time pH <6 during lactation, showed greater time pH <6 during the last week before calving (140 ± 232.1 min/d) compared with the pH + group (5 ± 15.6 min/d), based on the Tukey post hoc test (P < 0.05).Clinically diseased cows did not differ prepartum compared with both groups (1 ± 4.2 min/d).The pH − group, clustered based on β 0 and SqEr, did not show differences in β 0 and SqEr, respectively, during the dry period (data not shown).
In this study, 5 cows with complete pH data suffered from clinical diseases [displaced abomasum (n = 2), hypocalcemia (n = 2), severe lameness (n = 1)].As clinical diseases may largely affect blood metabolites, clinically diseased cows were removed before performing the k-means clustering, and those cows were considered as a separate group.As such, all clinically diseased cows were combined in the group "clinically diseased." Table 2 shows the metabolic profile of the clusters.Cows in the pH − group for β 0 clustering had higher NEFA concentrations compared with those in the pH + group (during the first 9 d in lactation).When clustering was based on SqEr, no differences in NEFA concentrations were noticed after post hoc pairwise comparisons.Clusters based on time pH <6 did not differ in blood NEFA concentrations.No differences were observed for BHB concentration.Glucose, insulin, and IGF-1 did not differ between the pH − and pH + clusters, except for insulin prepartum for clustering based on time pH <6.The concentration of IGF-1 was lower for diseased cows than for cows of the pH − cluster based on SqEr.Differences for fructosamine between the pH − and pH + groups were observed only when clustering was based on SqEr, whereas cows in group pH + were characterized by higher concentrations on d 9 (Supplemental Table S5; https: / / doi .org/ 10 .6084/m9 .figshare.19336457;Heirbaut et al., 2022).
The output of the mixed models studying the associations between pH metrics and total DMI is given in Table 3. Coefficients of the parameters health status and DIM differed from zero in the 3 models.Time pH <6 was not related to total DMI, whereas SqEr was negatively associated with it (P < 0.05).A positive interaction effect was observed between β 0 and DIM (P < 0.001).Squared error was not associated with milk production (P = 0.79), whereas β 0 had a negative interaction with DIM (P < 0.01).Time pH <6 (0-100 and 100-360 min/d) had negative estimates for milk production (P < 0.001 and P < 0.05, respectively).Milk fat content was not related to any of the pH metrics.However, protein tended to relate to time pH <6 (except for the interval 0-100 min/d; P < 0.1) and SqEr (P < 0.05).A positive estimate for the interaction effect was found between β 0 and DIM (P < 0.05).The fat: protein ratio was not related to pH metrics.Only SqEr was associated with BHB concentration in the milk, with a negative main effect (P < 0.05) and a positive interaction term SqEr × DIM (P < 0.05).
Time < pH 6 was not associated with ruminating time, but a positive association was found with β 0 (P < 0.001).The negative effect of SqEr on ruminating time (P < 0.001) decreased or reversed with increasing DIM (positive estimate SqEr × DIM, P < 0.001; Table 4).Similar findings were observed for the total number of ruminating chews, except that time pH <6 (100-360 min/d) was significant in the model (P < 0.05).In contrast to the results for total DMI, no association was   found between total eating time and SqEr (P = 0.19).Also, the estimate of β 0 was not significant (P = 0.81), but time pH <6 (range 100-360 min/d) was negatively associated with eating time (P < 0.01).The number of lying bouts and motion index (both ln-transformed) were only associated with β 0 .The model of lying bouts had a negative estimate for β 0 × DIM (P < 0.01), and motion index had a positive β 0 estimate (P < 0.05).The number of steps (ln-transformed) was only associated with time pH <6 (360-1,440 min/d): a positive main effect of time pH <6 (360-1,440 min/d; P < 0.001) and a negative interaction of time pH <6 (360-1,440 min/d) with DIM (P < 0.05).All pH metrics were associated with lying time and interacted with DIM.Time pH <6 (0-100 and 360-1,440 min/d) had negative main effects (P < 0.05 and P < 0.01, respectively); however, this effect was counteracted by increasing DIM (interaction of time pH <6 and DIM; P < 0.01).β 0 and SqEr also interacted with DIM, resulting in negative (P < 0.001) and positive (P < 0.01) interaction estimates, respectively.

Dynamic Linear Model
The DLM was fitted to all observations.The standardized forecast errors followed an approximately normal distribution.The RMSE for training data was 0.12.When fitting the DLM on test data, a slightly higher RMSE of 0.14 was obtained.

Prediction Models
Given the significant effect of pH deviations (SqEr) on total DMI (Table 3) and the importance of feed intake during the transition period, several prediction models were constructed using daily feed intake and pH (model) data, with a total of 1,123 observations for the dry period and lactation.The aim of these models was to evaluate which type of pH metrics (including data from the DLM) could give the best insight into variation of feed intake.The results of the 10-fold crossvalidation with 10 repeats are given in Supplemental Figures S1 (RMSE) and S2 (R 2 ; https: / / doi .org/ 10 .6084/m9 .figshare.19336457;Heirbaut et al., 2022).The differences in RMSE were rather small, with the lowest RMSE obtained for the random forest models using the full data set (mean RMSE = 3.02 kg/d).All other models performed worse (P < 0.001), with the worst performance, based on RMSE, obtained for the benchmark model with a mean RMSE of 3.40 kg/d.All models, except for the β 0 model (P = 0.41), performed better than the benchmark model (P < 0.05).Similar to the evaluation based on RMSE, the random forest models using full data showed the highest R 2 (mean R 2 Different letters indicate differences between groups pH − and pH + within each clustering approach. A,B Different letters indicate differences from the clinically diseased group within each clustering approach. 1 The top row within each parameter contains the prepartum concentration of the blood parameter (d −7 relative to expected calving date) ± SEM; the bottom row shows the mean ± SEM of the postpartum period (d 3, 6, 9, and 21 relative to calving) of the pH clusters (pH − and pH + ) based on time pH <6, β 0 , or squared error.
2 Day of blood sampling, pH cluster, and parity were used as fixed effects and forced into the models.Biologically 2-way interactions were evaluated and omitted if P > 0.10.Day of blood sampling was always significant (P DIM < 0.01), whereas parity was never significant in the models (P parity > 0.05). 3 Five cows (out of 38) could not be sampled prepartum because they calved more than 7 d before the expected calving date.
= 0.55).The benchmark model had the lowest R 2 with a value of 0.42.All other models performed better (P < 0.05).

DISCUSSION
This study focused on the description of rumen pH through different metrics and their association with animal characteristics during the transition period.The metrics were mutually correlated, particularly absolute u t , SqEr, and β 0 parameters, which is not surprising because they all reflect the (diurnal) variation in pH profile.However, there are also some differences between those metrics: the DLM and the static model, and their metrics u t and SqEr, attempt to describe (deviations from) normal diurnal pH patterns, whereas the β 0 parameter of the logistic curve describes the overall daily variation in rumen pH without accounting for normal diurnal pH fluctuations.There was a strong correlation between the absolute u t and the SqEr parameter, which is reasonable because both rely on the same harmonic functions.Whereas the static model already had some dynamic property by means of refitting new data every 24 h, the DLM updated every hour.Correlations were observed between time pH <6 and the absolute u t , SqEr, and β 0 parameters.Pearson correlations based on aggregated data at the cow level of time pH <6 (ln-transformed) with the absolute u t , SqEr, and β 0 were 0.36, 0.44, and −0.37, respectively.Those correlations suggest an increase in pH variation when the cow is characterized by a longer period of low pH.In contrast, based on a principal components analysis including rumen pH variables and milk fatty acid proportions, Colman et al. (2012) concluded that β 0 and time pH <6 were uncorrelated because each of the parameters correlated highly with either the first or the second principal component (62.1% of the variance).However, no correlation matrix was reported in that study.Despite the correlations, clustering based on time pH <6 did not result in significant differences for β 0 and vice versa.
Subacute ruminal acidosis is generally associated with depressed milk yield (Abdela, 2016).In this study, milk yield was associated with β 0 and time pH <6.Denwood et al. (2018) found that cows with lower pH deviations had a higher milk yield.In the current study, however, the coefficient of SqEr was not significantly different from zero in the mixed effects model.This discrepancy could be related to parity effects, which were not included in the models of Denwood et al. (2018).It could also be related to lactation stage, because Denwood et al. (2018) did not focus on early lactation.
During the transition period, DMI is generally considered a critical factor determining lactation success (LeBlanc, 2010).High pH variation, in terms of SqEr, was characterized by a depressed total DMI, whereas high time at pH < 6 was not.This suggests that the diurnal pH pattern is more associated with feed intake than the time pH <6 during the transition period.This finding is in line with the results of Denwood et al. (2018).Nevertheless, several former studies indicated a reduced total DMI for SARA cows (or goats), categorized based on time when pH was below a specific threshold (Krajcarski-Hunt et al., 2002;Danscher et al., 2015;Giger-Reverdin, 2018).Differences in the experimental design could explain this discrepancy, with SARA induction protocols applied in the former studies but not in the current one.Feed intake data (combined with rumination, water intake, and animal performance characteristics) have also been used to describe temporal pH patterns (Mensching et al., 2020).Conversely, in our study, pH patterns were used to predict feed intake on a daily basis.However, the predictive ability, assessed by RMSE, of these models was only slightly superior to that of the benchmark model (the best random forest model had an R 2 of 0.55 and RMSE of 3.02 kg/d, whereas the benchmark model had a cross-validated R 2 of 0.42 and RMSE of 3.43 kg/d).A possible reason for the limited improvement in predictive ability could be related to the importance of feeding behavior (i.e., the eating events) rather than total DMI in determining pH dynamics (Dijkstra et al., 2020).Also, rumination behavior is an important aspect to consider when modeling pH dynamics (Mensching et al., 2020).Both β 0 (P < 0.001) and SqEr (P < 0.05) were associated with daily rumination time.However, time below pH 6 was not significant in the model, whereas in the study of DeVries et al. (2009), a shorter rumination time was observed for cows at risk for SARA.Nevertheless, in that study, differences were accompanied by dietary variations, with cows at risk of SARA receiving a lower forage: concentrate ratio.
Numerous studies have shown that the disease state of an animal can influence activity levels or pattern or both.Often the transition period is used as a monitoring window, because most diseases occur during this period (Hendriks et al., 2020;Lomb et al., 2020;Steele et al., 2020;Stevenson et al., 2020).In this study, time pH <6 was associated with the number of steps and lying time, but not with the motion index or number of lying bouts.In the research of DeVries et al. (2009), SARA bouts were accompanied by changes in behavior.However, no overall differences were found between cows at high risk of SARA and cows at low risk of SARA.Therefore, it could be more appropriate to study intra-animal changes in activity in relation to SARA events rather than the overall activity level of an animal considered at risk over a longer period.In this regard, Wagner et al. (2020) used machine learning to predict SARA events by detecting abnormal activity patterns, albeit still with limited prediction potential.Generally, clinically diseased cows had a low rumen pH, combined with a clearly depressed feed intake.Whether this low pH is a cause or consequence of the disease remains unclear.Although no prepartum differences in time pH <6 were observed in diseased cows, their feed intake was depressed during the 2 wk before calving compared with cows in cluster pH − or pH + .Depressed feed intake is a general sign of illness; for example, Van Winden et al. (2003) observed a lower total DMI before diagnosis of displaced abomasum.Also, in the study of Huzzey et al. (2007), cows at risk of metritis could be identified using total DMI.Lukas et al. (2008) used changes in DMI as an indicator of health status.Obviously, predicting disease events by monitoring deviation of an expected pattern (of feed intake or rumen pH) should be a future research focus, as suggested by Denwood et al. (2018).However, for this modeling purpose, a much larger data set would be required, including many more clinically diseased cows than the 5 cows in this study.
Because SqEr was negatively associated with feed intake (P < 0.05) but not associated with milk production (P = 0.79) and β 0 was positively associated with milk production (P < 0.001) and interacted with DIM for feed intake (P < 0.001), differences in circulating NEFA concentrations in the blood could be expected for clusters based on β 0 and SqEr due to differences in energy status and mobilization of body reserves.Only for clustering based on β 0 were significant differences found between the pH − and pH + groups.Post hoc testing did not reveal any differences between pH − and pH + for the SqEr clustering; however, the average NEFA concentration was above the 0.60 mmol/L threshold for the SqEr pH − group but not for the SqEr pH + group.In the literature, no studies were found describing the association between metabolic health and pH variability, because most studies focus only on the use of single threshold values to classify SARA cows.Stefańska et al. (2017) attributed the lower NEFA concentrations for acidotic cows compared with healthy cows to increased availability of propionate and glucose.Similar results were found by Guo et al. (2013).Although in our study, differences in NEFA were not observed for clusters based on time pH <6, this does not necessarily contradict the literature because the former studies included more extreme pH profiles than observed in our study (Guo et al., 2013;Stefańska et al., 2017).When performing post hoc tests, BHB did not differ when clustering using SqEr; however, there was a tendency for the fixed effect of the pH − group (P = 0.054), which could indicate insufficient availability of oxaloacetate to oxidize all acetyl-CoA produced during the β oxidation of NEFA (White, 2015).Group pH − based on SqEr was also characterized by a lower fructosamine concentration on d 9 (interaction effect DIM × group), which could indirectly reflect a fatty liver because fructosamine has been suggested as a parameter to monitor glucose status over a prolonged period (Mostafavi et al., 2015).In contrast, glucose concentration did not differ.However, glucose is considered to have important shortterm variations, whereas fructosamine is more stable (Mostafavi et al., 2015).Although the hormones insulin and IGF-1 are important indicators of metabolic health status, they did not differ among pH clusters (except insulin at d −7 when clustering using time pH <6).In line with this, it should also be noted that none of the groupings had an average BHB concentration >1.2 mmol/L serum, which is often used as a threshold for hyperketonemia (as reviewed by Benedet et al., 2019).

CONCLUSIONS
In this study, we used different metrics to describe (diurnal) pH patterns during the transition period and showed that metrics focusing on pH variation could be more important under practical dairy management systems compared with a single pH threshold value.In line with previous results during the transition period, pH patterns were relevant.In particular, DMI and rumination time, both considered important parameters during the transition period, were associated with SqEr and β 0 but not with time pH <6.Moreover, time pH <6 was also not associated with NEFA concentration.Finally, based on our findings and the relation with metabolic health, monitoring rumen pH could be beneficial in regard to disease identification.However, more data from different settings need to be collected to evaluate this further.
Ilke van Hese [Animal Sciences Unit, Flanders Research Institute for Agriculture, Fisheries and Food (ILVO), Melle, Belgium; Department of Reproduction, Obstetrics and Herd Health Faculty of Veterinary Medicine, Ghent University, Merelbeke, Belgium] and Mirjan Thys (Animal Sciences Unit, ILVO, Melle, Belgium) for assistance during sampling.The authors have not stated any conflicts of interest.
19336457; Heirbaut et al., 2022): (1) benchmark model (only using DIM, parity, and sine cosine-transformed day of the year); (2) basic pH model (features of the benchmark model combined with mean and SD rumen pH, area < pH 6, time pH <6, and the cumulative mean of time pH <6 over the different days); (3) β 0 model (features of the benchmark model combined with β 0 and cumulative mean β 0 over the different days); (4) SqEr model (features of the benchmark model, SqEr, and cumulative Heirbaut et al.: RETICULORUMINAL pH METRICS DURING THE TRANSITION PERIOD mean SqEr over the different days); a-c Means with different superscripts within a row are different (P ≤ 0.05). 1 Clustering procedures using time pH <6, β 0 (the slope of a logistic curve that describes the accumulated time below a certain pH value), or SqEr (the squared error of a harmonic model describing the diurnal pH patterns of each cow separately) as clustering parameter.2Variable for which ANOVA comparison is performed, including the clustering parameter.3Largest SEM of the groups is reported.4Parameter used to perform k-means clustering. 5Number of cows within each cluster.

Figure 1 .
Figure 1.Venn diagram of the k-means clustering solutions based on time pH <6, β 0 (the slope of a logistic curve that describes the accumulated time below a certain pH value), or SqEr (the squared error of a harmonic model describing the diurnal pH patterns of each cow separately).The number of cows in each segment is reported as the number of cows in cluster pH − followed by the number of cows in cluster pH + in parentheses.
Heirbaut et al.: RETICULORUMINAL pH METRICS DURING THE TRANSITION PERIOD

Table 1 .
Heirbaut et al.: RETICULORUMINAL pH METRICS DURING THE TRANSITION PERIOD Mean pH characteristics during lactation (d 1 to 21) of the 2 clusters (pH − and pH + ) based on 3 clustering procedures 1

Table 2 .
Heirbaut et al.:RETICULORUMINAL pH METRICS DURING THE TRANSITION PERIOD Summary of metabolic parameters of successive samplings1,2