Modeling the effects of heat stress in animal performance and enteric methane emissions in lactating dairy cows

Heat stress (HS) negatively affects dry matter intake (DMI), milk yield (MY), feed efficiency (FE), and free water intake (FWI) in dairy cows, with detrimental consequences to animal welfare, health, and profitability of dairy farms. Absolute enteric methane (CH 4 ) emission, yield (CH 4 /DMI), and intensity (CH 4 /MY) may also be affected. Therefore, the goal of this study was to model the changes in dairy cow productivity, water intake, and absolute CH 4 emissions, yield, and intensity with the progression (days of exposure) of a cyclical HS period in lactating dairy cows. Heat stress was induced by increasing the average temperature by 15°C (from 19°C in the thermoneutral period to 34°C) while keeping relative humidity constant at 20% (tem-perature-humidity index peaks of approximately 83) in climate-controlled chambers for up to 20 d. A database composed of individual records (n = 1,675) of DMI and MY from 82 heat-stressed lactating dairy cows housed in environmental chambers from 6 studies was used. Free water intake was also estimated based on DMI, dry matter, crude protein, sodium, and potassium content of the diets, and ambient temperature. Absolute CH 4 emissions was estimated based on DMI, fatty acids, and dietary digestible neutral detergent fiber content of the diets. Generalized additive mixed-effects models were used to describe the relationships of DMI, MY, FE, and absolute CH 4 emissions, yield, and intensity with HS. Dry matter intake and absolute CH 4 emissions and yield reduced with the progression of HS up to 9 d, when it started to increase again up to 20 d. Milk yield and FE reduced with the progression of HS up to 20 d. Free water intake (kg/d) decreased during the exposure to HS mainly because of a reduction in DMI; however, when expressed in kg/kg of DMI it increased modestly. Methane intensity also reduced initially up to d 5 during HS exposure but then started to increase again following the DMI and MY pattern up to d 20. However, the reductions in CH 4 emissions (absolute, yield, and intensity) occurred at the expense of decreases in DMI, MY, and FE, which are not desirable. This study provides quantitative predictions of the changes in animal performance (DMI, MY, FE, FWI) and CH 4 emissions (absolute, yield, and intensity) with the progression of HS in lactating dairy cows. The models developed in this study could be used as a tool to help dairy nutri-tionists to decide when and how to adopt strategies to mitigate the negative effects of HS on animal health and performance and related environmental costs. Thus, more precise and accurate on-farm management decisions could be taken with the use of these models. However, application of the developed models outside of the ranges of temperature-humidity index and period of HS exposure included in this study is not recommended. Also, validation of predictive capacity of the models to predict CH 4 emissions and FWI using data from in vivo studies where these variables are measured in heat-stressed lactating dairy cows is required before these models can be used.


INTRODUCTION
Global climate change is affecting temperature, precipitation and water availability, which directly affects agriculture and livestock productivity Gisbert-Queral et al., 2021). In the United States, it is estimated that $1.2 billion is lost annually due to heat stress (HS) in dairy production systems (Key et al., 2014). Heat stress negatively affects DMI, milk yield (MY), and feed efficiency (FE; MY/DMI) in dairy cows (Becker et al., 2020). It can also have variable effects on free water intake (FWI) depending on the production level, and has detrimental consequences to animal welfare, health, and profitability of dairy farms, especially high-producing dairy cows (Polsky and von Keyserlingk, 2017;Collier et al., 2019). Climate change is likely to exacerbate these drawbacks with the projected higher frequency and intensity of heat wave events worldwide (IPCC, 2022).
Enteric methane (CH 4 ) emissions, including yield (CH 4 /DMI), and intensity (CH 4 /MY) in ruminants is affected by changes in DMI, NDF digestibility (NAS-EM, 2021), mean digesta retention time (Huhtanen et al., 2016), and animal productivity (Niu et al., 2018). Therefore, the expected decline in DMI and animal productivity and changes in water intake caused by HS may also affect environmental costs of production in cattle. With the current interest in reducing CH 4 emissions in ruminant production systems to limit global warming (Arndt et al., 2022), models that predict CH 4 emissions have become an important tool to evaluate mitigation strategies (Niu et al., 2018), when other technologies are not available to measure individual enteric CH 4 emissions. However, to our knowledge there are no models available to predict CH 4 emissions in cattle that account for the effect of HS. In addition, despite of the importance of CH 4 emissions to the entire dairy industry, there are no studies evaluating the effect of HS on CH 4 emissions in lactating dairy cows. As such, mathematical models might be a tool to forecast CH 4 emissions, whereas in vivo studies are not available. With the current trend of increased global average temperatures, higher daily maximums, and more frequent heat waves (Key et al., 2014), measuring accurately how much CH 4 is being emitted during HS exposure is essential.
Although the ability to understand and mitigate HS in dairy production systems has improved over the decades (Kadzere et al., 2002;Amaral et al., 2009;Baumgard and Rhoads, 2012;Broucek et al., 2020), there is paucity of predictive models that can accurately represent the changes in production responses with the progression (days of exposure) of HS. Such models could support more precise on-farm decisions based on the expected changes in production responses of cows exposed to HS bringing not only animal health and welfare, but also economic and environmental benefits. With the predictions from these additional models, adjustments in the diet composition or in the cooling system can be made to mitigate the negative effects of HS on animal performance and environmental costs of production with the progres-sion of HS. Recently, Benni et al. (2020) developed a model to characterize susceptibility of cows to HS based on the relationship between milk production and temperature-humidity index (THI) to define targeted treatments to individual cows according to their characteristics. However, their model does not predict DMI, FE, FWI, absolute CH 4 emissions, yield, and intensity. Therefore, the goal of this study was to model the changes in animal productivity, FWI and CH 4 emission, yield and intensity with the progression of a cyclical HS period [days of exposure to an average temperature increased by 15°C (from 19°C in the thermoneutral to 34°C), while keeping relative humidity constant at 20% (THI peaks of approximately 83) in climate-controlled chambers for up to 20 d] in lactating dairy cows.

MATERIALS AND METHODS
No Animal Care Committee approval was necessary for the purposes of this study as all the information required was obtained from pre-existing data sets.

Database
The database used in this study was composed of one data set of 1,675 individual records (raw data) of DMI and MY from 82 multiparous lactating dairy cows from 6 studies conducted in The University of Arizona Shwartz et al., 2009;Skrzypek et al., 2010;Wheelock et al., 2010;Baumgard et al., 2011). Cows were housed in 12-h light and dark cycles in climate-controlled chambers at 2 experimental periods with exception of cows from the study of Baumgard et al. (2011), which had 14-h light and 10-h dark cycles. All cows were fed ad libitum. In the first period, DMI and MY were recorded in thermoneutral conditions with an average temperature of 19°C (varying between 18 to 20°C) and relative humidity constant at 20% (average THI = 63.3; varying between 62 and 64) for up to 9 d. In the second period, cows were exposed to a cyclical HS (average temperature increased by 15°C; from 19°C in the thermoneutral period to 34°C) with daily temperature varying between 29.4 to 38.9°C and relative humidity constant at 20% (THI varying between 73 and 83) for up to 20 d. During HS period the temperature remained around 29.4°C between 0000 and 0700 h (THI = 73); thereafter, the conditions became increasingly warmer until peaking at a temperature of 38.9°C between 1300 and 1500 h (THI = 83). After the peak, temperatures gradually declined until 29.4°C at 2300 h. A description of studies is available in Table 1. The data are presented in Supplemental Figures S1 to S10 (https: / / data .mendeley .com/ datasets/ 8c93y2k2y4/ 1; Souza, 2023) and summary statistics of animal performance, CH 4 emissions (absolute, yield, and intensity), FWI, animal and dietary characteristics are given in Table 2. Free water intake and absolute CH 4 emissions were not measured in the studies included in the database and therefore were estimated.   Baumgard et al. (2011), and Rhoads et al. (2010). Data from pair-fed cows are not included. Absolute enteric methane emission and free water intake were predicted using the equations described by NASEM (2021).
Free water intake was predicted using the equation recommended for lactating dairy cows in the NASEM (2021) as follows: FWI (kg/d) = −91.1 + (2.93 × DMI) + (0.61 × DM) + (0.062 × NaK) + (2.49 × CP) + (0.76 × TMP), where DM is the dry matter content (%), NaK is the sodium plus potassium concentration (mEq/kg DM), CP is the crude protein content (% DM), TMP is daily mean ambient temperature (°C). Sodium and K were estimated based on diet composition reported in the experiments and the feed composition tabular values reported in the NASEM (2021).
Absolute CH 4 emissions (g/d) was predicted using the equation described in NASEM (2021) for lactating dairy cows as follows: where FA is the dietary fatty acid content (% DM) and dNDF is dietary digestible neutral detergent fiber content (% DM) obtained after a 48-h in vitro incubation. Dietary fatty acids and dNDF were estimated using the diet composition reported in the experiments and the feed composition tabular values reported in the NASEM (2021). In the equation described above the unit Mcal/d was converted to g/d considering that 1 g of CH 4 releases 0.01425 Mcal of energy when combusted as described in the supplemental materials (https: / / data .mendeley .com/ datasets/ 8c93y2k2y4/ 1; Souza, 2023).

Modeling the Effect of HS Length in Lactating Dairy Cows
A semiparametric statistical modeling approach was used to represent the relationship between each of the following animal responses: DMI, MY, FE, absolute CH 4 emissions, yield, intensity, or FWI (kg/d and kg/ kg DMI) and HS length. Specifically, generalized additive mixed models (Lin and Zhang, 1999) were used with errors assumed to be normally distributed and with an identity link function. Models were fitted to the data with the mixed GAM computation vehicle with automatic smoothness estimation from mgcv package (Wood, 2023) in R version 4.1.1 (R Core Team, 2021).
Smooth functions were constructed with thin plate splines as described in Wood (2006) and variance com-ponents and random effects estimated with maximum likelihood (Pinheiro et al., 2021).
Two models were developed to describe the relationship of the response variables with the progression (days of exposure) of HS. The model for the response variables DMI, FE, absolute CH 4 emissions, yield and intensity, or FWI was specified as follows: where y ijk is the kth record of the response variable during HS on the jth cow (j = 1, …, 81) in the ith study (i = 1, …, 6), β 0 is the intercept and β 1 is the parameter describing the relationship between the response variable at the HS period and its average at the thermoneutral period (baseX ijk ) for cow j in study i, s is a centered twice-differentiable smooth function of the variable day ijk describing the time from the onset of the HS, b i is the random effect associated with the ith study, b ij is the random effect associated with the jth cow within the ith study, and ε ijk is the error. Study random effects b i were assumed to be independent for different i. Random effects b ij were assumed to be independent for different i or j and to be independent of study random effects. The within group errors ε ijk were assumed to be independent for different i or j and to be independent of the random effects. Further, it was as-

( )
Similarly, the model for MY was specified as follows: where y ijk is the kth record on MY during HS on the jth cow (j = 1, …, 82) in the ith study (i = 1, …, 6), β 0 is the intercept and β 1 is the parameter describing the relationship between milk yield at the HS period and the average milk yield at the thermoneutral period (baseMY ijk ) for cow j in study i, β 2 and β 3 are parameters describing the relationship between MY and DMI. It is important to note that we allowed a potential quadratic effect to represent the decline in the rate of MY secretion at high DMI intakes, s is a centered twicedifferentiable smooth function of the variable day ijk describing the time from the onset of the stress. Further, b i is the random effect associated with the ith study, b ij is the random effect associated with the jth cow within the ith study and ε ijk the error. Mutual independence assumptions of errors and random effects were the same

( )
A leave-one-out cross-evaluation was used to evaluate model performance using the function training of the R package caret, where 75% of the data used for training and 25% for model evaluation with resampling iterations (k) = 10 (James et al., 2013). Root mean square error was used as model metrics to assess model performance.

Determining Changes in Response Variables During HS
The contributions of the additive smooth terms to the expectations were used to determine the changes in the response variables. Specifically, the fitted smooth curves were used to determine changes in DMI, MY, FE, absolute CH 4 emissions, yield and intensity, and FWI with the progression of HS. The instantaneous or marginal change (k) in a response variable with respect to an explanatory variable x (day after stress onset) was determined with the first derivative of the smooth curve: . [5] The average change between points a and b, k a b , ( ) in a response variable, was therefore computed by . [6]

Modeling the Effect of HS Length in Lactating Dairy Cows
In this modeling exercise we evaluated the effect of length of HS [days of exposure to an average temperature increased by 15°C (from 19°C in the thermoneutral period to 34°C) while keeping relative humidity constant at 20% (THI peaks of approximately 83) in climatecontrolled chambers for up to 20 d] on DMI, MY, FE, FWI, absolute CH 4 emissions, yield, and intensity in lactating dairy cows under controlled conditions. The THI is one of the major indicators used to evaluate HS in dairy operations (Polsky and von Keyserlingk, 2017). However, due to the lack of variability in the average THI across days and between the studies included in the data set the option was to not include this variable in the models, which would be a challenge from an application standpoint due to variability in THI between days, seasons, and geographic regions. Thus, application of the developed models outside of the ranges of THI and period of HS exposure included in this study is not recommended. In addition, caution should be taken when using the models to predict CH 4 emissions (absolute, yield, and intensity) and FWI because the data used to develop these models were based on estimated values.

Predicting Changes in DMI, MY, and FE During HS
The reduction in feed intake is a common response of heat-stressed dairy cows (Kadzere et al., 2002). The data of DMI is presented in Supplemental Figure S1. In this study, DMI during HS was, on average, 76.5% of the mean DMI from the thermoneutral period ( Table  2). The model for DMI fitted the data well (Supplemental Figure S2) and as expected, DMI decreased with HS progression (Figure 1a). A sharp, and roughly linear, decrease in DMI was observed in the first 5 d of exposure to HS. The greatest effect of HS on DMI occurred from 5 to 9 d after stress onset. From d 10, a slow increase in DMI was observed but it did not return to the levels observed before HS exposure ( Figure 1a). The data of MY is presented in Supplemental Figure  S1. Milk yield during HS was, on average, 78.1% of the mean MY from the thermoneutral period (Table 2), which was expected because of the reduction in DMI. Similar results were reported by Benni et al. (2020) who demonstrated that cows more susceptible to HS only produce 76% of their potential when exposed to HS. The model for MY predictions during HS exposure also fitted the data well (Supplemental Figure S2). A sharp and roughly linear decrease in MY was observed in the first week of the HS progression ( Figure 1b). Milk yield was maintained low during the second week and continued to decrease in the third week but at a slower rate (Figure 1b). The initial sharp decrease in DMI lasted about 5 d, whereas for MY, the initial sharp decrease lasted about 7 d (Figure 1a and b). Therefore, even with a small recovery in the feed intake after d 5, there was a 2-d lag for the establishment of a lesser reduction rate of MY. The 2-d lag could be explained by the time it takes for the nutrients be consumed, digested, metabolized, and used for MY (West, 2003).
With the additive mixed-effects model, we determined DMI in heat-stressed dairy cows using DMI measured in the thermoneutral period and the duration of the cyclical HS. Parameters estimated in the DMI model are given in Table 3. Predicted values for the smooth function determining changes in DMI with the stress progression are given in Table 4. Predictions of DMI during HS for the population (i.e., setting random effects to zero) may, therefore, be determined with DMI baseDMI day = + × + ( ) 6.37 0.46ˆ. s For example, the predicted DMI in the third day of HS of a cow that was consuming 18 kg of DM before the stress onset is DMI = 6.37 + 0.46 × 18 + 1.52 = 16.2 kg DM. Similarly, parameters estimated for the MY model are given in Table 3. Predicted values for the smooth function determining changes in MY with the stress progression are given in Table 4. Predictions of MY during HS for the population (i.e., setting random effects to zero) were determined with the following:  intake starts to return to the pre-exposure (or not differ from that of cows not exposed to HS), whereas milk is "permanently" compromised or it does not recover to the same rate. It has been suggested that in addition to the reductions in DMI during HS, liver dysfunction , hypoglycemia, and increased glucose disposal rates (Wheelock et al., 2010)   1 DMI, milk yield (MY), feed efficiency (FE), enteric methane (CH 4 ), methane yield (CH 4 /DMI), methane intensity (CH 4 /MY), free water intake (FWI), β 0 is the intercept, β 1 is the parameter describing the relationship between the response variable at the heat stress period and its average at the thermoneutral period, β 2 and β 3 are parameters describing the relationship between MY and DMI, σ c is the SD of the random cow effect, σ s is the SD of the random study effect, and σ ε is the SD of the error. 2 RMSE is root mean square error obtained in the leave-one-out cross-evaluation. reinforce the importance of providing heat abatement opportunities to lactating dairy cows. Indeed, Gisbert-Queral et al. (2021) showed much smaller effects of HS when mediating effects of farmer management (e.g., housing, fans, shading, and sprinklers) were incorporated in models predicting heat-induced MY losses across US states. There are several strategies to reduce the negative the effects of HS on animal performance such as the dietary inclusion of yeast (Schingoethe et al., 2004), fatty acids (Wang et al., 2010), the reduction of dietary fiber concentration in high-fiber diets (Kanjanapruthipong et al., 2010) or increasing diet NDF digestibility by replacing feed ingredients (Halachmi et al., 2004), AA (Pate et al., 2020;Dou et al., 2021), chromium (Shan et al., 2020), choline (Holdorf and White, 2021), and the provision of evaporative cooling systems (Correa-Calderon et al., 2004). Thus, the predictions made using smooth values in Table 4, may help nutritionists to decide when and how to apply the strategies available to attenuate the reductions in DMI and MY with the progression of HS, which might have great potential to prevent the side effects of HS exposure. However, the decision of when and how to intervene is variable and will depend of the degree of HS and expected economic return with each strategy available based on environmental and nutritional management conditions, which is outside of the scope of this study.
Feed efficiency has been used as a metric to evaluate the efficacy in which dairy cows convert the diet into milk, and therefore it has economic implications because as more energy from feed is converted into milk, the more profitable are the cows (Britt et al., 2003). It also has environmental implications because more efficient cows will reduce nutrient excretion in manure, and lower CH 4 and carbon dioxide emissions per kilogram of milk produced (Arndt et al., 2015). During HS conditions, in addition to the reduction in feed intake, cows reduce nutrient uptake for production purposes and direct it toward thermoregulation Renaudeau et al., 2012), which can explain the reductions in MY and FE during HS. Therefore, models that predict changes in FE taking into account the HS effects can be used to evaluate animal productivity and to develop management strategies to mitigate these side effects. The data of FE is presented in Supplemental Figure S3. Parameters estimated in the FE model are given in Table 3. Predicted values for the smooth function determining changes in FE with the stress progression are given in Table 4. Predictions of FE during HS for the population (i.e., setting random effects to zero) may, therefore, be determined with FE b aseFE s day . = 0.30 + 0.81 × + ( ) The model fit-ted the data reasonably well (Supplemental Figure S4). For example, the predicted FE in the 15th day after the onset of the HS of a cow that had a FE of 1.54 before it started is FE = 0.30 + 0.81 × 1.54 + (−0.11) = 1.44. For the FE model, we determined 3 average rates of changes ( Figure 1c Therefore, from d 1 to 7, there was a reduction of 0.01 points/d in FE, from d 7 to 11, an increase of 0.01 points/d, and a reduction of 0.03 points/d from d 11 to 20 (Figure 1c). These results are in line with previous studies where FE was reduced under more prolonged HS conditions compared with the studies used in this modeling exercise (Britt et al., 2003;Su et al., 2013). The FE predictions highlight the importance of the use of heat abatement and nutritional strategies to reduce the effects of HS in dairy cows. Longer studies are needed to model the long-term effects of HS on FE in dairy cows.

Predicting Changes in FWI During HS
The data of FWI in kg/d or kg/kg DMI are presented in Supplemental Figure S5. Free water intake is positively related to DMI and might increase or decrease during HS exposure, depending of the level of production (Collier et al., 2019). Parameters estimated in the FWI model are given in Table 3. Predicted values for the smooth function determining changes in absolute CH 4 intensity with the stress progression are given in Table 4. Predictions of FWI (kg/d) during HS for the population (i.e., setting random effects to zero) may, therefore, be determined with FWI = 50.2 + 0.48 × baseFWI + ŝ(day) with ŝ day ( ) given in the eighth column of Table 4. The model fitted the data well (Supplemental Figure S6). For example, the predicted FWI in the 10th day after the onset of the HS of a cow that had an FWI of 107.0 kg before it started is FWI = 50.2 + 0.48 × 107.0 + (−4.44) = 97.1 kg. Based on Figure 1d, we determined 2 intervals of changes in FWI (kg/d) as follows: From d 1 to 10, there was a reduction of 1.97 kg/d in FWI, and an increase of 0.59 kg/d from d 10 to 20. The reductions in FWI may be interpreted as a consequence of the decrease in DMI during HS. In addition, Collier et al. (2019) demonstrated that FWI behaves differently based on the level of production, where highproducing cows (MY >25 kg/d) decrease FWI and low-producing cows (MY ≤25 kg/d) increase FWI. On average, 50% of our MY records were from cows producing more than 25 kg/d of MY. However, a better way to understand variations in FWI with HS is to express it as a proportion of DMI. Parameters estimated in the FWI/DMI model are given in Table 3. Predicted values for the smooth function determining changes in FWI/DMI with the stress progression are given in Table 4. Predictions of FWI (kg/kg DMI) during HS for the population (i.e., setting random effects to zero) were determined with FWI = 1.29 + 1.03 × baseFWI_DMI + ŝ(day) with ŝ day ( ) given in the last column of Table 4. The model fitted the data reasonably well (Supplemental Figure S6). For example, the predicted FWI on the 10th day after the onset of the HS of a cow that had a FWI of 4.84 kg/kg DMI before it started is FWI = 1.29 + 1.03 × 4.84 + 0.29 = 6.57 kg/kg DM. Based on Figure 1e, we determined 4 intervals of changes in FWI (kg/kg DMI) as follows:  (Kadzere et al., 2002;NASEM, 2021). However, the increase in FWI was rather small. In a recent study, where FWI was measured in heat-stressed dairy cows, it was demonstrated that although cows increased the number of visits and the time spent at the drinker, FWI was not greatly affected by the increase in THI (McDonald et al., 2020), which corroborates with the predictions of the FWI model. The authors hypothesized that cows may increase the number of visits and the time spent on the drinker because of the cooling effect that water causes on their skin and in the area around the drinker, where passing air over water leads to evaporative cooling.

Predicting Changes in Absolute Methane Emission, Yield, and Intensity During HS
The data of absolute CH 4 emissions is presented in Supplemental Figure S7. Rumen methane production is primarily affected by the amount of feed consumed and digested (Niu et al., 2018;NASEM, 2021). Parameters estimated in the CH 4 model are given in Table 3. Predicted values for the smooth function determining changes in absolute CH 4 emissions with the stress progression are given in Table 4. Absolute CH 4 emissions (g/d) decreased with the progression of HS up to d 9 (Figure 1f). Predictions of absolute CH 4 emissions during HS for the population (i.e., setting random effects to zero) were determined with CH 4 = 91.5 + 0.48 × baseCH 4 + ŝ(day). The model fitted the data well (Supplemental Figure S8). For example, the predicted absolute CH 4 emitted in the third day of HS of a cow that was emitting 397 g before the stress onset is CH 4 = 91.5 + 0.48 × 397 + 31.4 = 313 g CH 4 . Based on Figure 1f, we determined 3 intervals of changes in absolute CH 4 emissions as follows: The reductions and increase in absolute CH 4 emissions followed the decline and increase in feed intake during HS exposure. The data of CH 4 yield (CH 4 /DMI) and intensity (CH 4 /MY) are presented in Supplemental Figure S9. Parameters estimated in the CH 4 yield model are given in Table 3. Predicted values for the smooth function determining changes in absolute CH 4 yield with the stress progression are given in Table 4. Methane yield and intensity decreased with the progression of HS and increased again after 10 and 5 d of HS exposure, respectively (Figure 1g and 1h) Table 4, respectively. Both models fitted the data reasonably well (Supplemental Figure S10). For example, the predicted CH 4 yield in the 10th day of HS of a cow that was emitting 17.8 g/kg of DM before the stress onset is Therefore, from d 1 to 10, there was a 0.14 g CH 4 /kg DMI reduction in CH 4 yield, and from d 10 to 20 there was an increase (compared with d 1-10 period) of 0.06 g CH 4 /kg DMI. The relationships between digestibility, passage rate, and DMI are complex and will control the amount of enteric CH 4 produced (Knapp et al., 2014). It has been shown that prolonged exposure to HS may reduce rumen cellulolytic activity (Yousri et al., 1977;Bernabucci et al., 1999;Kim et al., 2022), rumen motility, and fermentation (Schneider et al., 1988) to avoid heat generation from feed digestion, which was associated with the reduction in DMI may explain the reduction in CH 4 yield. In agreement with the model predictions of this study, Yadav et al. (2016) reported reductions in CH 4 yield in nonlactating crossbred cattle when average temperature was increased from 25 to 35°C in a climatic chamber. However, Bernabucci et al. (1999) showed a reduction in passage rate and an increase in nutrient total-tract digestibility in Friesian heifers exposed continuously for 40 d to an elevated THI of 84 (33°C and 60% relative humidity), which in theory could increase the amount of CH 4 generated per kilogram of DMI. Thus, it is also possible that the CH 4 yield predictions from this study are biased because DMI is being used to predict CH 4 , which may push down CH 4 yield when DMI is decreased. In addition, CH 4 yield and intensity are ratios which pose some challenges from a modeling stand point (Souza and White, 2021). Parameters estimated in the CH 4 intensity model are given in Table 3. Predicted values for the smooth function determining changes in absolute CH 4 intensity with the stress progression are given in Table 4. The predicted CH 4 intensity in the fifth day of HS of a cow that was emitting 12.1 g/kg of MY before the stress onset is Therefore, from d 1 to 5, there was a 0.15 g CH 4 /kg MY reduction in CH 4 intensity, from d 5 to 8 there was an increase of 0.09 g CH 4 /kg MY, from d 8 to 11 there was another increase of 0.02 g CH 4 /kg MY, and from d 11 to 20 there was a third increase of 0.20 g CH 4 /kg MY. The increase in CH 4 intensity from d 11 to 20 can be explained by the simultaneous increase in CH 4 emissions (g/d) controlled mainly by the increase in DMI, and the decrease in MY during this period.

Implications and Limitations
This study provides quantitative predictions of the changes in animal performance (DMI, MY, FE, FWI) and CH 4 emissions (absolute, yield, and intensity) with Souza et al.: MATHEMATICAL MODELS FOR PREVENTING HEAT STRESS EFFECTS the progression of HS in lactating dairy cows. The models developed in this study could be used as a tool to help dairy nutritionists to decide when and how to adopt strategies to mitigate the negative effects of HS on animal performance and related environmental costs. Thus, more precise and accurate on-farm management decisions may be taken with the use of the models developed in this study. However, it is important to highlight that it is outside of the scope of this study to provide recommendations of when and how intervene in the management to mitigate HS in dairy cows, because none of the existing strategies were evaluated.
Despite efforts to derive comprehensive models explaining the changes in animal performance, CH 4 emissions, and FWI with the progression of HS, there are several limitations in this study that should be considered before these models are used out of context. Although the results of leave-one-out cross-evaluation are presented, cross-evaluation is an imperfect means of model evaluation for predictive purposes (Souza and White, 2021), and should not be considered as ratification of any models derived in this work as acceptable. Therefore, users should have caution when applying these models on farm. In addition, although it was not evaluated in this study, the data used for model derivation was obtained in climate-controlled chambers, which may differ from data collected under farm conditions. Also, the predictions of FWI and absolute CH 4 emissions, and the respective calculations of yield and intensity were based on estimates from the NASEM (2021), not observed values, which may introduce some biases in the model predictions. Thus, validation of predictive capacity of the models to predict FWI and CH 4 emissions using data from in vivo studies where these variables are measured in heat-stressed lactating dairy cows is required before these models can be used. However, the models developed in this study advance the understanding on the quantitative responses of dairy cows to the exposure to the HS, which is becoming more important with the current trend in increase in global average temperatures in the next decades. Future work should include: comparisons of extant models to the models developed in this study using external data sets; rederivation of current models using data from heatstressed dairy cows submitted to different levels of THI (broader range of temperature and humidity, including the succession of HS to simulate more closely real conditions) with or without access to evaporative cooling systems; inclusion of measured instead of estimated values of FWI and CH 4 ; longer period of exposure to HS to evaluate the long-term effect of HS; a dynamic mechanistic model integrating the HS consequences on performance, digestive and metabolic processes.

ACKNOWLEDGMENTS
Funding for this work was provided by the National Institute of Food and Agriculture of USDA (Washington, DC) through an AFRI grant (award 2021-68014-34141). The authors have not stated any conflicts of interest.