INTRODUCTION

In data analysis of animal breeding programs, random regression models (RRM) have been used for modeling characteristics measured repeatedly during the animal's life and that change gradually and continuously over time. According to ^{RESENDE et al. (2001}), RRM can more realistically express the phenomena associated with longitudinal data than models of repeatability and finite dimension, by specifying the environmental effects related to a particular stage of lactation and enable the estimation of genetic parameters at any point of the lactation curve, even in animals with incomplete lactations, allowing more frequent genetic evaluations.

To adjust the trajectory of milk production over time using RRM, different functions can be employed. Among the non-linear parametric functions, ^{Wilmink (Wilmink, 1987}) and Ali and Schaeffer (^{ALI & SCHAEFFER, 1987}) are the most important. According to ^{BROTHERSTONE et al. (2000}), the functions of Ali and Schaeffer and Wilmink are more appropriate for adjusting the lactation curve of primiparous Holstein cows, compared to orthogonal Legendre polynomials with the same number of parameters.

In random regression models using non-linear parametric functions and Legendre polynomials, the residual variances can be considered homogeneous or heterogeneous throughout lactation. According to ^{JAMROZIK & SCHAEFFER (1997}), taking homogeneous residual variance into consideration can overestimate the additive genetic variances.

The fact that variances are assumed as heterogeneous tends to improve the partition of the total variance between the variances attributed to random effects included in the analysis models. ^{LÓPEZ-ROMERO et al. (2003}) related that the heterogeneity of residual variances is associated with the lactation stage and is larger at the beginning and at the end of lactation, due to a combination of non-specific factors in the model, such as pregnant stage, dry period features, and body condition at calving.

This study aimed to compare the functions of Wilmink and Ali and Schaeffer with Legendre polynomials in random regression models with different residual variance structures, in the estimation of genetic parameters for test-day milk production of Holstein Friesians reared in Rio Grande do Sul.

MATERIALS AND METHODS:

There were analyzed 5,880 records from the test-day milk production from the first lactation of 907 Holstein cows, born during 2003 to 2009, which were daughters of 235 bulls. Data were provided by the Dairy Herd Analysis Service of the University of Passo Fundo.

In preparing the working file, the test-day milk production was grouped into 20 fortnightly classes of lactation, with Class 1 comprising lactations measured between 6 and 20 days, Class 2, between 21 and 35 days, and so on, until class 20 comprising lactation measured between days 291 to 305 of lactation. Number of observations per fortnight ranged from 146 to 387, with the last three fortnights showing the lowest number, 187, 165 and 146 observations, respectively.

Calving cows aged less than 22 and more than 48 months, as well as records of the test-day milk production showing 3.5 standard deviations for more or less relative to the population mean within the fortnight class were eliminated. Contemporaries groups (CG) comprised animals born in the same herd, year, and control month; and groups with less than five animals were eliminated totalizing 481groups.

The test-day milk production was analyzed using an single trait random regression animal model, considering as fixed effects the contemporaries group and the covariate cow age at calving in months (linear and quadratic). In all models, the mean trajectory of the population (fixed curve) was modeled using orthogonal Legendre polynomials of order three. In the matrix form, the model is given by: y = Xb + Za + Wc + e, where: y is the vector of observations, measured on Nd animals; b is the vector of fixed effects; a is the vector of random coefficients for the additive genetic effect; c is the vector of random coefficients for the permanent environment effect; and e is the vector of residual effects; X, Z and W are the incidence matrices for fixed effects, random genetic additive and permanent environment, respectively, for which, it is assumed:

, K_a and K_c, matrices of (co) variances between the coefficients of additive genetic random regression and permanent environment, respectively; A, the matrix of numerators' relationship between animals; I_Nd,the identity matrix of dimension N_d; and R, a diagonal matrix of residual variances.

The covariance function estimated for the additive genetic effects and permanent environment were modeled for the following functions: Wilmink (WF): WF.33 = a₀ + a₁t + a₂ exp^(-0,05t); Ali and Schaeffer (ASF): ASF.55 = a₀ + a₁u + a₂u² + a₃v + a₄v²; and orthogonal polynomials on a Legendre scale, of the same order of WF and ASF, that is, order three (LEG.33) and five (LEG.55), in which: u = t/305, v = ln(305/t), t = lactation fortnight and a_i = regression coefficients.

For all functions, analyses were performed using homogeneous and heterogeneous residual variance with 2, 3, 4, 5, 6, 7, 9, 12 e 20 classes. Then, for each function, the different residual variance structures were compared using the likelihood ratio test, to identify statistical differences (P≥0.05). Thus, for this study, for all functions, the following residual variance structures were used: homogeneous (HO) and heterogeneous with two classes (HE2) = 1^st to 2^nd and from 3^rd to 20^th fortnight; 4 classes (HE4) = 1^st, 2^nd, 3^rd to 17^th and 18^th to 20^th; 6 classes (HE6) = 1^st, 2^nd, 3^rd, 4^th, 5^th to 17^th and 18^th to 20^th; 7 classes (HE7) = 1^st, 2^nd, 3^rd, 4^th, 5^th to 17^th, 18^th to 19^th, and 20^th; and 9 classes (HE9) = 1^st, 2^nd, 3^rd, 4^th, 5^th to 13^th, 14^th, 15^th to 17^th, 18^th to 19^th, and 20^th lactation fortnight.

The classes of residual variances were formed to improve the model fit, considering the changes in the temporary environment variances that occur at different stages of lactation : crescent, peak and decline, as reported by ^{LOPEZ-ROMERO et al. (2003}).

Choosing the best random regression model was based on the following statistical criteria: logarithmic maximum likelihood function (LML); Akaike information criterion (AIC); Bayesian information criterion (BIC); number of parameters; and graphical visualization of the estimates of variances and genetic parameters. However, each statistical criterion pointed different models as the better, so the following index was used: I = |LML|+ AIC + BIC + number of parameters (adapted of ^{Lui et al., 2006}).

Covariance components and genetic parameters were estimated by restricted maximum likelihood (REML) using the statistical software WOMBAT (^{MEYER, 2006}).

RESULTS AND DISCUSSION:

Mean milk production on the test-day was 28.70kg with a standard deviation of 7.80kg and 27.18% coefficient of variation. Milk production in the first fortnight was 25.50kg, increasing to 30.80 in the fourth fortnight, and decreasing gradually in subsequent fortnights up to 25.25kg in the last fortnight.

Genetic analysis using the Ali and Schaeffer function with homogeneous and heterogeneous residual variance were performed, and convergence difficulties were verified using various algorithms. In addition, even in models that converged, the trajectory of estimates of variances and genetic parameters differed completely from the expected, so results where not considered useful and thus not presented in table 1.

Table 1 Models, number of parameters (p), logarithm of maximum likelihood function (LML), Akaike information criteria (AIC), Bayesian information criteria of Schwarz (BIC), and Index. 

*Models: LEG.kakp = Legendre polynomials; WF.kakp = Wilmink Function; ka and kp order of adjustment of the function of additive genetic covariance and of permanent environment, respectively; HO = residual homogeneous variance; HE residual heterogeneous variance; and ^**Index = |LML| + AIC + BIC + number of parameters.

Within each group, the functions of LML, AIC, and BIC diverged on the model choice, the LML always indicate the ones with more number of parameters WF.33_HE9, LEG.33_HE9, and LEG.55_HE9; the BIC chose the ones with less number of parameters WF.33_HO, LEG33.HE2, and LEG.55_HO; and the AIC chose the intermediates WF.33_HE7, LEG.33_HE7 and LEG.55_HE6. Therefore, to determine the best model, the index which gathers the number of parameters information, AIC, BIC, and LML was used.

The best models pointed by the index were WF.33_HE2, to the Wilmink functions group; LEG.33_HE2 to the Legendre Polynomials, order three group and LEG.55_HE4 Legendre Polynomials, order three group. The LEG.55_HE4 resulting in the best index value (Table 1). A similar result was observed by ^{BIGNARDI et al. (2011}) and ^{SANTOS et al. (2014}) for Holstein and Guzerá breeds, respectively.

^{BROTHERSTONE et al. (2000}) observed that the non-linear parametric functions were more appropriate for adjusting the lactation curve of primiparous Holstein cows, compared to orthogonal Legendre polynomials with the same number of parameters. In this study, superiority of the non-linear parametric functions with respect to the Legendre Polynomials was observed only when using the WF, as observed by ^{KHEIRABADI et al. (2014}) in Holstein cattle.

Overall, the WF.33_HE2, LEG.33_HE2, and LEG.55_HE4 models showed little difference in phenotypic variance partition in residual variance for permanent environment and genetic traits (Figure 1). The additive genetic variance estimated by the four models decreased from the first (@ 31.50) until the 17^th fortnight (@12.15), and increased at the end of lactation, where the WF.33_HE2 and LEG.33_HE2 reached values close to 19.00 and the LEG.55_HE4 reached 14.85 (Figure 1A). Similar behavior was observed by ^{BIGNARDI et al. (2011}) and ^{PEREIRA et al. (2010}) for Holstein and Gir breed cows, respectively.

Figure 1 Additive genetic variance (A), permanent environment (B), residual (C) and phenotypic (D) estimated using WF.33_HE2 (), LEG.33_HE2 (○) and LEG.55_HE4 (--). 

The permanent environment variance estimates were higher than those of additive genetic variance after the eighth fortnight for LEG.55_HE4, the eighth to 18^th fortnight for WF.33_HE2, and the eighth to the 16^th fortnight for LEG.33_HE2 (Figure 1B). With respect to residual variance differences between models at the beginning and the end of the lactation, the LEG.55_HE4 resulted in lower estimates in the range of 7.19 in the last 3 fortnights to 17.47 in the second fortnight of lactation (Figure 1C).

The phenotypic variance showed similar behavior to the additive genetic variance, and although the residual variance was lower than the genetic variance and permanent environment throughout lactation, this influenced the behavior of the phenotypic variance at the beginning and at the end of the curve (Figure 1D). According to ^{MEYER (2000}), the differences in the structures of the residual variance have higher reflection than the phenotypic variance in the components attributable to different causes. ^{VAN VLECK & HENDERSON (1961}) related that milk production at the beginning and at the end of lactation are more subject to temporary environmental variation than the production in the middle of lactation, which is more influenced by genetic differences and permanent environment between animals.

Heritability estimates by WF.33_HE2 and LEG.33_HE2 models were similar, ranging from 0.31 to 0.50, with a tendency to overestimate the beginning and the end of lactation (Figure 2), in agreement with those reported by ^{BREDA et al. (2010}) and ^{SANTOS et al. (2014}). The LEG.55_HE4 diverged from these, in early lactation with higher estimates and from the 16^th fortnight with lower estimates. The higher estimates observed in this study are probably due to the fact that the herds studied were not being selected. According to ^{BULMER (1971}), although individually the genes do not have an important effect on the infinitesimal model, in a selection program, the selection of parents that mate randomly generates an imbalance caused by the covariance between genes in different loci of the same gamete (gametic disequilibrium phase). This covariance being under negative directional selection causes a decrease in additive genetic variance, and thereby in the value of heritability.

Figure 2 Heritability estimated using WF.33_HE2 (), LEG.33_HE2 (○) and LEG.55_HE4 (--). 

The genetic and permanent environment correlations estimated by LEG.55_HE4 were not always smaller with the increase in distance among controls, as expected. The WF.33_HE2 and LEG.33_HE2 presented estimations of genetic and permanent environment correlations, which decreased with an increase in distance among controls; however, negative correlations between the beginning and the end of lactation were verified only by genetic correlations (Figure 3). This behavior is probably associated with the difficulty to model milk production at the beginning of lactation (^{BIGNARDI et al., 2011}), which can be explained by a reduced number of observations and by the proximity to post-partum and pregnancy beginning stress period (^{BIGNARDI et al., 2009}). Negative correlations estimated by random regression models using different functions were shown by ^{COBUCI et al. (2005}); ^{COSTA et al. (2008}) and ^{BIGNARDI et al. (2011}) in Holstein Friesian cows and by ^{PEREIRA et al. (2010}) in the Gir breed.

Figure 3 Estimates of genetic correlations (A) and permanent environment (B) among fortnights of control using WF.33_ HE2, LEG.33_HE2 and LEG.55_HE4. 

The LEG55_HE4, among the three better models shown by the index, has the highest number of parameters (14 vs 34) and resulted in lower estimations of residual variance at the beginning and at the end of lactation, but overestimated the heritability in the first fortnight and presented a greater difficulty to model genetic and permanent environment correlations among controls.

CONCLUSION:

Random regression models that used the Ali and Schaeffer function are not recommended for this data base. Random regression models that used the Wilmink and Legendre polynomials functions with two classes of residual variance appropriately described the genetic variation during the lactation period of Holstein Friesians reared in Rio Grande do Sul.