HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS		QUICK SEARCH:   [advanced] Author: Keyword(s): Year:  Vol:  Page:

This Article Abstract Full Text (PDF) Alert me when this article is cited Alert me if a correction is posted Services Similar articles in this journal Similar articles in PubMed Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via Google Scholar Google Scholar Articles by Martin, O. Articles by Sauvant, D. Search for Related Content PubMed PubMed Citation Articles by Martin, O. Articles by Sauvant, D.

Metaanalysis of Input/Output Kinetics in Lactating Dairy Cows

UMR INRA INA P-G, Physiologie de la Nutrition et Alimentation, 16 rue Claude Bernard, 75231 Paris Cedex 05, France

	ABSTRACT

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION CONCLUSIONS ACKNOWLEDGEMENTS REFERENCES

 
The kinetics of dry matter intake, raw milk yield, milk fat yield, and body weight of 140 treatment groups of cows reported in 37 articles published between 1959 and 1999 were collected in a database. The majority of these experiments were performed with Holstein-Friesian cows (primiparous: n = 51, multiparous: n = 89). Animals were fed with forages and concentrates, offered separately (n = 98) or as a total mixed ration (n = 42). Empty body weight was calculated assuming that an increase of 1 kg in dry matter intake corresponded to an increase of 4 kg in gastrointestinal tract content. Each kinetic, denoted Y, was fitted with the model proposed by Grossman (1999) for raw milk yield. This model allowed the partitioning of the kinetic into three stages (early, middle, and late stages of lactation), each characterized by its own linear tendency. Kinetics were therefore summarized with parameter estimates of the model: the plateau value of Y in middle stage, the slopes of Y in early and late stages, and the weeks of start and end of the middle stage. The following additional variables were also calculated: the total increase of Y in early stage, the extrapolated value of Y at calving, and the persistency of Y defined as the duration of the middle stage. Variability of these parameters for each kinetic across treatment groups was examined through mixed models including the year of publication as covariate, parity as a fixed between-group factor, and trial as a random factor. Correlations between parameters within a Y-kinetic and between different Y-kinetics were calculated by parity group and feed allocation system. Finally we propose predictive models which emphasized (1) the effect of empty body weight at calving on the depletion of empty body weight in early lactation, (2) the incidence of the kinetic of dry matter intake on the kinetic of raw milk yield in middle and late stages of lactation, and (3) the relationships between the kinetics of raw milk yield and milk fat yield in middle and late stages of lactation.

Abbreviation key: CO = concentrate, , EBW = empty body weight, , FO = forage, , MEI = metabolizable energy intake, , MFC = milk fat content, , MFY = milk fat yield, , MLC = milk lactose content, , MLY = milk lactose yield, , MPC = milk protein content, , MPY = milk protein yield, , Multi = multiparous cow, , Primi = primiparous cow, , RMSE = root mean square error, , RMY = raw milk yield, , YOP = year of publication

Key Words: metaanalysis • milk yield • dry matter intake • body weight

	INTRODUCTION

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION CONCLUSIONS ACKNOWLEDGEMENTS REFERENCES

 
The majority of articles published on dairy nutrition and feeding have dealt with short periods of experimentation or with designs where animals receive successive experimental diets (e.g., Latin square designs). These trials provide static data, which allows empirical models of the effect of a particular type of diet or allocation system on animal production responses to be calculated for a given period of lactation. However, lactation is not a set of successive independent periods; it needs to be considered as a dynamic process (Broster and Broster, 1984). This is of practical importance. Knowledge regarding the impact of feeding on the kinetics of lactation is necessary for the optimization of a feeding strategy subject to constraints of feed inventory and fluctuating input and output prices. It therefore appeared important to examine the kinetic patterns that characterize the course of the lactation (DMI, yield, BW) and to model their interactions. In particular, the purposes of this study were: 1) to evaluate the variability of the kinetic patterns, their evolution through years of experimentation, and their differences according to parity; and 2) to highlight the dependencies within a given kinetic or between different kinetics, which determine a particular profile.

To achieve this goal, we decided to build a database by pooling the results from published trials dealing with the dynamic features of feeding, lactation, and BW. This paper presents the major characteristics of this collection of curves and the main effects revealed by a first metaanalysis approach.

	MATERIALS AND METHODS

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION CONCLUSIONS ACKNOWLEDGEMENTS REFERENCES

 
Compilation of the Literature
Articles and experimental designs.
Experimental factors and results published in 37 articles covering input/output kinetics of lactating dairy cows were collected in the database. The articles, published between 1959 and 1999, reported trials performed in North America (n = 17), Europe (n = 17), Asia (n = 1), Africa (n = 1), and Oceania (n = 1; Table 1). The criteria used for inclusion were the existence of data concerning the kinetics of DMI (kg/d) and/or metabolizable energy intake (MEI, Mcal/d), raw milk yield (RMY, kg/d) and/or 4% FCM yield (kg/d), and BW (kg), all of which had to have been recorded for at least the first 8 wk of lactation. When available, the kinetics of milk fat yield (MFY, g/d) or content (MFC, g/kg) [MFC = MFY/RMY], milk protein yield (MPY, g/d) or content (MPC, g/kg) [MPC = MPY/RMY], and milk lactose yield (MLY, g/d) or content (MLC, g/kg) [MLC = MLY/RMY] were recorded (Table 1). The data sets used were weekly (32 articles) or monthly means (five articles) of daily records. Most of the data were published as graphs: consequently, they were stored in the database after digitalization. Each graph was scanned (IrfanView 3.36, Irfan Skiljan, http://www.irfanview.com/), and the picture file was incorporated as the background of an indexed graphical interface. Relative positions and scales of the picture axis were used to calibrate (X,Y) coordinates, and data points were translated into numerical values. These calculations were performed with Excel 2000 (Microsoft Corporation). This digitalization procedure was tested on simulated graphs. This method was found to be unbiased and accurate with less than 1% error.

Animals and housing.
The experiments were mainly conducted with Holstein-Friesian cows. However, the two trials reported by Wiktorsson (1971, 1973) were performed with cows of the Swedish Red and White breed. Ayrshire cows (pure- or cross-bred) comprised half of the livestock used in the experiment of Garnsworthy (1982). Moreover, Ayrshire and Guernsey breeds were involved in the trial of Everson (1976), and Jersey cows in the trial of Swanson (1962). Within a treatment group, animals were predominantly of the same parity: 51 treatments used primiparous cows (Primi) and 89 multiparous cows (Multi; including 10 heterogeneous treatment groups with between 25 to 34% Primi which were treated as Multi). Animals were mostly housed individually in tie stalls or stanchion barns.

Diets and feeding.
Cows from each treatment group were individually fed during one or several experimental subperiods according to a given nutritional plan. The dietary treatments were designed using combinations of the type of diet and its allocation system. The feeding regimens used three basic types of feeds: forage (FO), concentrate (CO), and TMR. The feed allocation systems were applied to either FO, CO or TMR, and there were four major types: 1) Restricted: allowance according to a predetermined flat rate or gradual increase up to a ceiling during the experimental period; 2) Ad libitum: feed allowance to appetite (seven papers reported a standard daily ort of 5 to 12% on a DM basis allowed); 3) Fitted: feed allowance adjusted at given intervals to current yield of milk (for production requirements) and/or BW (for maintenance requirements)—the rate of feeding per unit output (kg FCM, kg RMY, kg BW, or kg BW^0.75) was constant throughout the whole experimental period or gradually changed; and 4) Mixed: allowance of part of the ration in a restricted manner, the remainder being fed in a fitted manner. As a codifying rule, if a treatment group was fed more than 50% of the diet (in kg DM) ad libitum, fitted, or restricted, it was labeled as such. Since the fitted system was applied on the basis of expected nutrient requirements, the rate of feeding per unit output determined the allocation scheme. Thus, depending on the rate, the fitted system could be close to the ad libitum system (liberal allocation) or close to the restricted system (limited allocation). A total of 111 different feedstuffs were used in these experiments: 29 TMR, 38 FO, and 44 CO. The FO were either hays (mainly meadow grass or alfalfa), grass or corn silages. The CO were mixtures of cereals (mainly barley, corn or oats), cereal byproducts (mainly wheat bran, brewers grains or corn gluten meal), oil meals (mainly soybean meal, peanut meal, rapeseed meal, coconut meal or cotton meal), and agro-industrial byproducts (mainly sugar beet pulp or molasses). The TMR were mixtures of roughages (mainly corn silage, meadow or alfalfa hay, or grass silage) and energy feeds (commercial dairy compound or mixtures of barley, corn, sugar beet pulp, brewers grains, and oil meals). All the available data concerning the chemical composition of feeds were systematically recorded in the database.

Analysis and records.
In the majority of cases, feed allowance and intake were weighed daily. In most of the trials, samples of feed were collected weekly for pooled monthly chemical analyses, which were carried out according to conventional AOAC methods. Cows were generally machine-milked twice daily, except in Pinchasov (1982) and three treatment groups in DePeters (1985) for which cows were milked three times a day. Analysis of milk composition was generally performed monthly on weekly-pooled frozen samples. MFC was determined using the Babcock method (in the United States), the Gerber method (in Europe), or by infrared analysis (in the 1990’s). MPC was assayed using Kjeldahl determination (N x 6.38), orange G dye binding method (Udy, 1956), amido-black colorimetry (Posthumus, 1960) or infrared analysis. The MLC content was determined according to Somogyi (1952), Rook and Line (1961), or infrared analysis. The energy values of feeds were calculated or estimated by the authors according to the method generally used in the corresponding country. Animals were usually weighed twice per month on 2 consecutive d at the same hour.

Data Analysis
Additional variables.
The percentages of FO and CO in diets allowed the calculation of the DMI of FO (kg/d), and CO (kg/d). DMI of the TMR was also noted (kg/d). The BW was used as such and also corrected for expected changes in gastrointestinal fill related to DMI, where a mean digestive tract content correction factor of 4 kg/kg DMI was used, assuming an increase in the weight of the reticulo-rumen content of 3 kg/kg DMI (Rémond, 1988) and an additional increase in intestinal content of 1 kg/kg DMI (Chilliard et al., 1987) [Empty BW (EBW, kg) = BW (kg) – 4 x DMI (kg)]. This correction was made to reduce the variability in the evolution of BW, especially in early lactation when DMI increased, and to estimate the real loss and gain of body reserves.

Treatment group general outcomes.
Letting Y be the response variable (DMI, FO DMI, CO DMI, TMR DMI, MEI, RMY, FCM, MFY, MFC, MPY, MPC, MLY, MLC, BW, or EBW), the mean values of Y between wk 1 and 8 of lactation (Y[1-8]) were calculated for each treatment group. This measure was regarded as an overall level indicator of the Y-kinetic, from a static standpoint.

Two sequential mixed models were fitted with Proc MIXED of SAS (1992) to determine the effects of parity and year of publication (YOP) on Y[1-8]. The first model focuses on the effect of parity on the slope of the regression:

[1]

where,

P_j=fixed effect of parity group j (j = 1, 2+),

t_k=random effect of trial k (various number of trial according to Y[1-8]),

B_j=regression coefficient of Y[1-8] on YOP for parity group j,

e_ijk= error term N(0, ²_e).

A significant effect of P_j indicated that the intercepts of the slope of the linear regression of Y[1-8] on YOP were different from zero for at least one parity group. Therefore, a second model was fitted as follows:

[2]

where symbols are as previously defined with the following changes:

µ=overall mean,

ß₀=overall regression coefficient of Y[1-8] on YOP,

ß_j=deviation of regression coefficient of Y[1-8] on YOP from the overall regression coefficient for parity group j.

A significant ß_j effect indicated that the slope of Y[1-8] on YOP was different for parity groups. Because of the unequal number of lactation (cow-year) between treatment groups and since the standard errors of Y by week were generally not known, least-squares were weighted with where N_i was the number of lactations of treatment group i.

This first statistical model was aimed at giving an overview of the database without respect of dynamic aspects of lactation.

Modeling kinetics.
Among the group of curves associated with lactation in the dairy cow, the kinetics of RMY have been most often described and modeled with time dependent functions (see reviews of Rowlands et al., 1982; Masselin et al., 1987; Papajcsik and Bodero, 1988; Rook et al., 1993). Among them, the gamma function y(t) = a•t^b•e^–ct proposed by Wood (1967, 1968, 1969, 1977) has been extensively used to summarize the pattern of lactational performance. In this function, y(t) is the average daily milk yield in wk t of lactation, and a, b, and c are positive parameters that determine the scale and shape of the curve. While improving the pioneering work of Wood, several authors have developed alternative forms of this model (Kuck et al, 1976; Cobby and Le Du, 1978; Dhanoa, 1981; Schneeberger, 1981). More recently, Goodall and Sprevak extended the fitting properties of Wood’s model with stochastic (1984) and recursive (1985) procedures. Since the 1980’s, numerous authors have suggested other mathematical approaches for fitting the milk yield x time data, either empirical (Emmans et al., 1983; Grossman and Koops, 1988; Morant and Gnanasakthy, 1989; Rook et al., 1993; Grossman et al., 1999), mechanistic (Neal and Thornley, 1983; Dijkstra et al., 1997), nonparametric (Elston et al., 1989) or autoregressive (Dhanoa and Le Du, 1982; Deluyker et al., 1990). Much less research has been undertaken on specific mathematical models to fit other curves characterizing the lactation profile. Nevertheless, several authors have proposed an application of their milk model to milk constituent yields (Schneeberger, 1979; Morant and Gnanasakthy, 1989) or contents (Wood, 1976), or to BW changes (Wood, 1979; Emmans et al., 1983). Stanton et al. (1992) proposed a statistical mixed model for milk, fat, and protein lactation curves. When fitting intake kinetics, Ostergaard (1979) proposed, in the context of ad libitum intake, to use the function y(t) = a – b•e^–ct – d•t, where y(t) is the average daily DMI in wk t, and a, b, and c are positive parameters which determine the scale and shape of the curve. In the present paper, we wanted to use 1) the same model to fit the kinetics of all variables Y, and 2) a model which parameters quantified a specific and concrete aspect of the global pattern of the curve. This model was aimed at summarizing all curves with the same analytical rule in equivalent vectors of parameters. Observed kinetics of Y over the recorded experimental period [t0; t3] (t₀ 1 and t₃ 45) of each treatment group were then summarized with the following triphasic model (only works concerning DMI, RMY, EBW, and MFY are considered in the present paper).

[3]

where,

Y*=predicted value of DMI, RMY, MFY, and EBW,

t=time, in weeks of lactation from calving,

[t₀; t₁]=early stage of lactation,

[t₁; t₂]=middle stage of lactation,

[t₂; t₃]=late stage of lactation,

k_E=slope of Y as a function of time in early stage of lactation (k_E 0 for Y = DMI, RMY, and MFY; k_E 0 for Y = EBW),

k_L=slope of Y as a function of time in late stage of lactation (k_L 0 for Y = DMI, RMY, and MFY; k_E 0 for Y = EBW),

Y_C=extrapolated value of Y for t = 0, calving time,

Y_M=plateau value of Y in middle stage lactation.

Using the above model, each stage was summarized by its linear trend: increasing or decreasing for early and late stages, constant for the middle stage. In order to build a smooth function from the linear equations system [3], we used the smoothing logistic transition between intersecting straight lines proposed by Koops and Grossman (1993) and employed by Grossman (1999) in the case of Y = RMY. Thus, the model finally used for each Y variable was (Figure 1):

[4]

View larger version (12K): [in this window] [in a new window]   Figure 1. Model proposed by Grossman et al.(1999) (smoothing logistic transition between intersecting straight lines): 25 where YM = plateau value of Y in middle stage, kE = slope of Y in early stage, t1 = upper bound of early stage, r1 and r2 = durations of smoothing transition between linear phases set to 0.15 wk, kL = slope of Y in late stage, and t2 = upper bound of middle stage. Additional parameters: YC = extrapolated value of Y at calving time (t = 0), YE = total variation in Y during early stage, YP = duration of middle stage.

[5a]

[5b]

where cov(k_E, t₁) is the asymptotic covariance of parameter estimates k_E and t₁. The extrapolated value of Y at calving (Y_C, in unit of Y) and its standard error s_C were estimated with the following formulas:

[6a]

[6b]

where cov(Y_M, t₁) and cov(Y_M, k_E) are the asymptotic covariances of parameter estimates Y_M and t₁, and Y_M and k_E respectively. Grossman (1999) proposed to define the persistency of RMY as the duration of the middle stage. We employed this terminology for all kinetics. The persistency of Y (Y_P, in wk) and its standard error s_P were estimated with the following formulas:

[7a]

[7b]

where cov(t₁, t₂) is the asymptotic covariance of parameter estimates t₁ and t₂. Estimates of parameters were used for further analysis if the following conditions were respected:

View larger version (29K): [in this window] [in a new window]   Figure 2a. Effect of increasing values of parameter r1 on the shape of the model proposed by Grossman et al. (1999): 26 (Model parameters: YM = 30.0, kE = 2.0, kL = –0.5, t1 = 5, t2 = 15, r1 = {0.15, 1, 2, 3, ..., 9}, r2 = 0.15). Figure 2b. Effect of increasing values of parameter r2 on the shape of the model proposed by Grossman et al. (1999): 27 (Model parameters: YM = 30.0, kE = 2.0, kL = –0.5, t1 = 5, t2 = 15, r1 = 0.15, r2 = {0.15, 1, 2, 3, ..., 9}).

t₂:t₂ t₀ and t₂ t₃ (middle stage detected)

k_E:t₁ t₀ (early stage detected)

k_L:t₃ – t₂ > 3 (slope estimated over at least 3 wk)

Y_M:t₃ > t₁ (end of early stage detected)

Y_E and Y_C:Y_M > Y_C (valid estimate of a positive value of Y_E)

Y_P:t2 t3 (non null persistency)

Finally a database was created with, for each Y-kinetic of each treatment group, the vector of valid estimates of eight parameters Y_C, Y_E, Y_M, k_E, k_L, t₁, t₂, Y_P and their standard errors s_C, s_E, s_M, s_kE, s_kL, s₁, s₂, and s_P.

This modeling approach allows one to decompose and to quantify the pattern of each Y-kinetic by the use of two slopes (k_E and k_L), two times (t₁ and t₂), one duration (Y_P), two levels (Y_M and Y_C), and one change (Y_E).

Two sequential mixed models were fitted with Proc MIXED of SAS (1992) to determine the effects of parity and YOP on parameter estimates. The first model focuses on the effect of parity on the slope of regression:

[8]

where,

(_p)_ijk=estimate of the parameter _p = {Y_C, Y_E, Y_M, k_E, k_L, t₁, t₂, Y_P} for the Y-kinetic of treatment group i in trial k (various number of treatment groups and trials according to Y-kinetic and parameter p) for parity group j (j = 1, 2+),

P_j=fixed effect of parity group j,

t_k=random effect of trial k,

B_j=regression coefficient of _p on YOP for parity group j,

e_ijk=error term N(0, 2_e).

As in model [1], a significant effect of Bj indicated that the slope of the linear regression of _P on YOP was different from zero for at least one parity group. Therefore, a second model was fitted as follows:

[9]

where symbols are as previously defined with the following changes:

µ=overall mean,

ß₀=overall regression coefficient of _p on YOP,

ß_j=deviation of regression coefficient of p on YOP from the overall regression coefficient for parity group j.

As in model [2], a significant ßj effect indicated that the slope of p on YOP was different for parity groups.

To account for unequal accuracy of estimates, least-squares means were weighted as follow: for parameter _p with standard error s_p, w₁ is defined as the inverse of the squared of s_p, and as its mean value over available data. Thus, w₂ = w₁/ was used as weight (St-Pierre, 2001). This statistical model was aimed at quantifying the variability of the pattern of Y-kinetics across YOP and parity.

Three levels of correlations between parameter estimates were examined for kinetics of DMI, EBW, RMY and MFY:

1. Between parameter estimates of the same Y-kinetic within treatment group (e.g. DMI_kE and DMI_kL for the first treatment group of Johnson (1983)). These correlations straight out of NLIN were implied by the mathematical properties of model [4].
   
2. Between parameter estimates of the same Y-kinetic across treatment groups (e.g., EBW_C and EBW_E). These correlations provided information about the structural dependencies within the pattern of each Y-kinetic.
   
3. Between parameter estimates between Y-kinetics across treatment groups (parameters can be equivalent, e.g. DMI_M and RMY_M, or different, e.g. DMI_kE and RMY_kL). These correlations provided information about the relationships between the patterns of different Y-kinetics.
   

Asymptotic correlations of level (1) are model-dependent and therefore not different between Y-kinetics. The sign of parameters k_E (negative for EBW, positive for DMI, RMY, and MFY) and k_L (positive for EBW, negative for DMI, RMY, and MFY) determine the sign of all correlations except that between t₁ and t₂. This latter correlation is always negative and induced by the structural constraint t₂ t₁. Non null correlations exist within stages (within linear segments): between k_E and t₁ in early stage, between t₁ and Y_M and between Y_M and t₂ in middle stage, and between t₂ and k_L in late stage. Correlations between parameters of different stage were approximately null (in particular between k_E and Y_M, Y_M, and k_L, and k_E and k_L).

These results justify the corrections made with covariance components for unbiased estimates of parameter YE [5a] and standard error of parameters Y_E [5b], Y_C [6b], and Y_P [7b].

Correlations of levels (2) and (3) were all calculated within parity groups and within allocation systems A, F, and R. Only correlations with the same sign and an absolute value greater than 0.5 in each parity group and allocation system were retained. Correlations within Y-kinetic between parameters Y_P and Y_t2 and between parameters Y_kE and Y_E were removed since they were implied by construction ([6a] and [7a]). Finally the conditions were matched for two couples of parameters within the EBW-kinetic and four couples of parameters between kinetics of DMI, RMY and MFY. These six correlations by parity group and allocation system are presented in Table 6.

[10]

where,

_Y)_ijk=estimate of the parameter _Y = {Y_C, Y_E, Y_M, k_E, k_L, t₁, t₂, Y_P} for the Y-kinetic of treatment group i in trial k (various number of treatment groups and trial according to couple of parameters (_Y, _X) and Y-kinetic) for parity group j (j = 1, 2+),

P_j=fixed effect of parity group j,

t_k=random effect of trial k;

B_j=regression coefficient of _Y on _X for parity group j;

(_X)_ijk=estimate of the parameter _Y = {Y_C, Y_E, Y_M, k_E, k_L, t₁, t₂, Y_P} for the Y-kinetic of treatment group i for parity group j,

e_ijk=error term N(0, ²_e).

As in model [1] and [5], a significant effect of Bj indicated that the slope of the linear regression of _Y on _X was different from zero for at least one parity group. Therefore, a second model was fitted as follows:

[11]

where symbols are as previously defined with the following changes:

µ=overall mean,

ß₀=overall regression coefficient of _Y on _X;

ß_j=deviation of regression coefficient of _Y on _X from the overall regression coefficient for parity group j.

As in model [2] and [6], a significant ßj effect indicated that the slope of _Y on _X was different for parity groups. If not significant, a third model was fitted as follows:

[12]

where symbols are as previously defined.

A significant P_j effect indicated that the intercept of the regression of _Y on _X was different for parity groups. If not significant, a fourth model was fitted as follows:

[13]

where symbols are as previously defined.

To account for unequal accuracies of estimates, least-squares were weighted as follow: for parameters _X and _Y with respective standard errors s_X and s_Y, w₁ is defined as the inverse of the product s_X · s_Y, and as its mean value over available data. Thus, w₂ = w₁/ was used as weight (St-Pierre, 2001). By this way, both error on _X and _Y were used to estimate the regression coefficients.

These four models were aimed at quantifying the relationships between:

1. the loss of EBW in early lactation and the EBW at calving,
   
2. the plateau values of RMY and DMI in middle stage of lactation,
   
3. the slopes of decrease of RMY and DMI in late stage of lactation,
   
4. the plateau values of MFY and RMY in middle stage of lactation,
   
5. the slopes of decrease of MFY and RMY in late stage of lactation.
   

The choice of independent (Y) vs dependent variable (X) was set in the sense of prediction of performances. Dependent variables were thus, either state at calving (EBW_C) for prediction of EBW loss, characteristics of intake pattern (DMI_M and DMI_kL) for prediction of equivalent characteristics of RMY, or characteristics of milk yield pattern (RMY_M and RMY_kL) for prediction of equivalent characteristics of MFY pattern.

In statistical models [1], [2], [8], [9], [10], [11], [12], and [13], trials were considered as blocks and their effect was incorporated as random intercept. When significant, random slopes (interaction of trial by dependent variable) were incorporated in the reduced models [10], [11], [12], and [13]. This methodology takes into account the fact that observations within a given trial have more in common that observations across trials (St-Pierre, 2001).

	RESULTS

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION CONCLUSIONS ACKNOWLEDGEMENTS REFERENCES

 
Database
Diets.
Table 2 describes the average chemical composition of diets fed to the 140 treatment groups. The availability of data was variable amongst articles. The major sources of variability appeared to be associated with DM, ash, phosphorus, ether extract, and crude fiber content (Coefficient of variation greater than 20%). Metabolizable energy and CP contents, ranging from 1.5 to 3.0 Mcal/kg DM and from 10.1 to 20.5% reflect experimental designs investigated.

Feeding systems.
A summary of feeding systems used across the 140 treatment groups is given in Table 3. From 1970 to 1999, TMR were only fed ad libitum. Restricted allowances of FO or CO were used before 1984. Fitted or ad libitum programs of feeding FO or CO were investigated from 1965 to 1999. Figure 3 shows the average fractions of diet fed A, F, and R according to the decade the article was published. The major part of DMI was mainly fed in a restricted way before 1970 (> 50% DMI) and was primarily fed ad libitum after 1980 (from 53% to 100% DMI amongst treatment groups). Before 1980, the underlying research was aimed at measuring the quantitative effects of intake on lactational performance (Table 1). Since 1980, ad libitum systems have been used predominantly, reflecting practical on-farm conditions: qualitative effects, mainly that of feed type, were therefore studied. Thus, data gathered across treatment groups were highly unbalanced with respect to feeding system.

View larger version (16K): [in this window] [in a new window]   Figure 3. Number of treatment groups (histogram) and average fraction of diet (%DMI) fed, ad libitum (), fitted to current performances () or restricted (•) according to decade of article publication.

View larger version (11K): [in this window] [in a new window]   Figure 4. Box plot of R2 of model of Grossman et al. (1999) fitted to treatment group kinetics of DMI, raw milk yield (RMY), milk fat yield (MFY), and empty BW (EBW). Bottom and top edges of the box are located at the sample 25th and 75th percentiles; center horizontal line is at the sample median; central filled circle is at the sample mean; asterisks are values more extreme than three interquartile ranges (distance between the 25th and 75th sample percentiles).

Weighted least squares means of parameters were different between parity group except for t₁ (6.3 ± 0.5 wk for pooled parity groups) and DMI_E (6.18 ± 0.48 kg for pooled parity groups). Primi had a lower level of DMI at calving than Multi (DMI_C: –2.44 kg) and a lower slope of increase in early stage (Primi: k_E = 0.94 ± 0.16 kg/wk, Multi: k_E = 1.16 ± 0.15 kg/wk) which induced a lower level in middle stage (Primi: DMI_M = 15.69 ± 0.33 kg, Multi: DMI_M = 18.63 ± 0.41 kg). It should be mentioned that the range of k_E estimates was higher for restricted schemes of feeding than for others. The end of the plateau stage was reached significantly later by Primi (t₂ = 26.4 ± 1.7 wk, DMI_P = 19.8 ± 1.5 wk) than Multi (t₂ = 12.6 ± 1.4 wk, DMI_P = 6.9 ± 1.2 wk).

No correlation was retained within the kinetic of DMI. Only kinetics of the ad libitum allocation system could express animal requirements. Patterns of the fitted and restricted allocation systems are mainly determined by experimental designs.

Empty body weight.
Goodness-of-fit was not different between parity groups, and the overall average RMSE was 4.2 ± 1.8 kg. The distribution of the 128 R² calculated is given in Figure 4 (R² = 93.1 ± 7.4%).

From 1959 to 1999, the annual increases in EBW_M and EBW_C were significant (P < 0.01) and not different between parity groups (EBW_M: + 1.5 kg/yr and EBW_C: + 2.1 kg/yr for pooled parity groups). The effect of YOP was higher for EBW_C than for EBW_M because of the associated reduction in k_E while t₁ was not affected. Thus, Primi and Multi were heavier at calving and their total loss of EBW in early stage was greater (EBW_E: – 0.8 kg/yr for pooled parity groups). Time t₂ was also influenced by the YOP for both parity groups (t₂: –0.2 wk/yr for pooled parity groups). No significant effect of the YOP on parameters k_L and EBW_P was detected.

Weighted least squares means of parameters were different between parity group except for the slopes k_E and k_L (–9.4 ± 1.3 kg/wk and 2.4 ± 0.2 kg/wk respectively for pooled parity groups) and for the time t₂ (12.4 ± 1.3 wk for pooled parity groups). Multi had a higher level of EBW_C (+ 80.2 kg) and a higher level of EBW_M (+ 69.8 kg) induced by a greater EBW_E in early stage (– 14.2 kg) achieved 1.1 wk later than Primi (Primi: t₁ = 3.3 ± 0.3, Multi: t₁ = 4.4 ± 0.2). Duration of the middle stage (plateau at minimum EBW before regrowth) was not different between parity groups (EBW_P = 8.7 ± 1.0 wk).

The only correlations between parameters within Y-kinetics deemed coherent within parity and allocation system (Table 6) were found between EBW_E and EBW_C, and between EBW_M and EBW_C (respectively = –0.73 and = +0.89 for pooled data). These relationships stressed the primary influence of animal size at calving on the early decrease of BW.

The relationship between EBW_E and EBW_C is presented in Figure 5. Weighting factors applied were globally homogeneous across treatment groups. Slopes of EBW_E on EBW_C were different from zero within parity (P < 0.001, model [10]). The difference in the slopes between parity groups was not significant (P = 0.688, model [11]) and the difference in the intercepts between parity group was significant (P < 0.001, model [12]).

View larger version (15K): [in this window] [in a new window]   Figure 5. Plot of EBWE (total loss of empty BW in early lactation) on EBWC (empty BW at calving) for multiparous (•) and primiparous (). Data are linked by trial and sized with respect to weight used in weighted least-squares analysis.

[14]

where = 0.0 for Primi and = 19.6 for Multi. This model was fitted for 104 treatment groups spread across 34 trials, and had a residual error of 9.0 kg (residual standard error of random intercept: 13.0 kg). Regression coefficient of model [14] gives an estimate of the within-trial effect of animal weight at calving on early depletion.

Raw milk yield.
Goodness-of-fit was not different between parity groups and the overall average RMSE was 0.66 ± 0.34 kg/d. The distribution of the 119 R² _{calculated is given in Figure 4 (R}² _{= 96.1 ± 4.5%). These results were better than those obtained with the alternative form of the lactation model of Wood (Wood, 1967) proposed by Dhanoa (Dhanoa, 1981) using the same dataset (RMSE = 0.79 ± 0.36 kg/d, R}² _{= 93.8 ± 6.9%).}

From 1959 to 1999, the annual increases in RMY_M and RMY_C were significant (P < 0.05) and not different between parity groups (RMY_M: + 0.40 kg/yr and RMY_C: + 0.28 kg/yr for pooled parity groups). The effect of YOP was significant on t₂ and RMY_P (t₂: +0.1 wk/yr and DMI_P: + 0.2 wk/yr for pooled parity groups). No significant effect of the YOP was detected on other parameters.

Weighted least squares means of parameters were different between parity group except for k_E and t₁ (k_E = 2.52 ± 0.40 kg/wk and t₁ = 4.1 ± 0.2 wk for pooled parity groups). Multi had a higher level of RMY_C (+ 7.03 kg) and a higher level of RMY_M (+ 7.93 kg). Multi had a shorter persistency (Primi: RMY_P = 8.1 ± 0.7 wk, Multi: RMY_P = 4.1 ± 0.4 wk) obtained with an earlier time t₂ (Primi: t₂ = 7.5 ± 0.5 wk, Multi: t₂ = 9.0 ± 0.7 wk). The weekly rate of decrease in RMY in late stage was significantly higher for Primi (k_L = – 0.29 ± 0.03 kg/wk) than Multi (k_L = – 0.56 ± 0.02 kg/wk).

No correlation was retained within the kinetic of RMY. Nevertheless, the plateau value and the slope of decrease in RMY during late stage were found to be correlated with equivalent parameters for DMI. The relationship between RMY_M and DMI_M is presented in Figure 6. Since the accuracy of parameter estimates was not equivalent across treatment groups, the weighting factors applied were heterogeneous. The difference in the slopes of RMY_M on DMI_M between parity group was significant (P < 0.001, model [11]). The difference in intercepts between parity groups was significant (P < 0.01, model [11]). Random effects on slope, on intercept, and random interaction between random slope and intercept were found to be significant (respectively P = 0.068, P = 0.066, and P = 0.071) and were therefore incorporated in retained model [10].

View larger version (14K): [in this window] [in a new window]   Figure 6. Plot of RMYM (plateau value of raw milk yield in middle stage of lactation) on DMIM (plateau value of DMI in middle stage of lactation) for multiparous (•) and primiparous (). Data are linked by trial and sized with respect to weight used in weighted least squares-analysis.

[15]

and the equation for Primi is:

[16]

This model was fitted for 119 treatment groups spread across 31 trials and had a residual error of 0.32 kg (residual standard error of random intercept: 13.7 kg). Regression coefficients of models [15] and [16] give an estimate of the within-trial increase in plateau value of milk yield associated with an increase in plateau value of DMI (about 1.5 and 0.4 kg RMY/kg DMI for Multi and Primi respectively).

The relationship between RMY_kL and DMI_kL is presented in Figure 7. Since the accuracy of parameter estimates was not equivalent across treatment groups, the weighting factors applied were heterogeneous. The difference in the slopes of RMY_kL on DMI_kL between parity group was significant (P < 0.001, model (11]) but only the slope within the Multi group was different from zero (P < 0.001, model [101]). The difference in the intercepts of RMY_kL on DMI_kL between parity group was also significant (P < 0.05, model [11]).

View larger version (13K): [in this window] [in a new window]   Figure 7. Plot of RMYkL (slope of weekly decrease in raw milk yield in late stage of lactation) on DMIkL (slope of weekly decrease in DMI in late stage of lactation) for multiparous (•) and primiparous (). Data are linked by trial and sized with respect to weight used in weighted least-squares analysis.

[17]

and the equation for Primi is (model [10]):

[18]

This model was fitted for 68 treatment groups spread across 21 trials and had a residual error of 0.06 kg/wk (residual standard error of random intercept: 0.06 kg/wk). The models [17] and [18] allow direct quantification of the incidence of the kinetic of intake on the kinetic of milk yield in late stage of lactation. This result emphasizes the potential effect that feeding management has on milk yield in late lactation.

Milk fat yield.
Goodness-of-fit was significantly better for Primi (20 treatment groups, RMSE = 27 ± 10 g/d) than Multi (45 treatment groups, RMSE = 40 ± 18 g/d). The distribution of the 65 R² calculated is given in Figure 4 (R² = 91.4 ± 9.0%).

The effect of YOP was significant only for parameter MFY_M (P < 0.1) and only for Primi (MFY_M: + 15 g/yr). A valid detection of an early stage with increasing MFY (t₁ t₀) was only detected for 23 treatment groups, and a valid estimate of MFY_E in early stage was only calculated for 4 treatment groups. Due to the number of parameter estimates, models [8] and [9] were not applied to parameters MFY_C, and MFY_E.

Weighted least squares means of parameters were different between parity groups for all parameters except k_E and t₁ (k_E = 11 ± 3 g/wk and t₁ = 4.0 ± 1.6 wk for pooled parity groups). Multi had a higher level of MFY_M (+ 301 g/d), an earlier time t₂ (–3.2 wk), and a higher slope of decrease in late stage (Primi: k_L = –11 ± 1 g/wk, Multi: k_L = –18 ± 1 g/wk). No correlation was retained within the kinetic of MFY. Nevertheless, the plateau value and the slope of decrease in MFY during late stage were found to be correlated with equivalent parameters for RMY.

The relationship between MFY_M and RMY_M is presented in Figure 8. Since the accuracy of parameter estimates was not equivalent across treatment groups, the weighting factors applied were heterogeneous. The difference in slopes between parity groups was significant (P < 0.001, model [11]), and the difference in intercepts between parity groups was also significant (P < 0.05, model [11]). Random effects on slope, on intercept, and random interaction between random slope and intercept were found to be significant (P = 0.058, P = 0.049, and P = 0.057, respectively) and were therefore incorporated in model [10].

View larger version (12K): [in this window] [in a new window]   Figure 8. Plot of MFYM (plateau value of milk fat yield in middle stage of lactation) on RMYM (plateau value of raw milk yield in middle stage of lactation) for multiparous (•) and primiparous (). Data are linked by trial and sized with respect to weight used in weighted least-squares analysis.

[19]

and the equation for Primi is:

[20]

This model was fitted for 65 treatment groups spread across 18 trials and has a residual error of 6 g/d (residual standard error of random intercept: 16 g/d). Predictions of models [19] and [20] are consistent with observations of peak yield of milk and fat of Schutz et al. (1990) and Stanton et al. (1992).

The relationship between MFY_kL and RMY_kL is presented in Figure 9. Since the accuracy of parameter estimates was not equivalent across treatment groups, the weighting factors applied were heterogeneous. The difference in slopes between parity groups was not significant (P = 0.283, model [10]), and the difference in intercepts between parity groups was also not significant (P = 0.088, model [12]). The regression equation proposed for both parity groups is (model [13]):

[21]

View larger version (11K): [in this window] [in a new window]   Figure 9. Plot of MFYkL (slope of weekly decrease in milk fat yield in late stage of lactation) on RMYkL (slope of weekly decrease in raw milk yield in late stage of lactation) for multiparous (•) and primiparous (). Data are linked by trial and sized with respect to weight used in weighted least-squares analysis.

	DISCUSSION

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION CONCLUSIONS ACKNOWLEDGEMENTS REFERENCES

 
The challenge of metaanalysis is to identify the general features of a system through the study of pooled data from different origins. This approach was selected to examine lactation kinetics.

Each kinetic collected in the database was recorded in a particular experimental context which can be defined by animal factors (breed, genotype, potential, parity, age, health, etc.), feed factors (composition, digestibility, particle size, energy value, etc.), environmental conditions (temperature, calving season, etc.), farm management (housing, feeding, milking, weighing, etc.), experimental design (allocation system, etc.), and analysis techniques (chemical composition of diet, milk constituents content, etc.). Not all these factors were defined in the different experiments used. Moreover, the criteria used to select these data, which were based on the availability of the three pieces of information (DMI, RMY, BW), induced a bias in the sample. The database was used as a collection of various kinetics rather than an experimental metadesign. Hence, our purpose was to quantify the variability of these kinetics between treatment groups. This overall variability includes the between-trial variability (random effect), associated with experimental conditions, and the within-trial variability, associated with experimental factors. Only parity as fixed effect was incorporated in statistical models to deal with this latter source of variability.

The majority of the trials were carried out on Holstein-Friesian cows, but animals differed in size, intake, and milk yield potential. Values of DMI, RMY, and EBW ranged in middle stage of lactation between 10.1 to 24.7 kg, 15.5 to 44.4 kg, and 320 to 573 kg, respectively. The variability was mainly explained by parity and YOP effects. From 1969 to 1999, the DMI_M increased by 6.6 kg, while the RMY_M increased by 15.3 kg and the EBW_M increased by 51 kg. The increase with time reflects genetic selection and improved management. The corresponding annual increase in 44-wk milk yield was of 165 kg.

The correction of BW for gastrointestinal contents was aimed at calculating the unbiased changes in body reserves. Nevertheless, the correction factor of 4 was applied without taking into account the type of diet, especially the ratio of FO to CO. Several equations, already proposed to explain and predict the rumen volume, have been reviewed by Sauvant (1996) and provided very different estimations of rumen content for the same diet. If rumen load and therefore gut fill is more closely related to fiber or FO intake, as several data have suggested (Paquay et al., 1971), the EBW_E would have been overestimated in the present model for higher yielding animals. Moreover, corrections were not made for uterine involution in early lactation, gravid uterus growth in late lactation, and udder growth in early lactation. Gier and Marion (1968) reported a decrease of 8.5 kg in BW due to uterine involution during the first 3 wk of lactation. This could lead to an overestimation in the loss of body reserves in early lactation. In order to more precisely estimate the early depletion in BW due to mobilization, BW must be split into uterus weight, udder weight, total digestive tract weight, and EBW.

The model proposed by Grossman et al. (1999) summarizes the kinetics into main intersecting straight lines and, hence, into a vector of concrete parameters. We proposed to use this analytical approach on all kinetics of performance data during lactation. This allowed the quantification and the linking of patterns of intake, yield, and BW. In the case of milk yield, the overall goodness-of-fit was greater with this model than with the model of Wood-Dhanoa. Considering its goodness-of-fit, its use for each kind of kinetic, and the concrete signification of its five parameters, the model of Grossman could be a satisfying alternative tool for modeling the course of lactation.

The evolution in EBW was not found to be influenced by DMI, RMY or MFY in early lactation. The major effect that explained EBW loss was the level of EBW at calving (model [14]). This result agrees with the review of Broster and Broster (1984), who observed no clear buffering effect regarded as the "cow’s ability to counter nutritional insults by withdrawal of body reserves to support milk production". This result needs to be validated in a context of severe under-nutrition in early lactation. Finally, the knowledge of the BCS at calving could probably improve the prediction of body reserve depletion. The evolution in RMY and in MFY exhibited very different shapes according to parity, especially according to the level of yield in midlactation and the slope of decrease in late lactation. These results are consistent with results of Schutz et al. (1990) and Stanton et al. (1992). In particular, the current data confirmed the well-known influence of parity on lactation kinetics. Models [15], [16], [19], and [20] allowed the prediction of plateau values of RMY and MFY according to DMI and RMY respectively, by parity. Models [17] and [21] focused on dynamic aspect of the lactation through the relationships between slopes of decrease in intake and milk yield, and milk and fat yield. The choice of DMI as a predictor for RMY is coherent with a "push" viewpoint of lactation (Mertens, 1996) for which milk secretion is considered as a consequence of energy intake. From the opposite "pull" viewpoint, intake is a response to requirements. Only mechanistic models could integrate regulation systems gathering together "push" and "pull" concepts.

	CONCLUSIONS

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION CONCLUSIONS ACKNOWLEDGEMENTS REFERENCES

 
This metaanalysis has highlighted the importance of three main types of effect in the dynamics of lactation: general scale relationship within the couple (intake, yield), incidence of animal size on mobilization of body reserves in early lactation, and effect of lower plane of nutrition on decreasing milk yield and associated decrease in fat yield in late lactation. This work has also opened up the possibility of building an empirical model of lactation kinetics by use of parameter vectors of the model proposed by Grossman (1999). Further trials should be included in the database to validate the regression models obtained.

	ACKNOWLEDGEMENTS

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION CONCLUSIONS ACKNOWLEDGEMENTS REFERENCES

 
The authors thank Dr. P. Schmidely, Dr. V. Berthelot, and Dr. S. Robin for suggestions and assistance.

Received for publication October 3, 2001. Accepted for publication June 17, 2002.

	REFERENCES

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION CONCLUSIONS ACKNOWLEDGEMENTS REFERENCES

 


Broster, W. H., and V. J. Broster. 1984. Reviews of the progress of dairy science: Long term effects of plane of nutrition on the performance of the dairy cow. J. Dairy Res. 51:149–196.[Medline]

Broster, W. H., V. J. Broster, and T. Smith. 1969. Experiments on the nutrition of the dairy heifer. VIII. Effect on milk production of level of feeding at two stages of the lactation. J. Agric. Sci. 72:229–245.

Broster, W. H., V. J. Broster, T. Smith, and J. W. Siviter. 1975. Experiments on the nutrition of the dairy heifer. IX. Food utilization in lactation. J. Agric. Sci. 34:173–186.

Chilliard, Y., B. Rémond, J. Agabriel, J. Robelin, and R. Vérité. 1987. Variations du contenu digestif et des reserves corporelles au cours du cycle gestation-lactation. Bull. Tech. C.R.Z.V. Theix, INRA 70:117–131.

Cobby, J. M., and Y. L. P. Le Du. 1978. On fitting curves to lactation data. Anim. Prod. 26:127–133.

Cowan, E. D., J. Oliver, and R. C. Elliott. 1970. Complete diets for dairy cows. 1. The health, reproductive performance, voluntary food intake and milk yield of cows fed with diets of different roughage content throughout lactation. Rhod. J. Agric. Res. 8:15–22.

Decaen, C., and J. Adda. 1970. Evolution de la sécrétion des acides gras des matières grasses du lait au cours de la lactation de la vache. Ann. Biol. Anim. Biochim. Biophys. 10:659–677.

Deluyker, H. A., R. H. Shumway, W. E. Wecker, A. S. Azari, and L. D. Weaver. 1990. Modeling daily milk yield in Holstein cows using time series analysis. J. Dairy Sci. 73:539–548.[Abstract]

DePeters, E. J., J. F. Medrano, and D. L. Bath. 1989. A nutritional evaluation of mixed winter cereals with vetch utilized as silage or hay. J. Dairy Sci. 72:3247–3254.