Abstract: Heterogeneity of the residual variances was investigated for test day milk yields of Turkish Holsteins. A third order random regression model including the fixed, random additive genetic, permanent environmental and residual effects were used. Residual Variances (RV) of test day milk yields divided into different classes with 5 schemes as RV1, RV3, RV5, RV7 and RV10. Also, a structural model was fitted to residual variances of milk yields. The residual variances were estimated from 2.19-29.47 kg² and additive genetic variances were estimated 3.19-10.97 kg². The heritability estimates ranged from 0.12-0.81. The third order random regression models with RV7 and RV10 schemes were found to be better than others for evaluation of residual variances.

INTRODUCTION

In animal production, repeated measures such as milk yield, body weight and feed intake have been often analyzed with test day models. An improved type of test day models is a random regression model that allows defining individual test day variation in the shape of the trajectory with a set of random regression coefficients. In the random regression model, residual effect is generally assumed constant throughout the trajectory. However, recent studies have shown that the Residual Variance (RV) has changing over time (Jaffrezic et al., 2000) because of several environmental factors such as region, herd management, weather conditions, lactation number, age at calving, month of calving, days in milk, pregnancy status, medical treatments and milking times etc. (Perochon et al., 1996; Swalve, 2000). These factors make the residual variances larger and more variable at the beginning and at the end of the lactation than those at middle lactation.

Assuming a homogeneous residual variance would directly affect the genetic evaluation. As stated by Rekaya et al. (2000), the impact of heterogeneity of the residual variance on evaluation goes through the weight given to information in each part of the lactation. If homogeneous variance is assumed, information from parts of the lactation having lower residual variance would receive lower weight and higher residual variance would receive higher weight. In addition, Olori et al. (1999a) stated that constant residual variance assumption causes residual variances to be underestimated and heritability estimates to be overestimated at early stages of lactation.

Instead of constant residual assumption, lactation has been divided into different classes with assuming homogeneity of residual variance within classes and heterogeneity between them (Olori et al., 1999a; Rekaya et al., 1999; Lopez-Romero et al., 2003). This approach provides information on the expected pattern of residual variance, which changes over lactation.

But a caution should be taken when defining different arbitrary classes. If the identification of the classes in terms of residual variance is not accurate, the proposed model will not be correct (Lopez-Romero et al., 2003). Olori et al. (1999a, b) reported that a correct estimation of residual variance in each class depends on, which lactation stages are assumed to have the same residual variance. Olori et al. (1999a) investigated the effect of constant residual variance on estimation of genetic parameters using a 3rd order orthogonal polynomial random regression model. They considered 4 stages of lactation in terms of heterogeneous residual variances and found larger heritability estimates at the late stages of lactation. Fujii and Suzuki (2006) used 3 kinds of heterogeneous residual variance and homogeneous residual variance over the calving year. They estimated permanent environmental variances were larger than additive genetic variances and the pattern of the permanent and additive genetic variances were very similar in different models.

Moreover, Rekaya et al. (1999) and Brotherstone et al. (2000) assigned residual variances with 30 and 12 classes in a random regression model, respectively. In these studies, residual variances showed larger but heritability estimates showed lower values at the extremes of the lactation. Lopez-Romero et al. (2003) assumed residual variances as constant, 3 or 30 classes in a set of random regression models based on Legendre polynomials with varying order on additive genetic and permanent environmental effects. The assumptions on the RV pattern did not affect the estimates of the daily additive genetic variance and only affected the estimates of daily permanent environmental variance in their 1st part of the lactation.

It can be seen that there is no enough study on definition of classes for residual variance of test day milk yields in the literature. The objective of this study was to compare random regression models, which have different RV schemes in order to account for heterogeneous residual variances of test day milk yields in Turkish Holsteins.

MATERIALS AND METHODS

Data: First lactation test day milk yields of 612 Holstein Friesian cows, raised at the 4 state farms in Turkey, were used in this research. Milk yields were collected at successive monthly periods (TD1-10) from 1987 through 1993. Test day months are used as the time variable rather than days in milk. Data set have limited with following restrictions: Lactation length was limited to have at least 150 days and maximum 308 days long. Therefore, test day milk records <5 kg were excluded from the analyses. Age at first calving was also limited between 18 and 51 mo. In the final data set, a total of 5918 test records were analyzed.

Statistical analysis: In this study, data set were analyzed with a 3rd order random regression Legendre polynomial model. Because Legendre polynomials are orthogonal, normalized and resulted in a better convergence and give more accurate results as compared to conventional polynomials (Kirkpatrick et al., 1990). Third order of polynomial regression was preferred due to the best fit (Van Der Werf et al., 1998; Pool and Meuwissen, 1999; Olori et al., 1999b; Veerkamp and Thompson, 1999; Kettunen et al., 2000; Strabel et al., 2003; Takma and Akbas, 2007).

For the fitting of the model, the DXMRR option of the DFREML statistical package (Meyer, 1998) was used. The mentioned random regression model is as follows:

where:

The residual variances were assumed to be different in several stages of lactation according to the pattern of residual variances for 10 test day milk yields, which have estimated from a univariate analyses done before. The F test for the different residual variance classes within a RV scheme was also performed. According to the significance of F-test (p<0.05), 5 residual variance schemes were constructed as RV1, RV3, RV5, RV7, RV10. The classes and related test day records in these RV schemes were shown in Table 1. In the RV1 scheme, residual variances assumed constant, in the RV3 scheme 3 number of residual variance classes, in the RV5 scheme 5 number of residual variance classes, in the RV7 scheme 7 number of residual variance classes were taken, respectively. In the last scheme (RV10), all test day records were taken as different groups (Table 1).

In addition to different RV schemes, heterogeneity of the residual variances was also modeled with a structural model as In σ²_eij = P'_ijδ, where, P'_ij = (1 t_ij t²_ij) (Foulley et al., 1998). For the quadratic function of time, the structural model would be like this: In σ²_eij = a + bt_ij + ct²_ij.

The goodness of fits for the models with different error variances were examined using likelihood based criteria as Akaike’s information criterion-AIC (Akaike, 1973) and Bayesian information criterion-BIC (Schwarz, 1978). These criteria have been calculated as:

AIC = -2*LogL + 2*p, BIC = -2*LogL +p*log (N-r(x))

where:

Table 1:	The classes design of the models in different residual variance schemes

The model, which gives the lowest AIC and BIC values was chosen as the better approximating model (Apiolaza et al., 2000; Lopez-Romero and Carabano, 2003).

Furthermore, different error structures were compared by Likelihood Ratio Test (LRT) (Rao, 1973). Calculation of LRT for models i and j was obtained with the formula:

LRT_ij = 2*(LogL_i-LogL_j)

In the LRT, the LogL differences were tested using Chi-square (χ²) test with the degrees of freedom determined as the number of the parameter differences between the models (Huelsenbeck and Bull, 1996).

RESULTS

The results of model comparison criterion from different RV schemes and the structural model were presented in Table 2.

An increase in logL and decreases BIC and AIC values with increased numbers of parameters was noticed among different RV schemes. The structural RV model had lower BIC and AIC values than the RV1 and RV3 models and higher than the RV5, RV7 and RV10 models. Although the structural RV model had the same number of parameters with the RV3 model, it gave better performance than the RV3 model (Table 2). However, the RV10 model gave best fit according to the lower AIC values. Moreover, the differences in LogL values of the models with different RV schemes were found significant (p<0.05). The LRT values between the RV1, RV3, RV5, RV7 and RV10 models were found significant (p<0.05). The LRT values between the structural RV model and the RV1, RV3, RV5, RV7, RV10 models were found significant as well (p<0.05) (Table 2).

The residual variances of the random regression models with different RV schemes and the structural RV model were plotted in Table 3. They were changed between 2.19 and 29.47 kg². The predicted residual variances were higher at the beginning of the lactation for all models, describing heterogeneity of RV. The structural RV model and the RV3 and RV5 models have shown different tendency up to the TD4. But, residual variances of the RV7 and RV10 models have shown same tendency over all stages of the lactation.

Estimations of additive genetic and permanent environmental variances and heritability values for the models with different RV schemes were summarized in Table 3. For the additive genetic variances, all models had same tendency and high variances at the end of the lactation. But it can be seen that 2 groups of models were appeared on the shape of the genetic variances at the beginning and the middle of the lactation. First group consists of the RV1, RV3 models and the structural RV model. The second group consists of the RV5, RV7 and RV10 models.

The 1st group had slightly lower genetic variances at early lactation than the 2nd group. However, genetic variance was increased at the middle of the lactation for the 1st group while, it was decreased for the 2nd group. During the late part of lactation (8-10 TD) all models showed similar pattern for genetic variances (Table 3).

The pattern of the permanent environmental variances, which changes over lactation was the same for the all models. The permanent environmental variances were higher at the beginning of the lactation. They were decreased up to the TD3 and increased from the TD3 up to the TD6 and then again they were decreased.

Moreover, heritability estimates were almost same with estimates of additive genetic variances as expected. The heritability estimates for the RV1 and RV3 models inclined to opposite direction of the RV5, RV7 and RV10 models but heritability estimates have increased after TD8 and higher at late lactation for all models.

Besides, the mentioned models, the structural RV model gave overestimated heritability estimates as compared to the heritability estimates of these models (Table 3).

Table 2:	Maximum log likelihood values, BIC and AIC values for different RV schemes
*Consecutive comparisons for LRT values between RV1, RV3, RV5, RV7 and RV10 models are significant (p<0.05), +LRT values between the RV3, RV5, RV7, RV10 models and structural RV model are significant (p<0.05)

Table 3:	Estimations of additive Genetic (G) and Permanent Environmental (PE) variances and heritability estimates (h2) for the models with different RV schemes

DISCUSSION

The aim of this study was to compare random regression models in respect of heterogeneity of residual variances during lactation. Lactation was divided by different RV schemes based on classes assuming homogeneity within classes and heterogeneity between them. This method requires the accurate definition of arbitrary test days, in which the residual variance is assumed constant, whereas the change of the residual variance continues over time (Jaffrezic et al., 2000).

In this study, definition of arbitrary classes has been designed by determining change of residual variances over the lactation with a univariate analysis and analyzing the significance of differences between classes with the F-test. After that residual variances were designed in 5 schemes.

Consecutive comparison of LRT values for the 5 schemes (RV1, RV3, RV5, RV7 and RV10 models) were found significant (p<0.05). Also, LRT values between the RV3, RV5, RV7, RV10 models and structural RV model were significant (p<0.05). Therefore, heterogeneous residual variances were defined as substantial based on increasing values of likelihood criteria with increasing numbers of parameters. On the other hand, there is no discrepancy between the structural RV and the random regression model with constant RV. Although, the structural RV model has the same number of parameters with those of RV3 model, it gave better performance than RV3 model.

Estimates of residual variances of the models were higher at the beginning of lactation and declined in mid lactation as the number of RV classes increased. This might be clarified with increasing ability of models to account for residual variances. Olori et al. (1999a, b) and Jaffrezic et al. (2000) observed a similar decrease for residual variances in models with different orders and different RV schemes. On the other hand, the residual variances were found lower than those of obtained by Rekaya et al. (1999) and Brotherstone et al. (2000) at the end of the lactation. In addition, the structural RV model gave smoother residual variance estimation as found by Jaffrezic et al. (2000) when compared to models with different RV schemes.

Further, estimates of additive genetic variances taken from all of the models (3.19-10.97) were lower than the estimates of Olori et al. (1999a, b) and Lopez-Romero et al. (2003) and higher than the estimates from Rekaya et al. (1999), Jaffrezic et al. (2000), Fujii and Suzuki (2006) and almost similar with the estimates from Brotherstone et al. (2000). On account of pattern for genetic variances, all models showed similar genetic variance patterns during the late part of the lactation (8-10 TD). These findings were similar with Olori et al. (1999a), Jaffrezic et al. (2000) and Fujii and Suzuki (2006), when they modeled different random regression models.

Nevertheless, the RV5, RV7 and RV10 models, which have more RV classes had lower genetic variances at early lactation and higher genetic variances at the middle of the lactation than RV1, RV3 models and the structural RV model. These fluctuations may be explained by having different flexibility of models to describe variances.

The patterns of permanent environmental variances were substantially similar to additive genetic variances of RV1, RV3 and the structural RV models at the beginning and the middle of the lactation. Magnitude of the permanent environmental variances were lower than those informed by Olori et al. (1999a,b), Rekaya et al. (1999), Jaffrezic et al. (2000), Lopez-Romero et al. (2003) and Fujii and Suzuki (2006). However, at the late stage of the lactation permanent environmental variances were observed in a declined shape. Changes in the permanent environmental variances were opposite to those have reported in studies of Olori et al. (1999a); Rekaya et al. (1999), Jaffrezic et al. (2000) and Fujii and Suzuki (2006). Consequently, decreasing permanent environmental variances within total variances might suggest that permanent environmental factors had the lowest influence on test day milk yields at late stages of lactation for the studied data set.

Moreover, in this study estimated heritability values were higher than those reported by Olori et al. (1999a, b), Rekaya (1999), Jaffrezic et al. (2000), Brotherstone et al. (2000) and Lopez-Romero et al. (2003), especially at the late stages of the lactation. The structural RV model also gave overestimated heritability estimates as compared to the heritability estimates obtained from other models with different RV schemes due to less residual and consequently lower total variances.

CONCLUSION

Accounting for heterogeneity of residual variances is vital for the accurate model definition. This heterogeneity can be clarified via polynomial models with several RV schemes. In this research, accounting for the residual variances with different schemes caused the estimates at each stage of lactation to vary. Especially the random regression models with different RV classes were found incapable of describing variances at the beginning and at the end of the lactation. Contrary to the random regression models, the residual variances don't have to be always equal leading to structural RV models.

After whole evaluation of the models, it was concluded that the random regression RV7 and RV10 models gave better fit. The models, which have more RV classes are recommended to define residual variances through lactation because of their better performance. Also, the structural RV approach can be recommended so that the number of estimated parameters will be lower than the model with different RV classes. In addition, when the changes in residual variance are considered to be continuous over time, there will be no need to define arbitrary RV classes in this approach (Jaffrezic et al., 2000).

This will also, avoid the necessity to increase the number of parameters to be estimated when the number of records per individual and hence, residual variance categories is large (Olori et al., 1999b).

Assuming different RV schemes for the residual variances in random regression models might have significant effect on the variance components at any stages of the lactation. Therefore, it should be taken into consideration when modeling the effects and residual variances.