Estimation of 305-day milk yield from test-day records of Chinese Holstein cattle

ABSTRACT This study compared six models, namely the Gaines, Sikka, Nelder, Wood, Dhanoa and Hayashi models, for the estimation of 305 days milk yield in Chinese Holstein cattle. We compared their ability to reliably predict 305-day lactation yield from incomplete (3 or 6 test-day (TD)) records. Our findings revealed that the accuracies (ACC) were 0.6655–0.9948, 0.8652–0.9977 and 0.9169–0.9968, whereas the mean square errors (MSE) were 0.0121–2.4807, 0.0139–1.0716 and 0.0170–0.5528 when 3 TD records were used in the first, second and higher lactations, respectively; when 6 TD records were used, the ACC were 0.8800–0.9992, 0.8742–0.9998 and 0.7950–0.9996, whereas the MSE values were 0.0017–0.3348, 0.0011– 0.8605 and 0.0021–1.4869 in the first, second and higher lactations, respectively. All the models were fitted more accurately with 6 TD than 3 TD records. Further analysis revealed that the curves made by the Nelder, Wood and Dhanoa models were close to the actual curves. These three models can be used to predict the 305-day yield for management decisions in farms and for the genetic evaluation of Chinese Holstein cattle.


Introduction
The 305-day milk yield is defined as a cow's milk yield from day 1 to day 305 of the lactation period. This parameter is an important basis for the selection and elimination of cows during the production process and for individual genetic evaluation. The genetic evaluation of dairy bulls and cows based on test-day (TD) yield by TD models has been widely studied (Strabel et al. 2005;Bohmanova et al. 2008). However, the total 305day yield is still the current basis for the genetic evaluations of dairy cattle and the information is used by producers to make their management and breeding decisions (Dongre and Gandhi 2014). The first attempts to mathematically represent the lactation curve and estimate the 305-day milk yield were made by Brody et al. (1923). Hence, several researchers have attempted to develop various models to describe the lactation curve, which can be classified into two main groups: empirical conceptions and mechanistic conceptions. The empirical conceptions models include the parabolic exponential model (Sikka 1950), the incomplete gamma model (Wood 1967), the Wilmink model (Wilmink 1987), the polynomial regression model (Ali and Schaeffer 1987) and the inverse polynomial model (Nelder 1966). The mechanist conception models describe the lactation curve with a biologically based approach; examples of these models are the Dijkstra model (Dijkstra et al. 1997), the Rook model (Rook et al. 1993), the Pollott model (2000) and artificial neural networks (ANNs) or the MilkBot® model (Ehrlich 2011). ANNs allow for the early and accurate predictions of milk production in cows (Grzesiak et al. 2003;Hosseinia et al. 2007;Sharma et al. 2007;Njubi et al. 2010). However, these models have yet to be applied in China to predict the 305-day milk yield with Dairy Herd Improvement (DHI) data to date. Consequently, we selected the more popular models to describe the course of lactation, namely those of Gaines (1927), Sikka, Wood, Nelder, Dhanoa and Hayashi et al. (1986). This study aimed to estimate the parameters of the six nonlinear models and to compare their performances. The best model could be applied in practice as an alternative to the currently used test interval method (TIM) (Everett and Carter 1968) with minimum mean square error (MSE) and maximum accuracy (ACC), and the curve predicted was close to the actual curve in three groups, because TIM method has to all records of at least 200 days in length were regarded as complete lactations. This work is important for developing effective selection strategies and regulating management practices in a herd.

Data sources
The dataset A was provided by the Shandong Dairy Herd Improvement Centre. The initial dataset consisted of 700,715 TD records of 73,130 Holstein cows from June 2011 to December 2014. The data were edited based on the following criteria: (a) TD records corresponded to 5-305 days milk yield (DIM); (b) the number of TD records for each lactation should be 10 times once a month; (c) the first reported TD could not exceed 35 DIM; (d) the time interval between consecutive TD must be between 28 and 35 days and (e) the individual daily milk yield should be between 5 and 50 kg. The final dataset contained 57,790 records, which belonged to three groups representing the first (24,290 records), second (16,100 records) and higher (17,400 records) lactations. The main statistical characteristics of these groups with respect to the daily milk yield are shown in Table 1. The average TD yield was 23.52 kg, whereas the peak yield and peak day were 31.15 kg and 70.02 days, respectively.
The B dataset was used as the control and consisted of 2,550,931 TD milk yield records of 6654 Holstein cows provided by two farms, where milk was automatically weighed with a DeLaval machine. The data were edited based on the following criteria: (a) the individual daily milk yield should be between 5 and 50 kg for TD records corresponding to the range of 5-305 DIM; (b) the number of TD records for each lactation should be more than 250; (c) the first reported TD could not exceed 35 DIM and (d) the abnormal and missing data can be substituted with the arithmetic average value of the herd on the same daily milk yield in the farm. The final dataset contained 1619 cows with 544,765 records, which belonged to three groups representing the first (458 cows, 134,300 records), second (1069 cows, 332,408 records) and higher parity (92 cows, 78,057 records) cows. The abovementioned datasets were derived from a farm which followed total mixed ration and batter management in practical production.
The C dataset was randomly selected from 10 TD every once a month from the B dataset of each cow, which were used to estimate the 305-day milk yield based on the first 3 or 6 TD records and evaluate the 6 studied models. The final C dataset contained 1619 cows with 19,120 records, wherein 4580 records were from 458 cows in the first lactation, 11,940 records were from 1069 cows in the second lactation and 2600 records were from 92 cows in the higher lactations.

Lactation equations
Six models have been used to describe the lactation curve and represent empirical models of varying levels of complexity as compared in Table 2. The first model from Gaines was an early attempt to model the lactation curve. The second model is a parabolic exponential model by Sikka (1950). The third model by Nelder (1966) was derived from the Sikka model with an inverse exponential parabolic function. The fourth model by Wood (1967) proposed the widely applied gamma equation. The fifth model by Dhanoa used the differences in the coefficients of the Wood model. The sixth model for describing lactation curves has been proposed by Hayashi et al. (1986).

Statistical evaluation
The parameters of the nonlinear lactation curves were estimated by a Gauss-Newton iterative method with the nonlinear regression NLIN procedure of the SAS software. Different statistical tests were used for ranking and evaluating models. Occasionally, the results from the different tests seemed contradictory; thus, an overall assessment is needed in this situation. The MSE, goodness of fit (R 2 ) and accuracy (ACC) were used to evaluate the models in study.
where SST is the total sum of squares, which is defined as the sum of the squared differences of each observation, over all observations, from the overall mean. ESS is the sum of the squares of the differences of the predicted values and the mean value of the response variable. y i is the predicted value, y i is the actual value and STD yi is the standard deviation of the actual value. The R 2 value is an indicator which measures the proportion of total variation in the mean as explained by the lactation curve model. The coefficient of determination is always between 0 and 1; the fit of a model is satisfactory if R 2 is close to unity.
The MSE is a measure of the closeness of a fitted line to the actual points. For every actual point, we vertically obtained the distance from the points to the corresponding y value on the fitted curve; this value is known as the error and is squared. Subsequently, we added all these values for all data points and divided the sum by the number of points. The error is squared so that the negative values do not cancel out the positive values. The smaller is the MSE, the closer is the fit to the data.
The MSE and ACC values were the most important criteria to compare the suitability of the lactation curve models in terms of the lactation milk yield properties. Therefore, the best model is the one with the lowest MSE but the highest R 2 coefficient and ACC (Liu et al. 2000). Notes: 1, means first lactation; 2, means second lactation and 3+, means higher lactations. Table 2. Equations to describe the lactation curve of dairy cows.

Model Equation Gaines
y t = ae −bt Sikka y t = ae (bt−ct 2 ) Nelder y t = t/(a + bt + ct 2 ) Wood y t = at b e −ct Dhanoa y t = at bc e −ct Hayashi Notes: y is the milk yield (kg/d); t is the time of lactation (d); e is the base of natural logarithm and a, b and c are parameters that define the scale and shape of the curve.

Statistical evaluation
Six models were evaluated with the collected A dataset. Table 3 gives the various parameters and summarizes the key statistical measures used to compare the performance of the models. Each model was separately evaluated for the first, second and higher lactations. For first-lactation cows, the R 2 by the six models was between 0.9195 and 0.9212 and the MSE was 0.0019 for all. For second-lactation cows, the R 2 was between 0.9114 and 0.9186 and the MSE was 0.0034 for all models, except the Hayashi model, which had an MSE of 0.0037. For higher-lactation cows, the R 2 was between 0.8979 and 0.9062 and the MSE was 0.0033 for all models, except that of Hayashi. As shown in Table 3, the Nelder equation performed best, with the maximum R 2 and the minimum MSE. The performance of the Hayashi model was relatively lower than those of the other models to produce the maximum MSE and the minimum R 2 . The different models can be fitted to the data without difficulty through nonlinear regression. The fitting effect decreased with increasing lactation in Table 3.
For Chinese dairy cows, the six equations have minimum MSE and maximum R 2 in the first lactation in general, genetic effects play an important role in the first lactation; this model can better respond to the milk yield.

Incomplete lactations
Six models were evaluated with the B and C datasets. A comparison of their predictive ability allowed the identification of a model capable of describing and providing a better perspective on the shape of the lactation curve of Holstein dairy cows with 3 or 6 TD records. The results of the predicted ACC and MSE values are shown in Table 4. The actual 305-day milk yield of the first, second and higher lactation groups were 8453, 11,016 and 9069.79 kg, respectively. The predicted 305day milk yield averaged between 3287.51 and 14,476.79; 6215.01 and 11,643.74; and 3287.51 and 9337.07 kg in the first, second and higher lactations, respectively. When 6 and 3 TD records were used, the results in Table 4 showed that the ACC ranged between 0.7950 and 0.9998, 0.6655 and 0.9977, whereas the MSE ranged from 0.0011 to Table 3. Parameter estimates and other measures when models were fitted to the first, second and higher lactations cows in dataset A.  The predicted milk yield during the lactation with 6 TD records. b The predicted milk yield during the lactation with 3 TD records.
1.4869 and 0.0121 to 2.4807, respectively. The ACC values calculated with 6 TD records were 0.8800 and 0.9992, 0.8742 and 0.9998, 0.7950 and 0.9996, whereas the ACC values calculated with 3 TD records were 0.7936 and 0.9650, 0.8652 and 0.9976, 0.9169 and 0.9771 in the first, second and higher lactation groups, respectively. In all groups, the ACC was higher with 6 TD records than with 3 TD records. For example, the ACC for the Gaines model was reduced by 0.3099 with 3 TD records than 6 TD records in the first lactation. However, an exception to this trend was the Hayashi model, which had opposite results in higher lactations. Similar results were reported by Aziz et al. (2006). The MSE values calculated with 6 TD records were 0.0017 and 0.3348, 0.0011 and 0.8605, 0.0021 and 1.4869, whereas the MSE values calculated with 3 TD records were 0.0121 and 2.4807, 0.0139 and 1.0716, 0.0170 and 0.5528 in the first, second and higher lactation groups, respectively. With 3 TD records, the maximum ACC was 0.9977 (Wood model in second lactation) and the minimum MSE was 0.0121 (Wood model in first lactation). With 6 TD records, the maximum ACC was 0.9998 (Dhanoa model in second lactation) and the minimum MSE was 0.0011 (Wood model in second lactation). Figures 1 and 2 illustrate the TD yield in the groups as predicted by six models over a full cycle (305 d) with 3 or 6 TD records, respectively, which could be compared with the actual daily yield. The most common shape was a rapid increase after calving to a peak after a few weeks, followed by a gradual decline until the cow was dried off. Basing from Figures 1(A-C) Figure 1. Average daily milk yield, actual and predicted values from 3 TD records by 6 models (Gaines, Sikka, Nelder, Wood, Dhanoa and Hayashi). and 2(D-F), we compared six models to the actual lactation curve. The curve by the Gaines model was always initially straight in the first, second and higher lactations. The Gaines model was simplistic and did not provide a physiological basis for the lactation curve (Korkmaz et al. 2011). However, the accuracy of the 305-day milk yield predicted by the model was above 0.9 in Table 2.
The lactation curve of the Sikka model was identical with the actual curve before the peak milk yield but decreased after the peak (about 101-305 DIM) when 3 TD records were used, as shown in Figure 1(A-C); after the peak, the predicted value was lower than the actual value in all groups. When 6 TD records were used, the lactation curve of the predicted daily milk yield was higher than the actual curve before ∼180 DIM and lower after ∼180 DIM in the first lactation. In the second lactation, the curve of the predicted daily milk yield was lower than the actual curve. In higher lactations, the curve was higher than actual curve before 30 DIM and close to actual curve afterward as shown in Figure 2(D-F). Hossein-Zadeh (2015) reported slightly different results, such that the fitted milk yield was only reasonable during the first lactation.
For cows with 3 TD records, the lactation curves by the Nelder model fitted the lactation law before ∼100 DIM and became higher than the actual curve after 105 DIM in first and second lactations. However, in higher lactations, the curve of the Nelder model was lower than the actual curve after 105 DIM. Similarly, for cows with 6 TD records, the curve was higher than the actual curve after 205 DIM in the first Figure 2. Average daily milk yield, actual and predicted values from 6 TD records by six models (Gaines, Sikka, Nelder, Wood, Dhanoa and Hayashi). lactation; in the second lactation, the predicted curve was higher than the actual curve before 30 DIM and after 155 DIM, but the prediction was lower between 30 and 140 DIM. In higher lactations, the predicted curve was higher after 225 DIM. The Nelder model slightly overestimated (overprojected) the lactation yields, and the researchers predict that the ACC will increase with rising parity. The Nelder model was the best option for describing the lactation curve in Figures 1(A-C) and 2(E-F), with the maximum R 2 and the minimum MSE in Table 2. This high level of accuracy has also been reported in previous studies of fitting models to mean yields (Hossein-Zadeh 2015).
The curve predicted by the Wood model in all groups was close to the actual curve when 3 or 6 TD records were used. However, the curve in the first lactation was close to the actual curve before 175 DIM with 3 or 6 TD records but became lower thereafter. Similar findings were described by Strabel et al. (2004), Macciotta et al. (2005) and Cankaya et al. (2011). The Wood model provided a physiological basis for the lactation curve, and the fit of the model was satisfactory (Cunha et al. 2010). However, the Wood model has been shown to slightly underestimate the daily milk yield (Tozer and Huffaker 1999;Cilek and Keskin 2008;Banu et al. 2012;Dohare et al. 2014). The Wood model serves as a guide for further modelling and provides novel insights. In general, the Wood equation provided a similar goodness of fit as the Dhanoa model (Korkmaz et al. 2011), although several studies have shown differences in the general shape of the lactation curve (Landete-Castillejos and Gallego 2000; Orhan and Kaygisiz 2002).
The curve in the first lactation as predicted by the Hayashi model was appreciably higher than the actual curve before the peak but was higher after the peak, regardless of whether or not 3 or 6 TD records were used. During the second lactation, the curve predicted with 3 TD records was lower than the actual curve before 100 DIM but was higher after the peak (Figure 1 (B)). By contrast, the curve predicted with 6 TD records was lower before the peak yield of 130 DIM but finally became higher until the end of lactation (Figure 2(E)). Previous studies showed the reverse trend in the shapes of the lactation curve of Egyptian buffalo. However, the model cannot fit the cow lactation curve in higher lactations with 3 or 6 TD records. The overall performance of the Hayashi model was lower compared with that of the other models.

Conclusion
In the study, monthly records from the DHI were used to estimate the 305-day milk yield with six models (Gaines, Sikka, Nelder, Wood, Dhanoa and Hayashi) and controlled with actual records from farms. The predicted 305-day milk yield was close to the actual values with the Nelder, Wood, Dhanoa, Gaines, Sikka and Hayashi models, which are listed with decreasing accuracy. The curves made by the Nelder, Wood and Dhanoa models were close to the actual curve for all groups, whereas that of the Hayashi model only fit during the first lactation. These models can be used to predict the 305-day milk yield and further assist farm management decisions and the genetic evaluation of Chinese Holstein cattle.

Disclosure statement
No potential conflict of interest relevant to this article was reported.