Comparison of models for describing the lactation curve of Awassi , Morkaraman and Tushin sheep

The aim of this study was to identify a suitable mathematical model for describing the lactation curve of Awassi, Morkaraman and Tushin sheep breeds and to determine breed differences. Data on milk yield of 182 Awassi, 47 Morkaraman and 74 Tushin ewes were used. Eight empirical models from the literature were used to fit the standard lactation curves. Among them the Wood model (WD) appeared the most appropriate according to mean square prediction error (MSPE), coefficient of determination (R2), Durbin-Watson statistic (DW), and its applicability to the data for all three breeds. There were statistically significant (P<0.05) differences among Awassi, Morkaraman and Tushin breeds in accordance with a, b and c parameters and peak yield. The Awassi breed had the highest peak yield and the Morkaraman and Tushin breeds had statistically similar lower peak yields. There were no significant differences among the parameters of the WD model except for peak yield and peak time in accordance with parities. Breed and parity interaction was significant (P<0.05) only for peak yield.


Introduction
The milk produced by dairy animals is crucial for feeding people.Although mainly cow milk is globally consumed for this purpose, a small amount of sheep and goat milk is also used.Especially in the Mediterranean, Middle East and Eastern Europe, milk is an important product from sheep (POLLOTT and GOOTWINE 2000).Sheep milk is an important product for Turkey.Total milk production of Turkey is about 8.4 million tons per year.Approximately 7.8 % of total milk production is produced by sheep (ANONYMOUS 2003).Lactation curves have been heavily studied in cows for a long time.However, few studies are available in sheep.Lactation curves have a wide variety of applications, such as extension of incomplete records for use in genetic evaluations, formulation of rations and economic evaluation of different management schemes (GRONEWALD et al. 1995, SWALVE 1995, HORSTICK and DISTL 2002, GIPSON and GROSSMAN 1990, SAKUL and BOYLAN 1992, DAG et al. 2005, KESKIN and DAG 2006).A lactation curve is particularly useful to make a decision that is time-dependent.If we know when to expect an animal to reach peak yield, we can arrange the feeding strategy to allow economic management of feed to the extent that would satisfy the animal's requirement during various stages of lactation, reduce cost, and possibly maintain peak yield for as long as possible (BAFFOUR-AWUAH et al. 1996, AMIN 2003, GRZESIAK et al. 2003).Appropriate mathematical models may also be used to predict future milk yield for lactations currently in progress.
Several models of lactation curves have been developed.The incomplete gamma function, or Wood's model, developed by WOOD (1967) is the most well known and commonly employed mathematical model, specifically applied to cattle (SWALVE and GUO 1999).The vast majority of the other models are also based on Wood's model (COBBY andLE DU 1978, WILMINK 1987).They usually include additional parameters to improve fitness and have been developed from experimental data.For example, other alternatives comprise a four-parameter model by Morant (MORANT and GNANASAKTHY 1989) and a six parameter diphasic function by Grossman (GROSSMAN and KOOPS 1988).More proposed functional forms of lactation curves can be found in ROOK et al. (1993).These alternative models have traditionally been applied to lactation data for dairy cattle, but less often in sheep and goats (KOCAK and EKIZ 2008).To model the lactation curve of sheep most studies have been based on Wood's model.However, some other models are used to fit sheep lactation data.NELDER (1966) suggested an inverse polynomial model (NE) be fitted to lactation data for dairy cattle and his model has occasionally been used for sheep and goats.This function has a value of zero as its initial yield and is also intrinsically linear.This was followed by Wood's three parameter model (WD), which has been widely used, especially in dairy cows and in sheep and goats.Wood proposed an incomplete gamma function in which constants can be interpreted biologically (WOOD 1967).Of these constants, a is initial milk yield just after lambing, and b and c are slopes of the curve before and after, respectively, the peak yield of lactation.Peak yield is a(b/c)b exp(-b) and time of peak yield is b/c.SCOTT et al. (1996) reported that WD has a tendency to overestimate milk yield prior to peak yield and in late lactation.COBBY and LE DU (1978) proposed another model (CD) in which an exponential decline and a linear decline are combined.The peak yield is estimated as c-In(ac/b) and there is a linear decline after that point.DHANOA (1981) reparameterised Wood's model (DH) to eliminate the problems based on its parameters.Parameter b in this model represents the time required to reach peak yield.WILMINK (1987) introduced a lactation model (WI) in which a quadratic term was added and an adjustment made to the exponential term.In this model, a may be interpreted as the level at which yield begins, b as the decrease after peak yield is reached and c as the initial increase to peak.The factor w was set equal to 0.05 based on the literature.MORANT and GNANASAKTHY (1989) adopted a different approach (MG) by considering the proportional rate of change in milk yield during lactation.They addressed the problem of a high level of correlation between the parameters of the incomplete gamma function and proposed a modified polynomial model with four parameters, which can also be used in log linear form, to overcome these correlations.GIPSON and GROSSMAN (1989) introduced a diphasic model (GR) to overcome problems of autocorrelation detected in models based on the gamma function.A diphasic logistic function has some advantages over the gamma function, which are due to smaller and more random residuals, prediction of total milk yield and meaningful functions of easily interpretable parameters that have biological importance as observed by GROSSMAN and KOOPS (1988).A non-linear modification of Wood's equation was developed by CAPPIO-BORLINO et al. (1995) (CB), which is especially appropriate when milk production dramatically decreases instantly after the lactation peak.
The aim of this study was to identify a suitable mathematical model for describing the lactation curve of Awassi, Morkaraman and Tushin sheep and to determine breed differences.

Experimental materials and procedures
The experiment was conducted at the Research and Application Farm of the College of Agriculture, Ataturk University, Erzurum, Turkey.The experiment was carried out with the lactations of 182 Awassi, 47 Morkaraman and 74 Tushin ewes.Age at first lambing was approximately 24 months.All ewes lambed in March.All lambs suckled their dams freely until first milk recordings.Ewes were hand-milked twice daily (morning plus evening) and the first milk test was performed in the first month after lambing.Lambs were kept with their dams until 2.5 months of age, the time when they were weaned completely from milk.All breed groups were treated similarly with respect to diet and management.
Approximately 87 % of Turkey's sheep population (25.4 million heads) is fat-tailed breeds (ANONYMOUS 2003).The sheep population of the eastern and northeastern parts of Turkey consists predominantly of the fat-tailed Morkaraman breed.
These sheep are thought to have evolved through natural selection under harsh environmental conditions.For centuries, eastern Anatolian farmers have used the Morkaraman breed to produce milk, meat and wool.It is estimated that there are about 5.09 million heads (ANONYMOUS 2003).The lactation period is 140-150 days.The average total milk yield is about 50-65 kg.
Awassi is a fat-tailed breed reared extensively in southern part of Turkey and it is estimated that there are about one million Awassi sheep in Turkey.Awassi sheep breed has been imported from southern part to eastern part of Turkey to increase milk and meat production of Morkaraman.It has been raised at Research and Application Farm of Agricultural Faculty since 1974.The Awassi is kept for milk, meat and wool production.The lactation period is 170-200 days.The average total milk yield is about 100-150 kg.Tushin sheep are also fat-tailed breed and it is estimated that they are about 120 000 heads (ANONYMOUS 2003).The meat of Tushin sheep is known for its palatability and a nice aroma free of the specific acrid smell of mutton.Milk yield is not high, 50-60 kg per lactation.The lactation period is 100-140 days.
Eight empirical models from the literature (Table 1) were used to fit the standard lactation curves in this study, where Y t is average daily milk yield in the t-th week of lactation and a, b, c and d are parameters.Models were fit to individual lactation data for each sheep.
To describe the lactation curves the parameters were estimated by fitting eight functions to the data.All analyses were conducted using the Statistical Analysis System (SAS) nonlinear least-squares regression (NLIN) procedure (SAS 2000).The Gauss-Newton algorithm was chosen for iterations of nonlinear fit.The maximum number of iterations used was 100 and the convergence criterion was (SSEj-1-SSEj)/(SSEj+10-6)<10-8, where SSE is the residual sum of squares after fitting the function to the data, and j denotes the round of iteration, as defaults for SAS.When needed, extra iterations were performed to provide convergence.To compare the suitability of the models, the following process was used: (1) Parameter values for each model and their standard errors were derived from the iterative process.
(2) Goodness-of-fit for each model was evaluated using mean square prediction error (MSPE) and the adjusted coefficient of determination: where R 2 is the multiple coefficient of determination, n is the number of observations, and p is the number of parameters in the model.Note that R 2 is adjusted for the number of parameters in the model to make a fair comparison of models, R 2 adj will be reported as R 2 for simplicity.
(3) The Durbin-Watson statistic (DW) was used as a measure of first order positive autocorrelation to test whether the residuals were randomly distributed (DURBIN and WATSON 1951): ∑ n (e t -e t-1 ) 2 t=2 (2) where e t is the residual at time t, and e t-1 the residual at time t-1.The observed value of DW was evaluated against the tabulated critical value to test for positive autocorrelation.Negative autocorrelation was not tested, because a negative autocorrelation coefficient of −1 implies that residuals fluctuate in a strict up and down rate around the actual curve, which in the case of lactation curves is not a problem (FERNáNDEZ et al. 2002).

Results and discussion
The goodness of fit statistics for the models tested is presented in Table 2.As stated before, all models had high R 2 values, ranging from 0.703 to 0.937, except for NE, suggesting overall good fits to the data.In general terms, apart from NE, the models provided accurate fits to the data as judged by MSPE, adjusted R 2 , and DW.The best fit was obtained with MG and was followed very closely by GR and WD.Again, the worst fit was found for NE.The MG, GR, and WD models had the lowest MSPE and a higher DW coefficient when compared with the others models.The nonlinear MG model gave better fit to the lactation data, having the lowest values of MSPE, and the highest values for the R 2 and the DW coefficients, meaning the fewest positive autocorrelated cases (Table 2).On the other hand, POLLOTT and GOOTWINE (2000) reported that the Morant function always gave the lowest MSPE in Awassi sheep in Israel.
The NE model gave the poorest fit according to MSPE and the R 2 .Although the MG model gave the smallest MSPE and the highest R 2 , it turns out to be inappropriate for representing the lactation curves of all individual animals in all three breeds (Table 2).This makes the curve unreliable.The GR model proved unsatisfactory for individual animals due to over fitting and resultant unrealistically shaped curves in Morkaraman and Awassi after MG.There was very little difference between MG and WD models, with MG having marginally smaller MSPEs.In addition, it is said that early in lactation the diphasic model results in a poor fit, because the hyperbolic tangent requires symmetry in both phases, and, when only two phases are considered, a symmetric curve does not fit the possible steep rise that occurs early in lactation (FERNáNDEZ et al. 2002, KOCAK andEKIZ 2008).
Given the overall results of these analyses, the WD model was determined to provide the best fit to our data for all three breeds.Averaged predictions of the Wood Model against observed milk yields for the three breeds are given in Figure 1-3.A number of other researchers also have tried to test the reliability and adaptability of the Wood model to the lactation curve of sheep.The general conclusion was that the Wood model may be useful to estimate milk production in ewes with a good level of approximation (CAPPIO-BORLINO et al. 1989, CAPPIO-BORLINO et al. 1995, SAKUL and BOYLAN 1992, PORTOLONA et al. 1996).
∑ n e t 2 t=1 DW = In the Figures 1, 2, and 3 it can be seen that the Awassi breed produced milk more than the others in a longer lactation period.The general conclusion was that the Wood model may be useful to estimate milk production in ewes with a good level of approximation (CAPPIO-BORLINO et al. 1989, CAPPIO-BORLINO et al. 1995, SAKUL and BOYLAN 1992, PORTOLONA et al. 1996).Given the overall results of these analyses, the WD model was determined to provide the best fit to our data for all three breeds.Averaged predictions of the Wood Model against observed milk yields for the three breeds are given in Fig 1-3.A number of other researchers also have tried to test the reliability and adaptability of the Wood model to the lactation curve of sheep.The general conclusion was that the Wood model may be useful to estimate milk production in ewes with a good level of approximation (CAPPIO-BORLINO et al. 1989, CAPPIO-BORLINO et al. 1995, SAKUL and BOYLAN 1992, PORTOLONA et al. 1996).In the Figures 1, 2, and 3 it can be seen that the Awassi breed produced milk more than the others in a longer lactatio The results of the Analysis of Variances of Wood model parameters according to breed and lactation parity are presented in Table 3.As seen, there were statistically significant (P<0.05)differences among Awassi, Morkaraman and Tushin breeds in accordance with a, b and c parameters and peak yield.The a parameter revealed the highly significant effect of breed groups.Awassi and Tushin were higher than Morkaraman.The b constant which measured the average slope of the lactation curve indicated by period of optimum daily production was significantly affected by breed.It was the lowest in Tushin.The c constant which describes the declining extreme of the curve was significantly affected by breed.It was the lowest in Tushin.The Awassi breed had the highest peak yield and the Morkaraman and Tushin breeds had statistically similar lower peak yields.Although it is not statistically significant, Morkaraman breed had the highest peak time.There were no significant differences among the parameters of the WD model, except for peak yield and peak time in accordance with parities.Breed and parity interaction was significant (P<0.05)only for peak yield (Figure 4).A positive relationship between scale and parity was observed, that is, first parity peak yields were lower than for later parities, similar to the findings of GIPSON and GROSSMAN (1990).Similar conclusions were also reported by RUIZ et al. (2000), andFERNáNDEZ et al. (2002).
In this study an attempt was made to find the mathematical model that most satisfactorily represents the lactation curves of Awassi, Morkaraman and Tushin sheep.Eight non-linear models were compared.The Wood model was considered as most accurate to describe lactation curves, according to R 2 , MSPE and DW criteria, and its applicability to the data for all three breeds.Among the three breeds studied significant differences were observed for the parameters of the Wood model and peak yield.The peak yield and peak time were significantly different among parities.GIPSON and GROSSMAN (1990).Similar conclusions were also reported by RUIZ et al. (2000), andFERNÁNDEZ et al. (2002).In this study an attempt was made to find the mathematical model that most satisfactorily represents the lactation curves of Awassi, Morkaraman and Tushin sheep.Eight non-linear models were compared.The Wood model was considered as most accurate to describe lactation curves, according to R 2 , MSPE and DW criteria, and its applicability to the data for all three breeds.Among the three breeds studied significant differences were observed for the parameters of the Wood model and peak yield.The peak yield and peak time were significantly different among parities.

Figure 1
Figure 1 Observed and predicted yields using WD model for Awassi Beobachteter und vorhergesagter Ertrag beim WD-Modell für Awassi Figure 1 Observed and predicted yields using WD model for Awassi Beobachteter und vorhergesagter Ertrag beim WD Modell für Awassi Figure 1 Observed and predicted yields using WD model for Awassi Beobachteter und vorhergesagter Ertrag beim WD Modell für Awassi

Figure 4
Figure 4 Interaction between breed and parity for WD model Interaktion zwischen Rasse und Parität beim WD Model Figure 4 Interaction between breed and parity for WD modelInteraktion zwischen Rasse und Parität beim WD Model

Table 1
The models for describing lactation curves Modelle zur Beschreibung von Laktationskurven

Table 2
Comparison of the goodness-of-fit statistics of several models to describe the lactation curves of Awassi, Morkaraman and Tushin Vergleich der Einigkeitsübereinstimmung Statistiken mehrerer Modelle für die Definierung der Laktation kurve der Awassi, Morkaraman und Tushin

Table 3
Least square means and standard errors for lactation curve parameters, peak yield and peak time for WD model Geringste quadratische Abweichung und Standardfehler der Laktationskurvenparameter, Höchstleistung und Dauer beim WD-Modell