Advertisement

Updating predictions of dry matter intake of lactating dairy cows

Open ArchivePublished:July 17, 2019DOI:https://doi.org/10.3168/jds.2018-16176

      ABSTRACT

      Our objective was to model dry matter intake (DMI) by Holstein dairy cows based on milk energy (MilkE), body weight (BW), change in BW (ΔBW), body condition score (BCS), height, days in milk (DIM), and parity (primiparous and multiparous). Our database included 31,631 weekly observations on 2,791 cows enrolled in 52 studies from 8 states of the United States, mostly in the Upper Midwest. The means ± standard deviations of these variables were 24 ± 5 kg of DMI, 30 ± 6 Mcal of MilkE/d, 624 ± 83 kg of BW, 0.24 ± 1.50 kg of ΔBW/d, 3.0 ± 0.5 BCS, 149 ± 6 cm height, and 102 ± 45 DIM. Data analysis was performed using a mixed-effects model containing location, study within location, diet within study, and location and cow within study as random effects, whereas the fixed effects included the linear effects of the covariates described previously and all possible 2-way interactions between parity and the other covariates. A nonlinear (NLIN) mixed model analysis was developed using a 2-step approach for computational tractability. In the first step, we used a linear (LIN) model component of the NLIN model to predict DMI using only data from mid-lactation dairy cows (76–175 DIM) without including information on DIM. In the second step, a nonlinear adjustment for DIM using all data from 0 to 368 DIM was estimated. Additionally, this NLIN model was compared with an LIN model containing a fourth-order polynomial for DIM using data throughout the entire lactation (0–368 DIM) to assess the utility of an NLIN model for the prediction of DMI. In summary, a total of 8 candidate models were evaluated as follows: 4 ways to express energy required for maintenance (BW, BW0.75, BW adjusted for a BCS of 3, and BW0.75 adjusted for a BCS of 3) × 2 modeling strategies (LIN vs. NLIN). The candidate models were compared using a 5-fold across-studies cross-validation approach repeated 20 times with the best-fitting model chosen as the proposed model. The metrics used for evaluation were the mean bias, slope bias, concordance correlation coefficient (CCC), and root mean squared error of prediction (RMSEP). The proposed prediction equation was DMI (kg/d) = [(3.7 + parity × 5.7) + 0.305 × MilkE (Mcal/d) + 0.022 × BW (kg) + (−0.689 + parity × −1.87) × BCS] × [1 – (0.212 + parity × 0.136) × exp(−0.053 × DIM)] (mean bias = 0.021 kg, slope bias = 0.059, CCC = 0.72, and RMSEP = 2.89 kg), where parity is equal to 1 if the animal is multiparous and 0 otherwise. Finally, the proposed model was compared against the Nutrient Requirements of Dairy Cattle (2001) prediction equation for DMI using an independent data set of 9,050 weekly observations on 1,804 Holstein cows. The proposed model had smaller mean bias and RMSEP and higher CCC than the Nutrient Requirements of Dairy Cattle equation to predict DMI and has potential to improve diet formulation for lactating dairy cows.

      Key words

      INTRODUCTION

      Actual or accurately estimated DMI is essential for the formulation of diets to prevent underfeeding or overfeeding of nutrients and to promote efficient nutrient use (
      NRC
      Nutrient Requirements of Dairy Cattle.
      ). Current equations used to predict DMI of lactating dairy cows are based on existing databases that contain observed intake and predictor variables such as milk energy (MilkE), BW, and DIM. Consequently, the accuracy and precision of these predictions are dependent on the quality of the database and the use of appropriate statistical analysis to derive the prediction equations.
      Traditionally, the prediction of DMI in lactating North American dairy cows is based on milk production and the energy required for maintenance, which together represent the major energy expenditures of the lactating dairy cow (
      • Tempelman R.J.
      • Spurlock D.M.
      • Coffey M.
      • Veerkamp R.F.
      • Armentano L.E.
      • Weigel K.A.
      • de Haas Y.
      • Stamples C.R.
      • Connor E.E.
      • Lu Y.
      • VandeHaar M.J.
      Heterogeneity in genetic and nongenetic variation and energy sink relationships for residual feed intake across research stations and countries..
      ;
      • VandeHaar M.J.
      • Armentano L.E.
      • Weigel K.
      • Spurlock D.M.
      • Tempelman R.J.
      • Veerkamp R.
      Harnessing the genetics of the modern dairy cow to continue improvements in feed efficiency..
      ). Furthermore, due to the peculiarities of each stage of lactation (i.e., transition and early-lactation dairy cows;
      • Allen M.S.
      • Bradford B.J.
      • Oba M.
      Board-invited review: The hepatic oxidation theory of the control of feed intake and its application to ruminants..
      ), DIM must also be considered.
      The most recent edition of Nutrient Requirements of Dairy Cattle (
      NRC
      Nutrient Requirements of Dairy Cattle.
      ) includes an empirical equation to estimate DMI of lactating Holstein cows based on animal factors that could be easily measured or known, as follows:
      DMI (kg/d) = (0.372 × FCM + 0.0968 × BW0.75) × {1 − e[−0.192 × (WOL + 3.67)]},


      where FCM = 4% FCM (kg/d), BW is given in kilograms, and WOL = week of lactation (
      NRC
      Nutrient Requirements of Dairy Cattle.
      ). This prediction equation was based on data collected on 1,284 Holstein cows using the equation proposed by
      • Rayburn E.B.
      • Fox D.G.
      Variation in neutral detergent fiber intake of Holstein cows..
      and further modified by
      • Fox D.G.
      • Van Amburgh M.E.
      • Tylutki T.P.
      Predicting requirements for growth, maturity, and body reserves in dairy cattle..
      . Also, the nonlinear effect of week of lactation was based on
      • Roseler D.K.
      • Fox D.G.
      • Chase L.E.
      • Pell A.N.
      • Stone W.C.
      Development and evaluation of equations for prediction of feed intake for lactating Holstein dairy cows..
      .
      A more recent and much larger database of weekly DMI data from many research stations was used by
      • Lu Y.
      • Vandehaar M.J.
      • Spurlock D.M.
      • Weigel K.A.
      • Armentano L.E.
      • Staples C.R.
      • Connor E.E.
      • Wang Z.
      • Coffey M.
      • Veerkamp R.F.
      • de Haas Y.
      • Tempelman R.J.
      Modeling genetic and nongenetic variation of feed efficiency and its partial relationships between component traits as a function of management and environmental factors..
      to estimate the genetic and nongenetic components of regression of DMI on MilkE and maintenance as a function of various management and environmental factors. Similarly,
      • Tempelman R.J.
      • Spurlock D.M.
      • Coffey M.
      • Veerkamp R.F.
      • Armentano L.E.
      • Weigel K.A.
      • de Haas Y.
      • Stamples C.R.
      • Connor E.E.
      • Lu Y.
      • VandeHaar M.J.
      Heterogeneity in genetic and nongenetic variation and energy sink relationships for residual feed intake across research stations and countries..
      developed station-specific prediction equations for DMI but without distinguishing genetic from nongenetic components as in
      • Lu Y.
      • Vandehaar M.J.
      • Spurlock D.M.
      • Weigel K.A.
      • Armentano L.E.
      • Staples C.R.
      • Connor E.E.
      • Wang Z.
      • Coffey M.
      • Veerkamp R.F.
      • de Haas Y.
      • Tempelman R.J.
      Modeling genetic and nongenetic variation of feed efficiency and its partial relationships between component traits as a function of management and environmental factors..
      . Nevertheless, both studies used only data between 50 and 200 DIM (
      • Tempelman R.J.
      • Spurlock D.M.
      • Coffey M.
      • Veerkamp R.F.
      • Armentano L.E.
      • Weigel K.A.
      • de Haas Y.
      • Stamples C.R.
      • Connor E.E.
      • Lu Y.
      • VandeHaar M.J.
      Heterogeneity in genetic and nongenetic variation and energy sink relationships for residual feed intake across research stations and countries..
      ;
      • Lu Y.
      • Vandehaar M.J.
      • Spurlock D.M.
      • Weigel K.A.
      • Armentano L.E.
      • Staples C.R.
      • Connor E.E.
      • Wang Z.
      • Coffey M.
      • Veerkamp R.F.
      • de Haas Y.
      • Tempelman R.J.
      Modeling genetic and nongenetic variation of feed efficiency and its partial relationships between component traits as a function of management and environmental factors..
      ) and therefore might not be valid for predicting DMI outside that range of DIM. Additional data were available in those studies that could be used to better develop intake predictions across the lactation. In addition, many of the cows in these data sets also had measures of BCS. Because body fatness alters the concentration of circulating leptin, which in turn alters appetite, the inclusion of BCS should increase the accuracy of predicting voluntary DMI (
      • Block S.S.
      • Smith J.M.
      • Ehrhardt R.A.
      • Diaz M.C.
      • Rhoads R.P.
      • Van Amburgh M.E.
      • Boisclair Y.R.
      Nutritional and developmental regulation of plasma leptin in dairy cattle..
      ;
      • Allen M.S.
      Drives and limits to feed intake in ruminants..
      ).
      In addition to models based only on animal characteristics, several models (
      • Arnedal S.
      Predictions for voluntary dry matter intake in dairy cows.
      ;
      • Huhtanen P.
      • Rinne M.
      • Mäntysaari P.
      • Nousiainen J.
      Integration of the effects of animal and dietary factors on total dry matter intake of dairy cows fed silage-based diets..
      ;
      • Zom R.L.G.
      • Andre G.
      • van Vuuren A.M.
      Development of a model for the prediction of feed intake by dairy cows: 1. Prediction of feed intake..
      ;
      • Krizsan S.J.
      • Sairanen A.
      • Höjer A.
      • Huhtanen P.
      Evaluation of different feed intake models for dairy cows..
      ) have been proposed to predict DMI based on both animal and diet characteristics.
      • Jensen L.M.
      • Nielsen N.I.
      • Nadeau E.
      • Markussen B.
      • Nørgaard P.
      Evaluation of five models predicting feed intake by dairy cows fed total mixed rations..
      evaluated 5 DMI prediction models and showed that models containing both animal and nutritional effects would yield better predictions for DMI than models based only on animal effect. Although DMI prediction may be improved by using both animal and dietary factors as covariates, in this paper we focused only on animal characteristics. One reason to use only animal characteristics is that predictions of intake are sometimes desired when diet composition is unknown, and the effects of feed factors are largely accounted for when milk yield and BW change are used in empirical predictions of intake. In addition, effects of feed factors differ depending on stage of lactation (
      • Allen M.S.
      Drives and limits to feed intake in ruminants..
      ). Predictions of DMI are especially important in formulating diets. Our overall idea is that the starting prediction for DMI when formulating diets should first be a prediction based on animal factors or on measured values, and then feed factors should be used to assess the effect of ration changes on DMI during the ration formulation process. A new equation to predict changes in DMI based on the filling effects of rations is described in a companion paper (
      • Allen M.S.
      • Sousa D.O.
      • VandeHaar M.J.
      Equation to predict feed intake response by lactating cows to factors related to the filling effect of rations..
      ).
      The objective of this study was to derive a prediction equation for DMI based on animal factors that include not only BW, MY, milk composition, DIM, and parity as in
      NRC
      Nutrient Requirements of Dairy Cattle.
      but also BCS and change in BW, which were not available for
      NRC
      Nutrient Requirements of Dairy Cattle.
      . Our data comprised an updated version of the database used in
      • Tempelman R.J.
      • Spurlock D.M.
      • Coffey M.
      • Veerkamp R.F.
      • Armentano L.E.
      • Weigel K.A.
      • de Haas Y.
      • Stamples C.R.
      • Connor E.E.
      • Lu Y.
      • VandeHaar M.J.
      Heterogeneity in genetic and nongenetic variation and energy sink relationships for residual feed intake across research stations and countries..
      that included 2,791 cows between 1 and 368 DIM. We expect that this database would better represent the current genetics and management practices of the North American dairy industry, although we recognize that it contains mostly cows from the Upper Midwest. We hypothesized that an equation derived from this larger and newer database would increase the accuracy and precision of DMI predictions based on animal factors compared with the equation provided in
      NRC
      Nutrient Requirements of Dairy Cattle.
      when tested in an independent data set.

      MATERIALS AND METHODS

      Data

      The database consisted of individual weekly values of DMI, MilkE, BW, metabolic BW (BW0.75), BCS, change in BW (ΔBW), height, DIM, and parity (primiparous and multiparous) from 10 research stations across the eastern half of the United States. We calculated BW0.75 as BW to the power of 0.75 and ΔBW as the final BW minus the initial BW through the experimetal period divided by number of days. Additionally, we calculated BW adjusted to a BCS of 3 (BWBCS3; kg) and BW0.75 adjusted to a BCS of 3 (BW0.75BCS3; kg) as a function of BW and BCS: BWBCS3 = BW + BW × 0.084 × (3 − BCS), and BW0.75BCS3 = BWBCS30.75. To calculate the BWBCS3, it was assumed that a 1-unit shift in BCS causes a 8.4% change in BW, as determined by
      • de Souza R.A.
      • VandeHaar M.J.
      Determining the change in body weight per unit of body condition score in Holstein cows..
      using a database comprising 2,181 Holstein cows. The MilkE was calculated according to
      NRC
      Nutrient Requirements of Dairy Cattle.
      .
      Recording frequencies of the variables included in the database varied from daily, thrice weekly, twice weekly, weekly, and biweekly depending on the variable and the research station. To standardize variables as weekly means across all studies, each variable (i.e., MilkE, BW, BCS, ΔBW, and height) on each animal on the same treatment for each study was analyzed using a fifth-order polynomial linear regression model on days since the start of that particular treatment. Subsequently, the predicted values from these animal-specific analyses were used to calculate weekly means for each animal on a particular treatment.
      The final database comprised 40,681 weekly observations on 4,785 Holstein dairy cows from 105 studies (Table 1). About 54% (22,045 weekly observations) of these data was previously described by
      • Tempelman R.J.
      • Spurlock D.M.
      • Coffey M.
      • Veerkamp R.F.
      • Armentano L.E.
      • Weigel K.A.
      • de Haas Y.
      • Stamples C.R.
      • Connor E.E.
      • Lu Y.
      • VandeHaar M.J.
      Heterogeneity in genetic and nongenetic variation and energy sink relationships for residual feed intake across research stations and countries..
      , but we used only their data from US research stations (University of Florida, Gainesville; Iowa State University, Ames; Michigan State University, East Lansing; University of Wisconsin, Madison; United States Dairy Forage Research Center, Madison, WI; and USDA Animal Genomics and Improvement Laboratory, Beltsville, MD) to be consistent with other data added to this database. The remaining 46% (18,636 weekly observations) of the data came either from locations already used in
      • Tempelman R.J.
      • Spurlock D.M.
      • Coffey M.
      • Veerkamp R.F.
      • Armentano L.E.
      • Weigel K.A.
      • de Haas Y.
      • Stamples C.R.
      • Connor E.E.
      • Lu Y.
      • VandeHaar M.J.
      Heterogeneity in genetic and nongenetic variation and energy sink relationships for residual feed intake across research stations and countries..
      but collected in more recent studies or collected at DIM less than 50 or greater than 200 (50% of all data) at these same stations or from new locations (14% of the data) including the Purina Animal Nutrition Center (Gray Summit, MO), Cargill Research and Development Center (Minneapolis, MN), Miner Institute (Chazy, NY), and Virginia Tech University (Blacksburg, VA). Our final database thereby comprised records from cows between 1 and 368 DIM and collected in research stations across the United States (Florida, Iowa, Maryland, Michigan, Missouri, New York, Virginia, and Wisconsin) from 2007 to 2016. All the data were collected from lactating Holstein dairy cows offered a TMR once per day and milked 2 or 3 times per day.
      Table 1Frequency of number of studies, diets, cows, and lactations by state
      Station
      AGIL = USDA Animal Genomics Improvement Laboratory, Beltsville, Maryland; UF = University of Florida, Gainesville; ISU = Iowa State University, Ames; MSU = Michigan State University, East Lansing; CRDC = Cargill Research and Development Center, Minneapolis, Minnesota; PANC = Purina Animal Nutrition Center, Gray Summit, Missouri; Miner = Miner Institute, Chazy, New York; VT = Virginia Tech, Blacksburg; UW = University of Wisconsin–Madison.
      No. of studiesNo. of dietsNo. of cowsNo. of lactationsNo. of weekly recordsAverage no. of weekly records per lactation
      AGIL234195865,4319.3
      UF11363664203,3267.9
      ISU349521,00810,33510.3
      MSU21803395405,68910.5
      CRDC1222373611,3313.68
      PANC13676796114.3
      Miner1458583405.9
      VT6221101118427.6
      UW33951,0881,18312,10610.5
      Total792693,6364,33440,6819.4
      1 AGIL = USDA Animal Genomics Improvement Laboratory, Beltsville, Maryland; UF = University of Florida, Gainesville; ISU = Iowa State University, Ames; MSU = Michigan State University, East Lansing; CRDC = Cargill Research and Development Center, Minneapolis, Minnesota; PANC = Purina Animal Nutrition Center, Gray Summit, Missouri; Miner = Miner Institute, Chazy, New York; VT = Virginia Tech, Blacksburg; UW = University of Wisconsin–Madison.
      The final database was subdivided into 2 independent data sets. The first data set, consisting of 31,631 weekly observations from 3,143 lactations (1,462 primiparous and 1,681 multiparous) on 2,791 cows, was used in the model development process (referred to as the development dataset), whereas the second data set, consisting of 9,050 weekly observations from 1,968 lactations (908 primiparous and 1,060 multiparous) on 1,804 cows, was used in the independent evaluation process. The distribution of the number of observations in each data set throughout the lactation is presented in Figure 1, and their respective means, medians, quartiles, maximums, and minimums are presented in Figure 2. The criteria used to determine the separation of the database into 2 data sets (model development and independent evaluation data sets) were based on the number of weekly observations per cow within a study. Studies that provided 7 or more weekly observations per cow constituted the model development data set, whereas the remaining smaller studies constituted the independent evaluation data set. Basic summary statistics of the raw weekly data for each of the 2 data sets are provided in Table 2 and shown in Figure 2.
      Figure thumbnail gr1
      Figure 1Distribution of the number of observations by DIM in the development (Dev) and evaluation (Eval) data sets.
      Figure thumbnail gr2
      Figure 2Means (×), medians (solid lines), quartiles (middle 2 quartiles in shaded boxes), and accepted ranges (whiskers) for DMI, milk energy (MilkE) output, BW, metabolic BW (MBW), BW change, BCS, height, and DIM for model development (Dev; n = 31,631) and evaluation (Eval; n = 9,050). The development data set is in blue, and the evaluation data set is in red. Whiskers extend up from the top of the box to the largest data element that is less than 1.5 times the interquartile range (IQR) and down from the bottom of the box to the smallest data element that is larger than 1.5 times the IQR. Anything outside of this range is considered to be an outlier (blue dots).
      Table 2Number of observations, mean, and standard deviation of the covariates used to develop (development data set) and evaluate (evaluation data set) the proposed equations for DMI
      Variable
      MilkE = milk energy; BW0.75 = metabolic BW; ΔBW = change in BW.
      Development data setEvaluation data set
      No. of observationsMeanSDNo. of observationsMeanSD
      DMI, kg/d31,63124.44.569,05024.04.58
      MilkE, Mcal/d31,63129.86.289,05028.96.72
      BW, kg31,63162479.49,05063988.2
      BW0.75, kg0.7531,63112511.99,05012713.1
      ΔBW, kg/d28,7290.2371.247,309−0.0032.08
      BCS
      Scale of 1 to 5.
      31,6313.040.4579,0502.990.429
      Height, cm7,6911495.284,1291504.44
      DIM31,63110846.79,05011465.6
      1 MilkE = milk energy; BW0.75 = metabolic BW; ΔBW = change in BW.
      2 Scale of 1 to 5.

      Statistical Analyses

      All data analysis was performed using SAS version 9.4 (SAS Institute Inc., Cary, NC). For numerical stability, all covariates were standardized to have a mean of zero and a standard deviation of 1 using PROC STANDARD and the functions MEAN = 0 and STD = 1; nevertheless, all reported prediction equations in this paper are provided on the actual observed or measured scale. All covariates were jointly checked with PROC REG for multicollinearity using variance inflation factors (multicollinearity analysis).
      In this paper, we refer to the models that initiated the modeling processes as the starting models, the models that were selected by the modeling processes and were used in cross-validation as the candidate models, and the model selected for best fit based on cross-validation and further evaluated against the
      NRC
      Nutrient Requirements of Dairy Cattle.
      DMI model as the proposed model. These models are more carefully explained in the subsequent sections.

      Modeling Process

      The nonlinear (NLIN) modeling process was accomplished in 2 major phases. The first phase was the development of a linear (LIN) model to predict DMI in mid-lactation cows (76–175 DIM) without the inclusion of DIM as a covariate. The second phase was the inclusion of the nonlinear effect of DIM on DMI across the entire modeling data set, including early (1–75 DIM), mid (76–175 DIM), and late (176–368 DIM) lactation.
      For the first phase of the modeling process, we specified a starting model that contained the fixed effects of partial linear regressions on MilkE, the corresponding BW measure (BW, BW0.75, BWBCS3, or BW0.75BCS3), ΔBW, BCS, parity (primiparous and multiparous), and all possible 2-way interactions between parity and the other covariates. Location, study within location, diet within study and location, and cow within study and location were treated as random effects. Furthermore, the model also included diet-specific partial regressions of DMI on MilkE and on BW, recognizing that the partial regressions of DMI on MilkE and on BW could depend on diet compositions. The generic starting model is presented in Equation 1:
      DMIpcdsl = β0 + u0,d + (β1 + u1,d) × MilkE + (β2 + u2,d) × body weight + β3 × ΔBW + β4 × BCS + β5 × Parityp + β6 × MilkE × Parityp + β7 × body weight × Parityp + β8 × ΔBW × Parityp + β9 × BCS × Parityp + Cowc(Studys × Locationl) + Studys(Locationl) + εpcdsl,
      [1]


      where DMIpcdsl (p = parity, c = cow, d = diet, s = study, l = location) is the observed DMI; β0 + u0,d is the intercept specific to diet d (d = 1–153); β1 + u1,d is the partial regression coefficient of DMI on MilkE specific to diet d; β2 + u2,d is the partial regression coefficient of DMI on body weight (BW, BW0.75, BWBCS3, or BW0.75BCS3) specific to diet d; β3 is the partial regression coefficient of DMI on ΔBW; β4 is the partial regression coefficient of DMI on BCS; β5 is the partial regression coefficient of DMI on parity (p = primiparous or multiparous); and β6, β7, β8, and β9 are the effects for the 2-way interactions between the corresponding covariates and parity. All other terms were random effects, including Cowc~NIID(0,σcow2) for cows c = 1 to 2,805 and Studys~NIID(0,σstudy2) for studies s = 1 to 57 across the various locations, with with~NIID(0,σcow2) being the error term assuming independent and identically distributed random variables with a mean of 0, and variance σ2. We assumed a multivariate normal distribution for diet effects and diet-specific partial regressions on MilkE and BW as follows:
      [u0,du1,du2,d]~N([000],G=[σu02σu0,u1σu0,u2σu0,u1σu12σu1,u2σu1,u2σu1,u2σu22]).
      [2]


      Note that the specification in Equation 2 also allowed for nonzero covariances, σu0,u1, σu0,u2 and σu1,u2 (i.e., off-diagonal values of the variance-covariance matrix G) between diet-specific intercepts and slopes for MilkE, between diet-specific intercepts and slopes for BW, and between diet-specific slopes for MilkE and BW, respectively.
      The starting models described above Equation 1 were subjected to model selection, which was accomplished in 2 phases. During the first phase, each of the full models was subjected to a backward model selection method using PROC GLMSELECT in SAS with the option HIERARCHY = SINGLE, such that the model was forced to include the corresponding main effects if the interaction between any 2 effects was deemed to be significant. Because PROC GLMSELECT does not allow for the specification of random effects, the model was limited to the specification of fixed effects only in this first stage of model selection. During the second phase of model selection, the model chosen by PROC GLMSELECT was analyzed using PROC HPMIXED, in which the above-mentioned random effects were included. The fixed effects with the highest P-values were successively removed until only significant fixed effects (P < 0.05) and all specified random effects remained in the model, yielding a linear mixed model as a candidate model to predict DMI for mid-lactation dairy cows. With this model selection process, we created 4 candidate models (i.e., 4 different ways to express the covariate BW; namely, BW, BW0.75, BWBCS3, and BW0.75BCS3).
      The NLIN model involved the estimation of the nonlinear effect of DIM on DMI using PROC NLMIXED in SAS. The nonlinear adjustment was based on the adjustment factor proposed by
      • Roseler D.K.
      • Fox D.G.
      • Chase L.E.
      • Pell A.N.
      • Stone W.C.
      Development and evaluation of equations for prediction of feed intake for lactating Holstein dairy cows..
      ; that is, 1-exp-[(β1-β2×PMM)×WOL+β3], where PMM is peak milk month and WOL is week of lactation, with further modifications as follows: (1) effect of parity was included to allow for parity specific effects, (2) PMM was removed because information was not available in our data set, and (3) a term was included that determines the maximum possible discount (β1 in Equation 3). The generic model used during the second phase is presented in Equation 3:
      DMI=DMIˆc×{1[β1+β2×Parity]×exp[(β3+β4×Parity)×DIM]}+εc,
      [3]


      where DMI is the observed DMI, DMIcˆ is the predicted DMI from candidate linear submodel c (c = 1–4) based on the 4 different ways to express the covariate BW, β1 is the coefficient that determines the maximum possible discount for DIM for primiparous cows, β2 specifies the difference for this discount between multiparous and primiparous cows, β3 is the nonlinear regression coefficient on DIM for primiparous cows, and β4 specifies the difference for this coefficient between multiparous and primiparous cows. Furthermore, parity is the dummy variable (0 for primiparous, 1 for multiparous) for parity in Equation 3, whereas ~NIID(0,σe2) is the error term.
      To verify the importance of using a nonlinear adjustment for DIM in the prediction of DMI, we refitted the 4 linear candidate models (based on model selection) to the entire development data set by adding a fourth-order polynomial for DIM and all possible 2-way interactions with parity as fixed effects as part of an LIN modeling. Hence, our model comparison involved a total of 8 candidate models [4 ways to express BW × model (NLIN vs. LIN)].

      Cross-Validation

      The 8 candidate models were formally compared using a 5-fold across-study cross-validation repeated 20 times for a total of 100-fold replicates using the development data set of 31,631 weekly observations (development data set described in the “Data” section). That is, for each of the 20 cross-validation replicates, the data were partitioned across rather than within studies into 5 almost-equal-sized subsets, each involving all data from 20% of the studies, such that for 1 fold, 4 of these subsets were used for development and the remaining subset was used for evaluation. Hence, across 5 folds, each subset took a turn being the evaluation data set. The candidate models were compared for fit criteria according to
      • Tedeschi L.O.
      Assessment of the adequacy of mathematical models..
      , including mean bias, slope bias, concordance correlation coefficient (CCC), mean squared error of prediction (MSEP), root MSEP (RMSEP), and the decomposition of MSEP in mean bias, slope bias, and random error (as percentage of MSEP). To determine our proposed model, the fit criteria were analyzed using PROC GLIMMIX of SAS version 9.4 according to the following model (Equation 4):
      Yijl = µ + ri + BWj + Modelm + BW × Modeljm + eijl,
      [4]


      where Yijl is one of the fit statistics of interest; µ is the overall mean; ri is the random effect of fold (i = 1–100); BWj is the fixed effect of corresponding BW measure (j = BW, BW0.75, BWBCS3, or BW0.75BCS3), Modelm is the fixed effect of model (m = LIN and NLIN), BW × Modeljm is the interaction between corresponding BW and model, and eijl is the residual error. The candidate model that consistently had the best goodness of fit was considered the best-fitting model and, therefore, our proposed model. A generic SAS code for model evaluation (using only linear models) and a short explanation are available in Supplemental Material S1 (https://doi.org/10.3168/jds.2018-16176).

      Independent Evaluation Process

      As an independent evaluation process beyond the cross-validation strategy described previously, the proposed model was compared with the DMI prediction equation proposed by
      NRC
      Nutrient Requirements of Dairy Cattle.
      using the independent evaluation data set of 9,050 weekly observations involving smaller studies with 6 or fewer weekly observations per cow. That is, for each of 20 evaluation replicates, the data were partitioned across studies into 5 almost-equal-sized subsets, each subset involving all data from 20% of the studies. For both models (proposed and NRC), the predictions of DMI were computed using the same random effects model. For each fold, the fit statistics (CCC, MSEP, RMSEP, and decomposition of MSEP) were computed as described previously. The major difference between the previously described cross-validation study and this independent evaluation study was that the same entire development data set of 31,631 weekly observations involving all of the larger studies (≥7 records per cow) was used for each evaluation subset such that predictions for each evaluation subset were based on the same fixed effects and variance components as based on the proposed model. Furthermore, to facilitate stage-specific comparisons, the evaluation data set was partitioned further according to the stage of lactation (early = 1–75 DIM; mid = 76–175 DIM; late = 176–368 DIM). The various fit statistics were compared using PROC GLIMMIX of SAS version 9.4 according to the following model (Equation 5):
      Yijk = µ + ri + Modelj + Stagek + Model × Stagejk + eijk,
      [5]


      where Yijk is the fit statistic of interest; µ is the overall mean; ri is the random effect of fold (i = 1–100); Modelj is the fixed effect of model (j = proposed or NRC), Stagek is the fixed effect of stage of lactation (k = early, mid, or late), Model × Stagejk is the fixed effect of the interaction between model and stage, and eijk is the residual error.

      RESULTS

      Regarding the 4 different measures of expressing BW (BW, BW0.75, BWBCS3, or BW0.75BCS3), we did not observe any difference in their ability to predict DMI based on the cross-validation (P > 0.9; results not shown). We did not observe differences between the 4 ways to express BW, possibly because all measures were highly correlated with each other (r = 0.96–0.99) and all cows in our data set were Holsteins. Thus, for the remainder of this paper, we show and discuss only data for the models developed using BW due to its simplicity and applicability.
      During the model development process, the main effect of ΔBW and the interactions of parity with MilkE, BW, and ΔBW were removed from the starting model because they were not significant (P > 0.05). With respect to the nonlinear effect of DIM in the NLIN model, parity was an important (P < 0.01) source of variation for the maximum discount on DMI due to DIM (i.e., β2 in Equation 3), but the effect of parity on the estimated nonlinear coefficient for DIM (i.e., β4 in Equation 3) was not significant (P = 0.68) and therefore was removed.
      Table 3Fit statistics for the cross-validation across studies of the candidate models (linear and nonlinear) using the modeling data set
      Fit statistic
      CCC = concordance correlation coefficient; MSEP = mean squared error of prediction; RMSEP = root MSEP.
      Candidate model
      The 2 modeling strategies, where the LIN model is a linear model including a fourth-order polynomial of DIM and the NLIN model contains the nonlinear effect of DIM. NS = not significantly different from zero.
      SEM
      SEM common to both LIN and NLIN.
      P-value
      For the mean comparison of LIN and NLIN.
      LINNLIN
      CCC0.7510.8000.006<0.01
      Mean bias−0.530.008NS0.059<0.01
      Slope bias−0.050NS−0.033NS0.0090.03
      MSEP10.96.830.143<0.01
      RMSEP3.202.610.073<0.01
      Decomposition of MSEP, %
       Mean bias13.66.900.616<0.01
       Slope bias1.962.280.1840.09
       Random error84.490.80.647<0.01
      1 CCC = concordance correlation coefficient; MSEP = mean squared error of prediction; RMSEP = root MSEP.
      2 The 2 modeling strategies, where the LIN model is a linear model including a fourth-order polynomial of DIM and the NLIN model contains the nonlinear effect of DIM. NS = not significantly different from zero.
      3 SEM common to both LIN and NLIN.
      4 For the mean comparison of LIN and NLIN.
      The NLIN model had better goodness-of-fit properties than the LIN model for all fit statistics, and therefore we selected a nonlinear relationship for the effect of DIM on DMI (Table 3). Our proposed NLIN model is shown in Equation 6:
      DMI (kg/d) = [(3.7 + Parity × 5.7) + 0.305 × MilkE (Mcal/d) + 0.022 × BW (kg) + (−0.689 + Parity × −1.87) × BCS] × [1 − (0.212 + Parity × 0.136) × e(−0.053 × DIM)],
      [6]


      where Parity is equal to 1 if the animal is multiparous and 0 if primiparous.
      The diagonal values of the estimated G matrix, which represent the estimated variance components for the diet-specific intercept, MilkE, and BW were σˆu02 = 0.56, σˆu12 = 0.056, and σˆu22 = 0.000247, respectively. The estimated DMI for multiparous animals adjusted to the mean values of MilkE, BW, BCS, and DIM as reported in Table 2 and using Equation 6 was 24.5 kg. Based on the multivariate normality assumption as invoked in Equation 3, this implies that roughly 95% of the diets lead to a mean DMI within 24.5±2×σˆu02 (or 23.0–26.0 kg) for multiparous animals adjusted to mean values of MilkE, BW, BCS, and DIM. Similarly, 95% of the diets would have a partial regression of DMI on MilkE of 0.305±2×σˆu12 (or −0.17 to 0.78 kg/Mcal; value of 0.305 from Equation 6), whereas 95% of the diets would have a partial regression of DMI on BW of 0.022±2×σˆu22 (or −0.01 to 0.053 kg/kg; value of 0.022 from Equation 6). This variability across diets naturally reflects, for example, differences in the filling effects and energy density across different diets whereas negative lower bounds likely reflect that these diet effects are either not precisely normally distributed or based on relatively small numbers of observations. The remaining variance component estimates were 3.39 for σcow2 4.41 for σstudy2, and 3.24 for σe2, implying that these were even larger sources of variability relative to diets within studies. Thus, as a proportion of the total variance, cow, diet within study, study, and residual variance accounted for 29, 5, 38, and 28% of the total variation, respectively, relative to cows of average MilkE and BW.
      For the comparison of the proposed model (Equation 6) and the
      NRC
      Nutrient Requirements of Dairy Cattle.
      DMI prediction equation, the effects of model, stage, and interaction between model and stage were highly significant (P < 0.01; Table 4) for all fit statistics. First, related to the effect of model, the proposed NLIN model had better goodness of fit than the
      NRC
      Nutrient Requirements of Dairy Cattle.
      model for all fit statistics (Table 4). Regarding the effect of stage, the prediction equations performed better for cows during mid lactation (CCC, slope bias, MSEP, RMSEP), followed by early lactation and last by late lactation. The interaction between model and stage of lactation was significant, where the goodness-of-fit properties were slightly better for our proposed model than the
      NRC
      Nutrient Requirements of Dairy Cattle.
      model during early lactation but much better during mid and late lactation (Table 4). The predicted versus observed plot for each stage of lactation is presented in Figures 3, 4, and 5, with data points for weekly DMI shown as measured and after adjustment for all random effects. The proposed model had slightly less mean slope bias than did the NRC model, but both models tended to overpredict intakes for cows at high intakes.
      Table 4Fit statistics for the across-study stage of lactation
      Early = 1–75 DIM; mid = 76–175 DIM; late = 176–368 DIM.
      (early, mid, and late)-specific cross-validation performance of the proposed and
      NRC
      Nutrient Requirements of Dairy Cattle.
      models using the independent evaluation data set
      Fit statistic
      CCC = concordance correlation coefficient; MSEP = mean squared error of prediction; RMSEP = root MSEP.
      Proposed model
      NS = not significantly different from zero.
      NRC
      Nutrient Requirements of Dairy Cattle.
      model
      SEM
      SEM common to all models.
      P-value
      EarlyMidLateEarlyMidLateModelStageModel × stage
      CCC0.717
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      0.735
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      0.688
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      0.702
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      0.681
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      0.645
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      0.010<0.01<0.010.01
      Mean bias−0.281
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      0.058
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      0.293
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      −0.801
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      1.54
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      1.45
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      0.012<0.01<0.01<0.01
      Slope bias
      Slope of the residuals for observed versus predicted intake.
      −0.131
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      −0.028
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      −0.068
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      −0.137
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      −0.117
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      −0.116
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      0.014<0.01<0.01<0.01
      MSEP10.2
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      7.13
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      8.42
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      11.8
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      10.5
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      11.3
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      0.339<0.01<0.01<0.01
      RMSEP3.19
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      2.67
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      2.90
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      3.42
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      3.23
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      3.32
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      0.048<0.01<0.01<0.01
      Decomposition of MSEP, %
       Mean bias7.15
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      3.45
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      8.02
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      12.0
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      25.9
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      24.2
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      1.57<0.01<0.01<0.01
       Slope bias2.73
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      0.74
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      3.91
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      5.24
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      2.28
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      1.32
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      0.5640.13<0.01<0.01
       Random error90.2
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      95.8
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      88.1
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      82.7
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      71.8
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      74.5
      Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      1.60<0.01<0.01<0.01
      A,B Different uppercase superscripts within fit statistic represent statistical differences for model. P < 0.05.
      a–c Different lowercase superscripts within fit statistic represent statistical differences for stage. P < 0.05.
      1 Early = 1–75 DIM; mid = 76–175 DIM; late = 176–368 DIM.
      2 CCC = concordance correlation coefficient; MSEP = mean squared error of prediction; RMSEP = root MSEP.
      3 NS = not significantly different from zero.
      4 SEM common to all models.
      5 Slope of the residuals for observed versus predicted intake.
      Figure thumbnail gr3
      Figure 3Plots of the observed weekly DMI of cows in early lactation (1–75 DIM) versus that predicted by
      NRC
      Nutrient Requirements of Dairy Cattle.
      (blue dots, dotted line) and the proposed model (red triangles, dashed line) with individual weekly DMI points as measured (a) and after adjustment for random effects (b). Systematic prediction error can be observed by comparing models with the line of equivalence (solid line).
      Figure thumbnail gr4
      Figure 4Plots of the observed weekly DMI of cows in mid lactation (76–175 DIM) versus that predicted by
      NRC
      Nutrient Requirements of Dairy Cattle.
      (blue dots, dotted line) and the proposed model (red triangles, dashed line) with individual weekly DMI points as measured (a) and after adjustment for random effects (b). Systematic prediction error can be observed by comparing models with the line of equivalence (solid line).
      Figure thumbnail gr5
      Figure 5Plots of the observed weekly DMI of cows in late lactation (>175 DIM) versus that predicted by
      NRC
      Nutrient Requirements of Dairy Cattle.
      (blue dots, dotted line) and the proposed model (red triangles, dashed line) with individual weekly DMI points as measured (a) and after adjustment for random effects (b). Systematic prediction error can be observed by comparing models with the line of equivalence (solid line).
      On dairy farms, diets are fed to groups of animals that contain a mix of primiparous and multiparous cows. In this case, Equation 6 can be modified to represent the proportion between primiparous and multiparous in the group. For example, in our development data set the ratio between primiparous and multiparous was 44:66; therefore, to use the equation for the entire development data set, the estimated coefficient for parity will be multiplied by 0.66 (proportion of multiparous cows in the data set). As a result, the prediction equation for DMI for the average cow in the development data set is presented in Equation 7. Equation 7 was named the simplified proposed model:
      DMIˆ=(6.89+0.305×MilkE+0.022×BW1.74×BCS)×[10.288×e(0.053×DIM)].
      [7]


      The difference in the prediction of DMI using the proposed model (Equation 6) and the simplified proposed model (Equation 7) is presented in Table 5. We recommend using Equation 6 rather than Equation 7. However, when parity information is not known, Equation 7 could be used as an alternative.
      Table 5Fit statistics for the cross-validation across studies of the proposed model and simplified model using the evaluation data set
      Fit statistic
      CCC = concordance correlation coefficient; MSEP = mean squared error of prediction; RMSEP = root MSEP.
      Model
      The proposed model with parity effect and the simplified model using the proportion of the primiparous:multiparous (44:56) cows in the development data set. NS = not significantly different from zero.
      SEM
      SEM common to both proposed and simplified.
      P-value
      For the mean comparison of proposed and simplified.
      ProposedSimplified
      CCC0.7210.7140.0060.52
      Mean bias0.0610.0650.0040.47
      Slope bias−0.039NS−0.042NS0.0090.01
      MSEP8.369.050.2410.01
      RMSEP2.893.010.0380.01
      Decomposition of MSEP, %
       Mean bias5.706.320.6670.02
       Slope bias2.112.150.2180.85
       Random error92.291.50.6460.02
      1 CCC = concordance correlation coefficient; MSEP = mean squared error of prediction; RMSEP = root MSEP.
      2 The proposed model with parity effect and the simplified model using the proportion of the primiparous:multiparous (44:56) cows in the development data set. NS = not significantly different from zero.
      3 SEM common to both proposed and simplified.
      4 For the mean comparison of proposed and simplified.

      DISCUSSION

      This study shows that the
      NRC
      Nutrient Requirements of Dairy Cattle.
      model does a reasonably good job of predicting intake in Holstein cows, but the proposed model had better goodness of fit than the NRC model in the independent evaluation analysis, especially during mid and late lactation. We suspect that this is because we had data from more cows and included parity and BCS in our prediction. As can be seen in Table 4 and in Figures 3, 4, and 5, both the NRC and proposed models show some slight systematic error with negative slope bias. Negative slope bias was worst in early lactation and indicates that the model overpredicts intake for cows with the greatest intakes, possibly because the filling effects of diets limit intake most in those cows with the greatest requirements (
      • Jensen L.M.
      • Nielsen N.I.
      • Nadeau E.
      • Markussen B.
      • Nørgaard P.
      Evaluation of five models predicting feed intake by dairy cows fed total mixed rations..
      ;
      • Allen M.S.
      • Sousa D.O.
      • VandeHaar M.J.
      Equation to predict feed intake response by lactating cows to factors related to the filling effect of rations..
      ).
      The DIM adjustment for the prediction of DMI is necessary due to the peculiarities of different postpartum stages over and beyond the effects of MilkE or BW (
      • Dann H.M.
      • Varga G.A.
      • Putnam D.E.
      Improving energy supply to late gestation and early postpartum dairy cows..
      ;
      • Drackley J.K.
      Biology of dairy cows during the transition period: The final frontier?.
      ;
      • Chan P.S.
      • West J.W.
      • Bernard J.K.
      Effect of prepartum dietary calcium on intake and serum and urinary mineral concentrations of cows..
      ;
      • Janovick N.A.
      • Boisclair Y.R.
      • Drackley L.K.
      Prepartum dietary energy intake affects metabolism and health during the periparturient period in primiparous and multiparous Holstein cows..
      ;
      • Allen M.S.
      Drives and limits to feed intake in ruminants..
      ). For the proposed model, as in
      • Roseler D.K.
      • Fox D.G.
      • Chase L.E.
      • Pell A.N.
      • Stone W.C.
      Development and evaluation of equations for prediction of feed intake for lactating Holstein dairy cows..
      , the greatest discounts in DMI are early in lactation and are not important after 70 DIM (Table 6). Multiparous animals had proportionately greater deviance in the effect of DIM on DMI than primiparous cows in early lactation (Table 6). As an example, at the onset of lactation (1 DIM), multiparous cows had a discount of 21% on the predicted DMI based on MilkE, BW, and BCS, whereas primiparous cows had a discount of 13%.
      Table 6Nonlinear adjustment of DIM on the predicted DMI estimated by the proposed (primiparous and multiparous) and
      NRC
      Nutrient Requirements of Dairy Cattle.
      models
      DIMProposed modelNRC model
      PrimiparousMultiparous
      100.8750.7950.624
      200.9270.8800.714
      300.9570.9290.783
      400.9750.9580.835
      500.9850.9760.875
      600.9910.9860.905
      700.9950.9920.928
      800.9970.9950.945
      900.9980.9970.958
      1000.9990.9980.968
      1100.9990.9990.976
      1201.0000.9990.982
      1301.0001.0000.986
      1401.0001.0000.989
      1501.0001.0000.992
      The
      NRC
      Nutrient Requirements of Dairy Cattle.
      applies a more aggressive discount for DIM on the DMI than the proposed model as shown in Table 6. The difference in the effect of DIM on DMI observed between the 2 equations may be attributable to differences in the covariates used to derive the equation. That is, whereas the proposed model included the effect of parity, BCS, and DIM expressed as days, the NRC model does not account for differences in parity or BCS and uses the week of lactation. Furthermore, the NRC prediction equation for DMI was derived mostly from data collected in the 1990s, such that it is likely that the data were collected on animals with lower milk production and lower body size than the animals in our database (mostly data collected since 2007). Additionally, since 2000, many researchers have focused on better understanding the mechanisms that control feed intake and the development of strategies to maximize intake during early lactation (
      • Drackley J.K.
      Biology of dairy cows during the transition period: The final frontier?.
      ;
      • Allen M.S.
      Effects of diet on short-term regulation of feed intake by lactating dairy cattle..
      ;
      • Grummer R.R.
      • Mashek D.G.
      • Hayirli A.
      Dry matter intake and energy balance in the transition period..
      ;
      • Allen M.S.
      • Bradford B.J.
      • Oba M.
      Board-invited review: The hepatic oxidation theory of the control of feed intake and its application to ruminants..
      ;
      • Allen M.S.
      Drives and limits to feed intake in ruminants..
      ). Although our proposed model accounts only for animal effects, we speculate that differences in management, environment, and diet composition (included as random effects in our model) could explain the difference observed between our proposed model and the NRC model with respect to the effect of DIM.
      The estimated partial regression coefficients of DMI on milk production and energy required for maintenance based on the proposed model (0.305 kg of DMI/Mcal of MilkE and 0.022 kg of DMI/kg of BW, respectively) differed from those of the
      NRC
      Nutrient Requirements of Dairy Cattle.
      model. For the NRC equation, DMI increased 0.372 kg of DMI/kg of FCM and 0.0968 kg of DMI/kg0.75 of BW, which convert to 0.330 kg of DMI/Mcal of MilkE and 0.014 kg of DMI/kg of BW for cows in the range of 500 to 700 kg of BW. Thus, in our proposed model, changes in DMI are less per unit change in MilkE output and more per unit change in BW than in the NRC model.
      Dry matter intake changes less per unit change in MilkE partly because the proposed equation has an intercept of approximately 1.7 kg/d and is adjusted for BCS, with low BCS increasing intake, especially for multiparous cows. Low BCS generally occurs when multiparous cows are at peak milk production, so this might decrease the coefficient for milk. In addition, our database includes only cows from the past approximately 10 yr, and they had greater milk production that those of the
      NRC
      Nutrient Requirements of Dairy Cattle.
      . These more recent cows also may have been fed diets that were higher in concentrates than the diets fed to the cows used in the NRC model so that less DMI was needed per unit of milk. Cows in our database also were of more modern genetics.
      • Moraes L.E.
      • Kebreab E.
      • Strathe A.B.
      • Dijkstra J.
      • France J.
      • Casper D.P.
      • Fadel J.G.
      Multivariate and univariate analysis of energy balance data from lactating dairy cows..
      showed that maintenance requirements have been increasing over the past 50 yr, and this may be part of the reason why the coefficient for BW is greater.
      A substantial difference between the proposed model and the
      NRC
      Nutrient Requirements of Dairy Cattle.
      DMI prediction equation is that our model includes the effect of BCS. The effect of BCS was highly significant (P < 0.01) and interacted with partiy, so that for each 1-unit decrease in BCS, intake was increased 0.7 kg/d for primiparous cows and 2.6 kg/d for multiparous cows. These results are consistent with the idea that thinner animals have greater appetites and consume more feed than fatter animals, at least partly because they have less circulating leptin (
      • Block S.S.
      • Smith J.M.
      • Ehrhardt R.A.
      • Diaz M.C.
      • Rhoads R.P.
      • Van Amburgh M.E.
      • Boisclair Y.R.
      Nutritional and developmental regulation of plasma leptin in dairy cattle..
      ;
      • Allen M.S.
      Drives and limits to feed intake in ruminants..
      ). The reason why BCS had less effect on DMI in primiparous cows may be that the cows are still growing and have flatter lactation curves than multiparous cows. We derived Equation 7 to illustrate how the proposed model (Equation 6) would be applied to a group of cows in which the proportion of primiparous and multiparous animals is known. Figure 6 illustrates the behavior of the proposed and
      NRC
      Nutrient Requirements of Dairy Cattle.
      models over a lactation curve of a dairy cow using the observed DMI, MilkE, and BW of the average cow in the evaluation data set.
      Figure thumbnail gr6
      Figure 6Lactation curve for the average Holstein dairy cow in the data set. The solid lines represent the observed milk energy (MilkE; blue with circles), observed DMI (red with diamonds), and observed BW (green with ×), and the dashed lines represent the DMI predicted by the
      NRC
      Nutrient Requirements of Dairy Cattle.
      (dotted line) and proposed (dashed line) models.
      Finally, we found important variability across diets in these partial coefficients. This is not surprising based on
      • Jensen L.M.
      • Nielsen N.I.
      • Nadeau E.
      • Markussen B.
      • Nørgaard P.
      Evaluation of five models predicting feed intake by dairy cows fed total mixed rations..
      and suggests that dietary characteristics can enhance the prediction of DMI when they are known. Dietary characteristics were not known for all cows in the current database, but the combination of diet and animal factors is addressed in a companion paper (
      • Allen M.S.
      • Sousa D.O.
      • VandeHaar M.J.
      Equation to predict feed intake response by lactating cows to factors related to the filling effect of rations..
      ).

      CONCLUSIONS

      The proposed model to predict DMI contains the effect of MilkE, BW, BCS, parity, and DIM. This model differs from the
      NRC
      Nutrient Requirements of Dairy Cattle.
      DMI prediction equation in that the NRC model contains only FCM, metabolic BW, and week of lactation. Using an independent evaluation data set, both models were similar for predicting DMI during early lactation (1–75 DIM), and the proposed model outperformed the NRC model during mid and late lactation (76–368 DIM). Finally, in regard to the estimated coefficients for milk and BW, the proposed model suggests less change in DMI for changes in milk production but greater changes in DMI for changes in cow BW. We recongnize that our equation is based on cows mostly from the eastern United States, but we suggest that our new equation can serve as a base prediction of voluntary DMI for lactating cows, especially when diet composition is unknown, and that feed factors can be used to adjust the predicted DMI during ration formulation, as described in
      • Allen M.S.
      • Sousa D.O.
      • VandeHaar M.J.
      Equation to predict feed intake response by lactating cows to factors related to the filling effect of rations..
      .

      ACKNOWLEDGMENTS

      This project was supported by Agriculture and Food Research Initiative Competitive Grant no. 2011-68004- 30340 from the USDA National Institute of Food and Agriculture (Washington, DC) and by Michigan State University AgBioResearch (East Lansing, MI). Rodrigo Araujo de Souza was supported by a PhD fellowship from National Council for Scientific and Technological Development (CNPq) from the Brazilian Ministry of Education (Brasilia, DF, Brazil). The authors thank L. Armentano, K. Weigel, D. Spurlock, C. Staples, E. Connor, H. Dann, W. Weiss, and M. Hanigan for providing some of the data used in this analysis.

      Supplementary Material

      REFERENCES

        • Allen M.S.
        Effects of diet on short-term regulation of feed intake by lactating dairy cattle..
        J. Dairy Sci. 2000; 83 (10908065): 1598-1624
        • Allen M.S.
        Drives and limits to feed intake in ruminants..
        Anim. Prod. Sci. 2014; 54: 1513-1524
        • Allen M.S.
        • Bradford B.J.
        • Oba M.
        Board-invited review: The hepatic oxidation theory of the control of feed intake and its application to ruminants..
        J. Anim. Sci. 2009; 87 (19648500): 3317-3334https://doi.org/10.2527/jas.2009-1779
        • Allen M.S.
        • Sousa D.O.
        • VandeHaar M.J.
        Equation to predict feed intake response by lactating cows to factors related to the filling effect of rations..
        J. Dairy Sci. 2019; 102: 7961-7969https://doi.org/10.3168/jds.2018-16166
        • Arnedal S.
        Predictions for voluntary dry matter intake in dairy cows.
        Department of Animal Nutrition and Management, Swedish University of Agricultural Sciences, Uppsala, Sweden2005 (PhD Diss)
        • Block S.S.
        • Smith J.M.
        • Ehrhardt R.A.
        • Diaz M.C.
        • Rhoads R.P.
        • Van Amburgh M.E.
        • Boisclair Y.R.
        Nutritional and developmental regulation of plasma leptin in dairy cattle..
        J. Dairy Sci. 2003; 86 (14594240): 3206-3214
        • Chan P.S.
        • West J.W.
        • Bernard J.K.
        Effect of prepartum dietary calcium on intake and serum and urinary mineral concentrations of cows..
        J. Dairy Sci. 2006; 89 (16428639): 704-713
        • Dann H.M.
        • Varga G.A.
        • Putnam D.E.
        Improving energy supply to late gestation and early postpartum dairy cows..
        J. Dairy Sci. 1999; 82 (10480103): 1765-1778
        • de Souza R.A.
        • VandeHaar M.J.
        Determining the change in body weight per unit of body condition score in Holstein cows..
        J. Dairy Sci. 2018; 101 (Abstr.): 153
        • Drackley J.K.
        Biology of dairy cows during the transition period: The final frontier?.
        J. Dairy Sci. 1999; 82 (10575597): 2259-2273
        • Fox D.G.
        • Van Amburgh M.E.
        • Tylutki T.P.
        Predicting requirements for growth, maturity, and body reserves in dairy cattle..
        J. Dairy Sci. 1999; 82 (10509256): 1968-1977
        • Grummer R.R.
        • Mashek D.G.
        • Hayirli A.
        Dry matter intake and energy balance in the transition period..
        Vet. Clin. North Am. Food Anim. Pract. 2004; 20: 447-470https://doi.org/10.1016/j.cvfa.2004.06.013
        • Huhtanen P.
        • Rinne M.
        • Mäntysaari P.
        • Nousiainen J.
        Integration of the effects of animal and dietary factors on total dry matter intake of dairy cows fed silage-based diets..
        Animal. 2011; 5 (22439992): 691-702https://doi.org/10.1017/S1751731110002363
        • Janovick N.A.
        • Boisclair Y.R.
        • Drackley L.K.
        Prepartum dietary energy intake affects metabolism and health during the periparturient period in primiparous and multiparous Holstein cows..
        J. Dairy Sci. 2011; 94 (21338804): 1385-1400
        • Jensen L.M.
        • Nielsen N.I.
        • Nadeau E.
        • Markussen B.
        • Nørgaard P.
        Evaluation of five models predicting feed intake by dairy cows fed total mixed rations..
        Livest. Sci. 2015; 176: 91-103
        • Krizsan S.J.
        • Sairanen A.
        • Höjer A.
        • Huhtanen P.
        Evaluation of different feed intake models for dairy cows..
        J. Dairy Sci. 2014; 97 (24508436): 2387-2397
        • Lu Y.
        • Vandehaar M.J.
        • Spurlock D.M.
        • Weigel K.A.
        • Armentano L.E.
        • Staples C.R.
        • Connor E.E.
        • Wang Z.
        • Coffey M.
        • Veerkamp R.F.
        • de Haas Y.
        • Tempelman R.J.
        Modeling genetic and nongenetic variation of feed efficiency and its partial relationships between component traits as a function of management and environmental factors..
        J. Dairy Sci. 2017; 100 (27865511): 412-427
        • Moraes L.E.
        • Kebreab E.
        • Strathe A.B.
        • Dijkstra J.
        • France J.
        • Casper D.P.
        • Fadel J.G.
        Multivariate and univariate analysis of energy balance data from lactating dairy cows..
        J. Dairy Sci. 2015; 98 (25892698): 4012-4029https://doi.org/10.3168/jds.2014-8995
        • NRC
        Nutrient Requirements of Dairy Cattle.
        7th rev. Natl. Acad. Press, Washington, DC2001 (ed)
        • Rayburn E.B.
        • Fox D.G.
        Variation in neutral detergent fiber intake of Holstein cows..
        J. Dairy Sci. 1993; 76: 544-554
        • Roseler D.K.
        • Fox D.G.
        • Chase L.E.
        • Pell A.N.
        • Stone W.C.
        Development and evaluation of equations for prediction of feed intake for lactating Holstein dairy cows..
        J. Dairy Sci. 1997; 80 (9178128): 878-893
        • Tedeschi L.O.
        Assessment of the adequacy of mathematical models..
        Agric. Syst. 2006; 89: 225-247
        • Tempelman R.J.
        • Spurlock D.M.
        • Coffey M.
        • Veerkamp R.F.
        • Armentano L.E.
        • Weigel K.A.
        • de Haas Y.
        • Stamples C.R.
        • Connor E.E.
        • Lu Y.
        • VandeHaar M.J.
        Heterogeneity in genetic and nongenetic variation and energy sink relationships for residual feed intake across research stations and countries..
        J. Dairy Sci. 2015; 98 (25582589): 2013-2026
        • VandeHaar M.J.
        • Armentano L.E.
        • Weigel K.
        • Spurlock D.M.
        • Tempelman R.J.
        • Veerkamp R.
        Harnessing the genetics of the modern dairy cow to continue improvements in feed efficiency..
        J. Dairy Sci. 2016; 99 (27085407): 4941-4954
        • Zom R.L.G.
        • Andre G.
        • van Vuuren A.M.
        Development of a model for the prediction of feed intake by dairy cows: 1. Prediction of feed intake..
        Livest. Sci. 2012; 143: 43-45https://doi.org/10.1016/j.livsci.2011.08.014

      Linked Article

      • Equation to predict feed intake response by lactating cows to factors related to the filling effect of rations
        Journal of Dairy ScienceVol. 102Issue 9
        • Preview
          Our objective was to predict the dry matter intake (DMI) response during ration formulation to factors related to the filling effect of rations and their interaction with milk yield (MY) by lactating cows past peak lactation. A data set was developed consisting of 134 treatment means from 34 experiments reported in 32 peer-reviewed articles published from 1990 through 2015. The data set included data for cows ranging from 60 to 309 d postpartum with mean DMI ranging from 17.6 to 30.6 kg/d and MY ranging from 20.3 to 51.1 kg/d.
        • Full-Text
        • PDF
        Open Archive