Chapter 4: Results

We estimate the single-breed model proposed in Equation 3.2.2. The applied estimation engine is gamm4j which we introduced in Section 3.4.4.2. For each breed, data is subsampled from our dataset presented in Table 3.4. The subsamples are generated by picking all samples from randomly selected farms until a sample threshold \( t_s \) is reached. \( t_s \) is determined by our available computational capacities. This resource limitation is discussed in Section 3.4.4.4. We deliberately select all available samples from each randomly chosen farm in order for the model to implicitly incorporate the structural characteristics of the dairy farm and cow herd. Initially, we conduct the analysis on subsamples spanning the entire study period which is our primary focus. Then, we divide the data into two distinct periods, prior to and subsequent to 2010 and briefly describe the results in Section 4.4. This data split is justified by previous work from Citation: Gisbert-Queral et al., 2021 Gisbert-Queral, M., Henningsen, A., Markussen, B., Niles, M., Kebreab, E., Rigden, A. & Mueller, N. (2021). Climate impacts and adaptation in US dairy systems 1981–2018. , 2(11). 894–901. https://doi.org/10.1038/s43016-021-00372-z  and also motivated by potential changes in the cows’ thermoregulatory mechanisms due to breeding, or the abolition of milk quotas (Table 2.3).

How to read the figures?

The following sections and subsections will systematically highlight the effect of THI on milk yield and ECM yield across breeds with one class of figures. On these figures, the abscissa represents the 3-day mean THI, while the ordinate indicates the milk performance metric, either milk yield or ECM yield. Each breed is denoted by a unique color and symbol, as outlined in Table 3.4. The main objective of this study is to identify the THI turning point for each breed. Consequently, each curve’s turning point, also referred as THI threshold or inflection point, is marked by a vertical line, accompanied by the corresponding THI value. In certain instances, a numerically derivable turning point is not achievable, as the slope of the curve does not transition from positive to negative derivatives, or vice versa. One category of figures compares breeds, including averages of the corresponding breed and parity, while another exclusively examines the marginal effects. The reader is strongly advised to refer to the THI table in Figure A.1 to facilitate the comprehension of this chapter.

How to read the tables?

In addition to the figures, we report tables with model specifics and principal results. The initial column in each table furnishes a reference to the model’s section within the appendix.

The first category of tables is devoted to technical model details, encapsulating the number of samples \( N \), random effect factor levels for farms \( \phi \), ZIP-month interactions \( \iota \), and animals \( \alpha \), as well as the number of basis functions \( K \) pertinent to the smooth terms associated with \( THI \) and days in milk \( DIM \), along with the actual estimated degrees of freedom \( EDF \). The table is complemented with the time required for the model estimation with gammJ. Entries necessitating a distinction between primiparous \( P \) and multiparous cows \( M \) are reported accordingly.

The second category of tables summarizes the content depicted in the figures, such as the THI turning points, if available, alongside the corresponding loss rates. The loss rates delineate the linear gradient between the THI turning point and the THI value belonging to the furthest rightmost datapoint available for the particular slope. Hence, they are a linear approximation of the non-linear marginal effect. If no turning point is available, no loss rate is provided. Given that the lactation curves \( f(DIM) \) are also smooth terms, the determination of lactation peaks utilizes the same numerical methodology as that employed for \( f(THI) \). Hence, we report the turning points for the lactation curves, albeit without further discussion. Furthermore, for each model, the marginal and conditional \( R^2 \) are reported according to the formulas in Citation: Nakagawa & , 2013 Nakagawa, S. & Schielzeth, H. (2013). A general and simple method for obtaining R2 from generalized linear mixed-effects models. , 4(2). 133–142. https://doi.org/10.1111/j.2041-210x.2012.00261.x [Equation 29 and Equation 30].

Model Details

The Appendix B provides comprehensive details for each estimated model, encompassing the model summary, diagnostic plots, and the marginal effects \( f(THI) \) over the entire available THI spectrum with the corresponding model data, as well as the lactation curves \( f(DIM) \).

Roadmap

First, an overview of the estimated models is presented in Section 4.1. Subsequently, Section 4.2 examines the impacts of THI on volumetric milk yield, while Section 4.3 addresses similar effects on the component-corrected ECM yield. The models in both sections consider the full period of study from 1982-2023. Then, a succinct description of the split-period results is provided in the subsequent Section 4.4. Finally, the results are synthesized and discussed in Section 4.5.

4.1 Model Overview

Milk Yield

Table 1.1 shows all single-breed models with the milk yield as dependent variable. \( t_s \) is set to 1’000’000 for all experiments. Only for the Jersey breed there are fewer samples available. The estimated degrees of freedom confirm that \( k \) is well-chosen for both, \( f(THI) \) and \( f(DIM) \).

ECM Yield

4.2 Milk Yield

Figure 4.1 summarizes the effect of the 3-day mean THI on the volumetric milk yield for both parities and all breeds at average 2023 levels. The subplots are in descending order with respect to the average yields with the high-yield Holstein yield on the left and the modest-yield Jersey breed on the right. Table 4.3 provides a summary of the THI turning points and the associated loss rates. For multiparous cows the peak THI points range from 53.6 for Brown Swiss to 58.12 for Jersey. The THI difference between the lowest and highest peak point is 4.56 THI points. For primiparous cows the the lowest THI peak point is 48.04 for Brown Swiss and 56.2 for Jersey. The THI difference between the lowest and highest peak point for primiparous cows is 8.16 points. In both cases, for Brown Swiss and Jersey, the curves for primi- and multiparous cows, the marginal losses are less accentuated than for the other breeds. For primiparous Simmental we observe a steep drop after a relatively stable period beyond the THI turning point. Figure 4.2 and Figure 4.3 only depict the marginal effects of THI without adding the average levels of year and parity. We observe a two-stage decrease for multiparous Holsteins: similarly as for primiparous Simmental, the milk yield decline starts at a THI of 54.71 at a moderate rate. Then, at an approximate value of 65, the loss rate is more steep.

4.3 ECM Yield

Figure 4.4 summarizes the effect of the 3-day mean THI on the component-corrected milk yield for both parities and all breeds at average 2023 levels. Hence, the milk yield is normalized with respect to fat and protein. First, by the definition of ECM, the average yield levels increase for all breeds¹. The jump for Jersey cows is expected since they are known as a high-fat and high-protein breed. Moreover, we are able to identify turning points for all breeds in the multiparous case. The lowest THI peak is at 51.05 for Holstein and the highest at 52.76 for Simmental. The difference between the lowest and highest THI peak point is 1.71. For primiparous cows, we cannot identify the THI turning points for Swiss Fleckvieh, Brown Swiss and Jersey. This is because their corresponding slopes are continuously decreasing. Therefore, from a numerical perspective, there are no turning points in these cases. However, in all three cases, fairly steep changes in the loss rates at THI values in the mid-fifties are identifiable. For the other three breeds with available turning points, these are Holstein (46.61), Simmental (48.88) and Original Braunvieh (49.71). Among them, the difference between the lowest and highest turning point is 3.41 THI points. Table 4.4 lists the individual THI turning points. Figure 4.5 and Figure 4.6 solely visualize the marginal effects of THI on ECM yield in separate subfigures for primi- and multiparous cows.

4.4 Split Period: Until 2010 - After 2010

Table 4.5 and Table 4.6 compare the THI turning points and loss rates for the split-period experiments. The associated figures are provided in Appendix C. We observe THI turning point shifts between the two periods. More specifically, the turning points are at lower THI values in the period 2011-2023. However, in many cases, leftward turning point shifts, come with lower linearized loss rates.

4.5 Discussion

General Observations

Figure 4.7 combines the findings of the full-period models presented in Section 4.2 and Section 4.3. Consistent with expectations, non-linear effects THI on milk yield and ECM yield have been identified. Specifically, the observed concave-downward responses of THI align with biological assumptions, indicating that the animals reach an optimal performance level at specific THI points, beyond which their productivity diminishes as THI levels further increase. For the majority of breeds, including both primiparous and multiparous cases, these THI turning points can be numerically determined. Additionally, clear differences in average performance levels across breeds are evident, with Holstein cattle exhibiting superior production, whereas dual-purpose breeds such as Simmental and Original Braunvieh exhibit lower performance metrics, contingent upon whether the yield is adjusted for energy content. Furthermore, it is possible to observe distinct THI turning points among different breeds. In general, the impact of THI on milk yield is moderate in relation to the average yields among different breeds. The effects of accumulated THI, commencing from the turning points, range from approximately 0 to 1.5 kg. None of the breeds appears to exhibit significantly superior performance over others due to THI effects. This observation is congruent with the findings of Citation: Ahmed et al., 2022 Ahmed, H., Tamminen, L. & Emanuelson, U. (2022). Temperature, productivity, and heat tolerance: Evidence from Swedish dairy production. , 175(1-2). 10. https://doi.org/10.1007/s10584-022-03461-5  that suggest diversification of a dairy cow portfolio does not improve resilience to heat. However, further research is required, employing simulations over extended time frames and including productivity parameters such as cost per cow and breed.

Low THI Values

An initial investigation reveals unexpectedly low THI values at the turning points. It is anticipated that the THI turning points would be situated in the high sixties or early seventies for both milk yield and ECM yield Citation: Armstrong, 1994 Armstrong, D. (1994). Heat Stress Interaction with Shade and Cooling. , 77(7). 2044–2050. https://doi.org/10.3168/jds.S0022-0302(94)77149-6 ; Citation: De Rensis et al., 2015 De Rensis, F., Garcia-Ispierto, I. & López-Gatius, F. (2015). Seasonal heat stress: Clinical implications and hormone treatments for the fertility of dairy cows. , 84(5). 659–666. https://doi.org/10.1016/j.theriogenology.2015.04.021 ; Citation: Vroege et al., 2023 Vroege, W., Dalhaus, T., Wauters, E. & Finger, R. (2023). Effects of extreme heat on milk quantity and quality. , 210. 103731. https://doi.org/10.1016/j.agsy.2023.103731  , which correspond to light heat stress and temperatures between 22°C and 33°C, depending on the relative humidity, according to Figure A.1. Contrary to this expectation, the peak THI values start in the low forties concerning ECM yield for primiparous cows and extend to the high sixties with regard to milk yield for multiparous cows. This implies dairy cows reach their peak performance with respect to THI independent of the lactation stage, at positive single-digit degrees Celsius in some cases. First, we did check our data pipeline multiple times for errors and also ran numerous sanity checks regarding the location-weather matching. Second, in Switzerland, the average farm are small-scale family operations. We speculated that this effect is due to the fact that in Switzerland the average farm is a small-scale family operation. When applying the same models to data solely from larger farms with more than 60 cows in 2023, supposedly well-managed, the results indicate an elevation in THI thresholds by approximately 1-2 points. Employing broken-stick fixed-effect models from econometrics, which account for animals and farms as unobserved heterogeneity and are estimated through pooled OLS regression with a search for the statistically optimal breaking point, yields similar THI threshold values. Furthermore, when using the same methodological framework but substituting mean THI values with maximum THI values, the threshold increases by several THI points; however, all remain significantly below the expected threshold of 68 THI points². These THI turning points have precedents: Citation: Hill & , 2015 Hill, D. & Wall, E. (2015). Dairy cattle in a temperate climate: The effects of weather on milk yield and composition depend on management. , 9(1). 138–149. https://doi.org/10.1017/S1751731114002456  report average THI values as low as 54.9 as threshold values for outdoor cows. Furthermore, Citation: Vinet et al., 2023 Vinet, A., Mattalia, S., Vallée, R., Bertrand, C., Cuyabano, B. & Boichard, D. (2023). Estimation of genotype by temperature-humidity index interactions on milk production and udder health traits in Montbeliarde cows. , 55. https://doi.org/10.1186/s12711-023-00779-1  identify optimal 3-day THI ranges between 50 and 55 for the Montbélliarde breed. Finally, our reference study Citation: Ahmed et al., 2022 Ahmed, H., Tamminen, L. & Emanuelson, U. (2022). Temperature, productivity, and heat tolerance: Evidence from Swedish dairy production. , 175(1-2). 10. https://doi.org/10.1007/s10584-022-03461-5  documents negative marginal effects on milk yield at temperatures of 4-5\(^{\circ}\)C based on a 7-day mean of daily maximum temperatures. Additionally, Citation: Schüller et al., 2013 Schüller, L., Burfeind, O. & Heuwieser, W. (2013). Short communication: Comparison of ambient temperature, relative humidity, and temperature-humidity index between on-farm measurements and official meteorological data. , 96(12). 7731–7738. https://doi.org/10.3168/jds.2013-6736  indicate that ambient temperatures recorded at official meteorological stations are notably lower when compared to temperatures within barns. They document temperature discrepancies of \(6.5 \pm 3.6^\circ\)C and corresponding variations of \(11 \pm 6.5\) THI points. Our data does not provide insights into the duration that cows spend indoors versus outdoors, nor which of the pasture-related policy frameworks outlined in Table 2.3 are applicable to specific farms.

THI Turning Points - Primiparous vs Multiparous Cows

We consistently observe lower turning points for primiparous cows than for multiparous cows. This contradicts the findings of Citation: Bernabucci et al., 2015 Bernabucci, U., Basiricò, L., Morera, P., Dipasquale, D., Vitali, A., Piccioli Cappelli, F. & Calamari, L. (2015). Effect of summer season on milk protein fractions in Holstein cows. , 98(3). 1815–1827. https://doi.org/10.3168/jds.2014-8788 ; Citation: West, 2003 West, J. (2003). Effects of Heat-Stress on Production in Dairy Cattle. , 86(6). 2131–2144. https://doi.org/10.3168/jds.S0022-0302(03)73803-X  and Citation: Becker et al., 2020 Becker, C., Collier, R. & Stone, A. (2020). Invited review: Physiological and behavioral effects of heat stress in dairy cows. , 103(8). 6751–6770. https://doi.org/10.3168/jds.2019-17929 , which identify an increased heat tolerance in primiparous cows. They state that primiparous cows, owing to their lower metabolic heat production and reduced milk yield, should exhibit an enhanced capacity for mitigating heat stress relative to multiparous cows. Nevertheless, Citation: Maggiolino et al., 2022 Maggiolino, A., Landi, V., Bartolomeo, N., Bernabucci, U., Santus, E., Bragaglio, A. & De Palo, P. (2022). Effect of Heat Waves on Some Italian Brown Swiss Dairy Cows’ Production Patterns. , 2. 800680. https://doi.org/10.3389/fanim.2021.800680  and Citation: Vinet et al., 2023 Vinet, A., Mattalia, S., Vallée, R., Bertrand, C., Cuyabano, B. & Boichard, D. (2023). Estimation of genotype by temperature-humidity index interactions on milk production and udder health traits in Montbeliarde cows. , 55. https://doi.org/10.1186/s12711-023-00779-1  report findings consistent with our study. One plausible explanation posits that primiparous cows are still undergoing maturation, and that their physiological systems, including thermoregulatory mechanisms, may not be as fully developed or efficient as those in older, multiparous cows.

THI Turning Points - Breeds

As anticipated, distinct THI turning points are observed for various breeds concerning both milk yield and ECM yield. Notably, Holstein cows present a lower THI threshold compared to Jersey and Simmental breeds, corroborating our hypotheses detailed in Section 1.3 at a first glance. In the context of milk yield, multiparous Jersey cows demonstrate the highest THI threshold at 58.12 points, accompanied by a relatively modest rate of loss when contrasted with other breeds. Upon evaluating the thresholds and the shape of the response curves, Jerseys are identified as possessing the highest resilience to heat stress. In contrast, Holstein and Simmental breeds exhibit the least tolerance to elevated temperatures. The Simmental example highlights that not only the THI turning points are relevant, but also the shapes of the curves. Hence, Brown Swiss exhibiting the lowest THI turning points, could still be qualified as more heat-robust than Simmental since the losses are considerably larger for Simmental with a rising THI values. Analyzing the response curve configurations between Holstein and Simmental breeds, Simmental exhibits superior robustness compared to Holstein, even though the differences are relatively small.

The THI thresholds are lower for ECM yield compared to milk yield. Consequently, the milk component yields in relation to THI induce a leftward shift. Nonetheless, this shift does not remain uniform across all breeds. While the disparity between the lowest and highest THI turning points for milk yield is 4.56 for primiparous cows and 8.16 THI units for multiparous cows, the differences for ECM yield are notably smaller. Apparent benefits for certain breeds with respect to milk volume dissipate upon adjusting the milk volume for fat and protein yields. Figure 4.8 and Figure 4.9 superpose the marginal effects of THI for each breed without confidence intervals to better illustrate this observation. The inverse relationship between THI and both fat and protein content is expected Citation: Moore et al., 2023 Moore, S., Costa, A., Penasa, M., Callegaro, S. & De Marchi, M. (2023). How heat stress conditions affect milk yield, composition, and price in Italian Holstein herds. , 106(6). 4042–4058. https://doi.org/10.3168/jds.2022-22640 . Citation: Vinet et al., 2023 Vinet, A., Mattalia, S., Vallée, R., Bertrand, C., Cuyabano, B. & Boichard, D. (2023). Estimation of genotype by temperature-humidity index interactions on milk production and udder health traits in Montbeliarde cows. , 55. https://doi.org/10.1186/s12711-023-00779-1  and Citation: Hill & , 2015 Hill, D. & Wall, E. (2015). Dairy cattle in a temperate climate: The effects of weather on milk yield and composition depend on management. , 9(1). 138–149. https://doi.org/10.1017/S1751731114002456  have determined that declines in fat content in relation to THI manifest at an earlier stage or are stronger than protein declines. The specific individual contributions of protein and fat per breed to the observed shift and dissipation require further investigation.

Confidence Intervals

The variability observed in the confidence intervals for each breed can be attributed to the employed subsampling strategy. The larger farms, the more samples per farm. Moreover, the available number of samples per farm also depends on the sampling strategy of the breeding associations. Table 4.1 and Table 4.2 support this statement. Consequently, our approach of incrementally incorporating farms into the dataset until a predetermined sample limit is reached results in a larger number of samples for these farms \(t_s\). As a result, this approach contributes to diminished THI diversity, given that all samples gathered on a particular test-day at a given farm are associated with a single THI data point. This raises the question whether alternative subsampling methodologies at the farm level might need to be employed. Nonetheless, a careful interpretation of confidence intervals is necessary in this non-experimental study design and is further elaborated in Chapter 5.

Split Period: Until 2010 - After 2010

A striking difference emerges when splitting the data into two time periods: the THI turning points tend to be lower in the more recent data period from 2011 to 2023 compared to the period from 1982 to 2023. The Jersey breed poses an exception for ECM yield.

The average daily milk production is substantially higher today for all breeds relative to 13 years ago, attributed to successful breeding programs. The reduced THI turning points suggest that selective breeding has made dairy cows less tolerant to heat. The principal objectives of breeding are to enhance dairy performance with equivalent or relatively diminished input. Furthermore, long-term genetic selection exerts an influence on the metabolic system. The current selection strategies may prioritize enhancements in dairy performance without necessarily fostering heat resilience. As breeding continues to result in increased yields, offsetting heat-induced losses despite heightened heat exposure for dairy cows in Switzerland, this development potentially prompts more concerns regarding animal welfare rather than economic implications Citation: König & , 2024 König, S. & Swalve, H. (2024). Breeding for metabolic robustness and longevity in dairy cows. InGross, J. (Eds.), . (pp. 531–553). Springer International Publishing. https://doi.org/10.1007/978-3-031-51788-4_23 . However, our results need further investigation with different cut-off years, multiple periods or different data subsampling strategies. Furthermore, shifts in THI turning points come along with changes in the loss rates. For example, simulations aggregating the daily effects with weather data over a long periods are necessary to further validate this negative change in heat-tolerance.

4.6 Bibliography