Estimating variance components for test day milk records by restricted maximum likelihood with a random regression animal model

Livestock Production Science 61 (1999) 53–63

Estimating variance components for test day milk records by restricted maximum likelihood with a random regression animal model a, a b a V.E. Olori *, W.G. Hill , B.J. McGuirk , S. Brotherstone a

Institute of Cell, Animal and Population Biology, University of Edinburgh, West Mains Road, Edinburgh, Scotland EH9 3 JT, UK b Genus Ltd., Vallum Farm, Military Road, Stamfordham, Newcastle-Upon-Tyne, NE18 0 LL, UK Received 20 August 1998; received in revised form 11 January 1999; accepted 22 January 1999

Abstract Variance components for weekly averages of daily milk yield were estimated by restricted maximum likelihood (REML) using a random regression animal model. Additive genetic (s 2a ) and permanent environmental (s 2pe ) covariances were modelled with orthogonal polynomial regressions of varying order, while residual variance (s e2 ) was assumed to be constant for yields in all or some weeks of lactation. The data comprised records of 488 first lactation Holstein–Friesian cows in one herd. The results indicate that the log likelihood increased as the order of polynomial regression in both the fixed and random part of the model increased from 3 to 5. However, only the first three eigenvalues of the additive covariance coefficient 2 matrix were greater than zero in all models. Estimates of s 2a and s pe did not show any trend due to the increase in the order of the covariance function. Additive genetic variance declined from about 9 kg 2 in week 4 to 6 kg 2 in week 10 and increased linearly afterwards to peak at about 16 kg 2 in week 35. Permanent environmental variance was relatively constant at about 12 kg 2 in the first 35 weeks of lactation for most models. Estimates of residual variance (s 2e ) in all stages of lactation declined as the order of fit of the other variance components increased and depended on the assumption about variation in measurement error across lactation. When assumed constant throughout lactation, s 2e was estimated as 3.2, 2.8 and 2.6 kg 2 for the quadratic, cubic and quartic models, respectively.  1999 Elsevier Science B.V. All rights reserved. Keywords: Covariance function; Dairy cattle; Random regression; REML; Test day records

1. Introduction Genetic evaluation of dairy cattle in most countries and the ranking of animals for selection is based almost entirely on best linear unbiased prediction (BLUP) of genetic merit (INTERBULL, 1992). For production traits, these are obtained by analysing *Corresponding author.

305-day or equivalent cumulative yield records predicted from a few test day yields, even though 305-day yields predicted from monthly test day records are biased (McDaniel, 1969; Anderson et al., 1989). The error of genetic evaluation may further increase, if 305-day yields are obtained by projecting partial lactations with factors (Wiggans and Van Vleck, 1979; Pander and Hill, 1993), which assume a constant shape of the lactation curve for all cows

0301-6226 / 99 / $ – see front matter  1999 Elsevier Science B.V. All rights reserved. PII: S0301-6226( 99 )00052-4

54

V.E. Olori et al. / Livestock Production Science 61 (1999) 53 – 63

contrary to reality (Shanks et al., 1981; Guo and Swalve, 1997; Olori et al., 1999). Genetic evaluation based directly on test day records can overcome the need to predict 305-day yields or project incomplete lactations (Swalve, 1995b, 1998). It will also enable us to avoid the use of expansion factors and the error in the estimation of breeding values that may occur due to the disproportionate weighting of 305-day yields projected from unequal number of test day records (VanRaden et al., 1991). Test day models allow the effect of short-term environmental factors specific to individual yields, such as the effect of gestation stage (Olori et al., 1997), to be accounted for directly by including them in the model. They can also facilitate a cheaper and more flexible recording scheme and perhaps allow the inclusion of owner sampled records in genetic evaluation. It is therefore important to explore the potential of any statistical and computing technique which allows a direct and more efficient utilisation, for genetic evaluation, of all available records. Analysis of test day records involves a range of assumptions about the structure of the covariances amongst records. At one extreme, test day records are analysed as repeated measures of the same trait with a constant variance (Ptak and Schaeffer, 1993). This assumes a genetic correlation of unity between yields at different stages of lactation, contrary to findings (Meyer et al., 1989; Pander et al., 1992; Swalve, 1995a). Conversely, test day records could be treated as different traits in a fully parameterised multivariate model with unstructured covariances (Meyer et al., 1989; Pander et al., 1992, 1993); but genetic analysis with this model is limited by the number of ‘traits’ or repeated measures per individual (Thompson and Hill, 1990; Wade et al., 1990; Wiggans and Goddard, 1996) and is computationally expensive. Jamrozik and Schaeffer (1997) proposed a random regression model for the analysis of test day records. This is a practical intermediate parameterisation approach, which assumes a structure for the covariance of repeated records (Kirkpatrick et al., 1994; Meyer, 1998a,b). The covariance structure can be characterised by a covariance function, estimated by fitting a set of orthogonal polynomials or other defined co-variables as random regressions on time in the analysis of repeated records (Meyer and Hill,

1997; Meyer, 1998a,b). This can allow, for example, the estimation of the covariance between daily milk yields within a given time scale, including days that were not sampled, thereby enabling a more efficient utilisation of test day records and providing the potential for a more accurate evaluation of individuals animals. Previous heritability estimates for test day milk yield obtained by random regression (Jamrozik and Schaeffer, 1997; Kettunen et al., 1997, 1998) were higher than expected, especially for yields at the beginning and end of lactation. Also the correlation between yields at the extremes of the lactation trajectory were negative contrary to the trend from multivariate analyses of test day yields (Meyer et al., 1989; Pander et al., 1992; Swalve, 1995a). These preliminary results (Jamrozik and Schaeffer, 1997; Kettunen et al., 1997, 1998) suggest a need for further studies on the derivation of covariance functions for dairy cattle production data. The aim of this study was to estimate coefficients of the covariance function and hence the variance components for weekly averages of daily milk yield by restricted maximum likelihood (REML), determine the effect of modelling the genetic and permanent environmental covariances with orthogonal polynomial regressions of varying order and compare the effect of assuming a constant residual variance with a value that may vary for different stages of lactation whilst modelling the other components with covariance functions. The analysis utilised a rather small but densely recorded data set.

2. Materials and method Records of daily milk production were obtained from a single multiple ovulation and embryo transfer (MOET) herd located in the North East of England. These comprised first lactation records of 488 Holstein–Friesian cows calving between 1990 and 1994. The animals were generally reared indoors and fed a ration made up of 60%whole crop silage and 40% concentrate (Strathie and McGuirk, 1995). The cows were daughters of 49 sires and 94 dams while the pedigree also included 52 grand sires and 58 grand dams with progeny records. Daily milk production from each cow was summarised weekly to obtain the

V.E. Olori et al. / Livestock Production Science 61 (1999) 53 – 63

mean daily milk yield in each week of lactation. A maximum of 37 records per cow, representing the mean yields between lactation weeks 4 and 40 (inclusive) were included in the analysis. Variance components were estimated by REML using a random regression animal model with the DXMRR statistical package (Meyer, 1998c). Detailed derivation of the log likelihood and REML estimation of covariance functions has been presented (Meyer and Hill, 1997; Meyer, 1997, 1998a). In this study, covariance functions for the random additive and permanent environmental components of variance were modelled with orthogonal polynomial regressions (on the Legendre scale) of order 2–4. The general model (Meyer, 1997) can be specified as; kA 21

y ij 5 F 1

k R 21

O a f (t * ) 1 O l im

m 50

m

ij

im

fm (t ij* ) 1 eij

m50

where y ij is the jth record of the ith animal; F represents the fixed effects of: 1. calving age (in months, 11 classes ranging from 23 to 33 months of age); 2. gestation stage (nine classes where class 1 represent records obtained when the animal was not pregnant, and classes 2–9 represent records obtained when the animal was 1–8 months pregnant, respectively (see Olori et al. (1997) for details); 3. production season (29 year by season of production classes with each season consisting of production over 2–3 calendar months. Data were obtained over a 4-year period); 4. lactation stage (modelled with polynomial regressions of the same order as the random regressions sub-model).

aim and lim represent the mth additive genetic and permanent environmental random regression coefficients for the ith animal, respectively; kA and k R represent the order of the fitted polynomial regression for the additive genetic and permanent environmental components of variance, respectively (always equal); t ij* is the jth standardised lactation week of the ith animal. This ranged between 21 and 1 on the Legendre scale; and fm (t ij * ) is the mth polynomial evaluated for t ij * .

55

Additive genetic, permanent environmental (co)variances and the fixed regression were modelled with the same order polynomial regression in all cases. The number of measurement error (ME) categories specified varied. For ME51, it was assumed that the residual variance was constant for all stages of lactation, while for ME537, it was allowed to differ each week of lactation. For ME54, the residual variance was assumed to be constant for yields within, but different between, the following groups of lactation week; 4–5, 6–14, 15–33 and 34–40, chosen to reflect different phases of the lactation. For ME510, the residual variance was assumed to be different weeks 4–9 of lactation, and assumed constant within but different between yields in the following groups of weeks; 10–14, 15–25, 26–33 and 34–40. Significant differences in the fit of models with the same order of fixed regressions, was tested using a x 2 -test of the likelihood (Kirkpatrick et al., 1990). The first three eigenfunctions of the covariance functions were derived and plotted and each eigenvalue was expressed as a percentage of the sum of all to determine its importance.

3. Results

3.1. Log likelihood The maximum of the log-likelihood function (ln L) increased as the order of the polynomial function used to model the additive genetic and permanent environmental covariances increased from 3 to 5. The changes in the likelihood peak as the order of fit increases are presented in Table 1 as deviations from the maximum of the quadratic regression model with ME51. In terms of these likelihoods, the model with the best fit was the quartic polynomial regression with 10 ME classes. The worst fit was obtained when ME was assumed constant for all stages of lactation and covariances were modelled with a quadratic function. The results in Table 1 indicate significant increases in the log likelihood function due to increases in the number of ME within each polynomial regression model. For example, ln L increased by 296, 162 and 130 when ME was increased from 1 to 10 for the quadratic, cubic and quartic models,

V.E. Olori et al. / Livestock Production Science 61 (1999) 53 – 63

56

Table 1 Maximum log likelihood and the change in log likelihood for different models Regression model

Number of ME classes

Number of parameters

Change in log likelihood

Quadratic Quadratic Quadratic Quadratic Cubic Cubic Cubic Quartic Quartic Quartic

1 4 10 37 1 4 10 1 4 10

13 16 33 49 21 24 30 31 34 40

– 251 296 440 512 637 674 808 915 938

respectively. Allowing residual variance to vary in each lactation week (ME537), instead of assuming it was constant throughout lactation (ME51), resulted in a change in the maximum likelihood of

about 440 for the quadratic model. The difference in the number of parameters estimated by these two models was 36.

3.2. Eigenvalues and eigenfunctions The first three eigenvalues of the coefficient matrix of the additive genetic covariance function accounted for about 99% of the sum of all eigenvalues for each model (Table 2). The leading eigenvalue alone accounted for about 93% of the sum of all eigenvalues on average. Eigenfunctions corresponding to the first three eigenvalues for the cubic model are plotted in Fig. 1. The leading eigenfunction was positive and similar in value (0.4–0.9) at all stages of lactation implying a positive genetic correlation across all stages. The second eigenvalue accounted for 6.1, 6.0 and 5.5% of the sum of the eigenvalues for the quadratic,

Table 2 First three eigenvalues and their proportions for the additive genetic covariance function Polynomial regression Quadratic (k53) Cubic (k54) Quartic (k55)

Eigenvalues

Proportion of total (3100)

First

Second

Third

Total

First

Second

Third

19.66 19.80 20.56

1.96 1.27 1.21

0.13 0.25 0.24

21.75 21.32 22.01

90.4 92.9 93.4

9.0 5.9 5.5

0.6 1.2 1.1

Fig. 1. Eigenfunctions of the additive genetic covariance function coefficient matrix (k54, ME510).

V.E. Olori et al. / Livestock Production Science 61 (1999) 53 – 63

cubic and quartic models, respectively. The corresponding eigenfunction was negative in early lactation and positive in late lactation for all models. This change in sign suggests association with a factor, or factors, having contrasting effects on milk yield in early and late stages of lactation. The third, fourth and fifth eigenvalues (where applicable) made up about 1% of the sum of the eigenvalues for each model. This suggests that the potential for overall change in the pattern of milk yield over lactation by selection based on factors represented by the third, fourth and fifth eigenfunctions is small as they account for a negligible proportion of the variation in additive genetic variance across lactation stage.

3.3. Variance components Table 3 shows a summary of the additive genetic 2 (s 2a ), permanent environment (s pe ) and phenotypic 2 (s p ) variances estimated by the different models. Fig. 2 shows the different components estimated with the cubic model (ME510) on a continuous scale. Estimates of the different components of variance depended on the lactation stage and the order (k) of the polynomial fitted in the random regression sub-model. Genetic variances of milk yield in week 4 (the earliest week in the data),

57

averaged over models with ME51, 4 and 10, were 6.8, 9.2 and 10.7 kg 2 for quadratic (k53) cubic (k54) and quartic (k55) models, respectively. The estimates declined at different rates between weeks 4 and 8 to more similar figures suggesting that the models were less robust in modelling the variance of yields in early lactation. This may be due to the more profound effect of constraining residual variance in the early stages of lactation and or the mathematical properties of the polynomial regressions. Table 3 indicates for the cubic model, little difference between ME classes in the mean estimates of s 2a and s 2pe at different stages of lactation. However, an examination of the estimates for individual weeks showed local variation between models, which differed only in the number of ME classes.

3.4. Measurement errors and residual variance The estimates of residual variance of yield at different stages of lactation depended on the lactation weeks assumed to have a constant error variance in the model fitted. When measurement error was assumed to be constant for all stages of lactation (ME51), estimates of s 2e were 3.2, 2.8 and 2.6 kg 2 for quadratic, cubic and quartic models, respectively, implying that residual variance declined as the order

Table 3 Mean (6S.D.) estimate of variance components for daily milk yield at different stages of lactation Variance components

Lactation stage (weeks)

Regression sub-models in fixed and random part Quadratic (k53) ME54

Cubic (k54) ME51

ME54

ME510

Quartic (k55) ME54

4–10 11–19 20–28 29–36 37–40

6.760.85 7.960.64 10.760.96 13.860.95 16.360.57

7.261.04 7.260.78 10.961.39 14.860.77 15.560.16

7.261.03 7.360.79 11.061.38 14.860.82 15.860.09

7.461.13 7.360.75 10.961.38 14.860.84 15.860.07

7.461.68 7.461.14 11.761.23 14.560.57 16.360.54

Permanent environment (s 2c )

4–10 11–19 20–28 29–36 37–40

11.760.93 10.860.22 11.660.20 12.460.48 14.760.89

11.460.38 12.260.19 12.260.06 12.760.22 14.461.01

11.560.46 12.060.16 12.060.06 12.660.31 14.661.06

11.760.57 12.160.17 12.060.06 12.560.30 14.561.07

11.360.48 12.360.07 11.960.19 12.160.43 14.561.27

Phenotypic (s 2p )

4–10 11–19 20–28 29–36 37–40

22.863.67 21.761.15 25.661.16 29.160.94 33.161.47

21.461.34 22.260.93 25.961.44 30.260.98 32.760.85

22.162.79 22.161.15 26.061.42 30.060.71 32.360.98

22.463.61 22.261.12 26.061.32 29.960.75 32.361.00

21.361.74 22.561.33 26.561.05 29.160.56 32.461.81

Genetic (s 2a )

V.E. Olori et al. / Livestock Production Science 61 (1999) 53 – 63

58

Fig. 2. Additive genetic (Va), permanent environmental (Vc) and residual (Ve) components of variance as well as the phenotypic variance of weekly averages of daily milk yield (k54, ME510).

of the polynomial sub-model increased. Table 4 shows a summary s 2e estimated by the different models when ME was assumed to vary for each week (ME537) or between 4 (ME54) and 10 (ME510) groups of lactation weeks, respectively.

Generally, s 2e was highest at the beginning of lactation and relatively constant in mid lactation. When residual variance was assumed to be constant for milk yield in lactation weeks 4 and 5 (ME54), the estimate was 8.6 kg 2 for the quadratic model.

Table 4 Residual variances (kg 2 ) estimated by random regression models with different order of fit by measurement error classes Lactation stage (weeks)

Variable errors a (ME537) Quadratic sub-model

Constant errors within 10 groups of lactation weeks (ME510) Quadratic sub-model

Cubic sub-model

Quartic sub-model

Quadratic sub-model

Cubic sub-model

Quartic sub-model

4 5 6 7 8 9 10–14 15–25 26–33 34–40

11.71 6.35 3.72 1.83 2.45 2.07 2.90 a 3.27 a 3.22 a 2.12 a

11.64 6.28 3.68 1.83 2.51 2.17 2.88 3.33 3.24 2.08

6.69 3.78 2.68 1.66 2.76 2.32 2.62 3.13 2.94 1.91

3.97 2.27 3.01 1.88 2.65 2.06 2.50 2.99 2.90 1.74

8.56 8.56 2.76 2.76 2.76 2.76 2.76 3.28 3.28 2.08

5.42 5.42 2.54 2.54 2.54 2.54 2.54 3.06 3.06 1.91

2.97 2.97 2.49 2.49 2.49 2.49 2.49 2.96 2.96 1.74

a

Means with S.D. about 0.8 for each class of weeks.

Constant errors within four groups of lactation weeks (ME54)

V.E. Olori et al. / Livestock Production Science 61 (1999) 53 – 63

However when s 2e was assumed to vary between these two weeks (ME510), the estimates were 11.6 and 6.3 kg 2 for yields in weeks 4 and 5, respectively. This shows that the estimate of s 2e depends on which stages of lactation are a priori assigned to the same measurement error class and hence constrained to have the same estimate of residual variance. Results not presented showed that accurate identification of lactation stages with a similar error variance was more important than the number of ME classes specified in the RR model.

3.5. Heritability Fig. 3 shows a plot of the heritability of milk yield at different stages of lactation estimated by random regression with different covariance functions (models) and Table 5 gives some estimates from the models with 10 ME classes. The cubic model gave moderate estimates of the heritability in most stages of lactation, ranging from 0.31 in lactation week 4 to a maximum of 0.51 in lactation week 35. Corresponding estimates were 0.22 and 0.50 for the

59

quadratic model and 0.41 and 0.52 for the quartic model. The h 2 estimates were highest in week 35 for all models. The quartic model yielded unexpectedly high h 2 (0.41) for yield in lactation week 4 and a large drop at week 10 (h 2 50.29). Estimates obtained by the quadratic model followed a different pattern being lowest in week 4 and highest in week 35.

3.6. Genetic and phenotypic correlation Genetic (r g ) and phenotypic (r p ) correlations of milk yield in selected weeks of lactation for the different models with ME510 are presented in Table 5. For any pair of lactation weeks, r g was consistently higher than r p and both declined as the interval between yields increased. For example, the genetic correlation between yield in week 4 and yields in weeks 8–36 ranged from 0.95 to 0.63 and the phenotypic correlations from 0.77 to 0.49 with the cubic (k54) random regression model. There were only slight differences between the models in both r g and r p and none of the models gave consistently higher estimates across lactation. Ge-

Fig. 3. Estimates of heritability of average daily milk yield by different models.

V.E. Olori et al. / Livestock Production Science 61 (1999) 53 – 63

60

Table 5 Heritability, maternal effect (c 2 ), genetic and phenotypic correlation amongst milk yield in selected weeks of lactation Lactation week (I)

Lactation week (j)

4 10 20 30 40

4 10 20 30 40

c2

Heritability Quadratic (k53)

Cubic (k54)

Quartic (k55)

Quadratic (k53)

Cubic (k54)

Quartic (k55)

0.22 0.34 0.39 0.46 0.49

0.31 0.31 0.37 0.48 0.47

0.42 0.29 0.40 0.49 0.49

0.42 0.52 0.47 0.43 0.45

0.43 0.56 0.50 0.42 0.47

0.43 0.58 0.48 0.41 0.46

Genetic correlations 4 4 4 4 4 4 4 10 10 10 20 20 30

10 15 20 25 30 35 40 20 30 40 30 40 40

0.94 0.85 0.76 0.69 0.64 0.61 0.58 0.93 0.84 0.77 0.98 0.92 0.98

0.89 0.78 0.73 0.72 0.69 0.65 0.54 0.94 0.86 0.73 0.97 0.90 0.97

netic correlation between yields in adjacent weeks were close to unity especially in mid lactation.

4. Discussion One advantage of a random regression model over a multiple trait model is the reduction in the number of parameters that need to be estimated (Meyer and Hill, 1997). In RR models based on REML, the number of parameters depends on the number of (co)variance components to be modelled and the order of the covariance function for each modelled variance component (Meyer, 1998a,b). For nested models with the same order of fixed polynomials but different numbers of random components, the likelihoods can be tested statistically using a x 2 -test (Kirkpatrick et al., 1990). Models with ME54, 10 and 37 can therefore be compared statistically with the model where ME51. Although theoretically, we expect the fit of a covariance function to improve with the order, the number of parameters and the accuracy of the models in partitioning total variance

Phenotypic correlations 0.90 0.78 0.75 0.72 0.70 0.65 0.60 0.94 0.86 0.78 0.97 0.93 0.99

0.69 0.63 0.52 0.50 0.46 0.44 0.40 0.78 0.70 0.59 0.84 0.72 0.84

0.68 0.60 0.53 0.53 0.52 0.50 0.42 0.79 0.67 0.60 0.82 0.72 0.82

0.70 0.63 0.59 0.59 0.55 0.53 0.47 0.78 0.67 0.58 0.82 0.73 0.83

into the various components is important in the choice of a random regression model. The magnitude of the first and second eigenvalues and the associated eigenfunctions in this study agree with previous values reported for test day milk yield in dairy cattle (Kirkpatrick et al., 1994). The results suggest that factors associated with the leading eigenvalues account for over 90% of the genetic variation in milk yield at all stages of lactation. The relatively constant values of the corresponding eigenfunction across lactation imply that milk yield is mainly controlled by genes with similar effects at all stages of lactation. The leading eigenfunction and corresponding eigenvalue can be associated with the constant coefficient of the covariance function and hence the mean additive genetic variance of milk yield at different stages of lactation. Genetic correlations obtained in this study were all positive, even between yields at the extremes of the lactation trajectory. The genetic correlation between yields in weeks 4 and 40 was over 0.50 for all models suggesting that selection for increased milk yield in early lactation will have a positive effect on

V.E. Olori et al. / Livestock Production Science 61 (1999) 53 – 63

yield in late lactation. This can be supported by the absence of trend in the first eigenfunction which captures over 90% of the variation in the additive genetic variance, indicating that selection on factors represented by this first canonical variate or almost equivalently mean yields at each stage, will cause a shift in the entire lactation curve. The change in sign of the second eigenfunction between early and late stages of lactation in this study tends to suggest an underlying factor with opposing effects on yields in early and late stages of lactation. Selection on such a factor, e.g. persistence (Gengler et al., 1995), can be used to change the shape of the lactation curve. Response to selection based on the second eigenvalue will be very small as it accounts for only about 6% of the variation in additive genetic variance. However the magnitude and trend of the first and second canonical variates suggests that improvement in milk production can be achieved by basing selection on mean daily milk yield and persistency. In this study, only the first three eigenvalues of the additive genetic coefficient matrix were greater than zero in all models. This implies that a relatively simple covariance function based on a coefficient matrix of rank 3 may model the additive (co)variances adequately. More sophisticated models could result in over parameterisation which may lead to high estimates of the additive variance at the extremes of the lactation trajectory (Jamrozik and Schaeffer, 1997; Kettunen et al., 1998). Estimates of the additive and permanent environmental components of variance for test day yields in a random regression model depend on the sub-model fitted as the covariance function. Different models, including standard lactation curve models and orthogonal polynomials of varying order, have been used to model the covariances and the fixed effect of lactation stage in RR models (Jamrozik and Schaeffer, 1997; Jamrozik et al., 1997, 1998; Kettunen et al., 1997, 1998). The general trend seem to indicate that the sub-model for one component has an effect on the optimal model and estimated parameters for the other components. This trend was more obvious in the estimation of the residual variance component in this study. Heritability estimates were lower in this study than those reported by Kettunen et al. (1997) and by

61

Jamrozik and Schaeffer (1997) when they assumed that permanent environmental variance was constant at all stages of lactation, but are higher than those reported by Kettunen et al. (1997) when they modelled the permanent environmental variance with a normalised third order (or quadratic) polynomial. When residual variance was assumed constant, its estimate was approximately the mean of the estimates obtained by allowing residual variance to vary. It may therefore be more appropriate to allow residual variances to vary than to constrain it to a single value for yields in some stages of lactation. Where possible, modelling residual variance with a continuous function similar to or less sophisticated than the model for the other variance components may be required for a more accurate partitioning of total variance into its components. This will also avoid the need to increase the number of parameters to be estimated when the number of records per individual and hence residual variance categories is large. The decline in residual variance observed as the order of the regression sub-model increased may be due to the increased ability of the higher order polynomials to model the lactation stage effect especially in the beginning of lactation. As the fit of the model increased, more variation in yield due to lactation stage is accounted for by the fixed regression thereby reducing residual variance. Guo and Swalve (1997) observed a similar increase in fit when they modelled lactation stage effect with several models and Olori et al. (1999) showed that a five parameter lactation curve model fitted yields in early lactation better than three parameter models. In the present study, increased ability to model yields in early lactation may explain the high phenotypic variance of the models with k53 compared to the higher order models in early lactation as shown in Table 3. Estimates of variance components from the analysis (k53) with 37 measurement error classes were similar to the estimates when only 10 measurement error classes were specified because the grouping of lactation weeks into 10 ME classes was based on the pattern of covariances observed when fitting 37. The accurate grouping of weeks with similar errors into the same ME class is more important than the number of groups specified. Usually, we may not

62

V.E. Olori et al. / Livestock Production Science 61 (1999) 53 – 63

have a prior knowledge of the test days over which residual variance is constant supporting the need to model the residual variances with a continuous function. Such a function should be simpler than the one used for the other variance components because covariances would be assumed to be zero (i.e. a line rather than a surface). As the estimates of variance components were obtained from one herd that is not representative of the UK dairy cattle population, they cannot be regarded as UK estimates of genetic parameters for test day milk yield of Holstein–Friesian cattle. This study has, however, pointed out how such analysis can proceed to provide useful results. It has shown that the critical issues in fitting a random regression model include the order of the polynomial used to model the lactation curve at both the fixed and random level, and the grouping of observations into measurement error classes.

Acknowledgements The authors are grateful to Genus Ltd. for providing the data and to the commonwealth Scholarship Commission in the UK and MDC for financial support for V.E.O. and S.B. respectively. We also thank Karin Meyer for help in installing and using the software.

References Anderson, S.M., Mao, I.L., Gill, J.L., 1989. Effect of frequency and spacing of sampling on accuracy and precision of estimating total lactation milk yield and characteristics of the lactation curve. J. Dairy Sci. 72, 2387–2394. Gengler, N., Keown, J.F., Van Vleck, L.D., 1995. Various persistency measures and relationship with with total, partial and peak milk yield. Br. J. Genet. 2, 237–243. Guo, Z., Swalve, H.H., 1997. Comparison of Different Lactation Curve Sub-models in Test Day Models. Paper Presented at the INTERBULL Open Meeting, 28–29 August 1997, Vienna, Austria. INTERBULL, 1992. Sire evaluation procedures for dairy production traits practised in various countries. International Bull Evaluation Bulletin No. 5. Uppsala, Sweden. Jamrozik, J., Schaeffer, L.R., 1997. Estimates of genetic parameters for a test day model with random regression for yield traits of first lactation Holsteins. J. Dairy Sci. 80, 762–770.

Jamrozik, J., Kistemaker, G.J., Dekkers, J.C.M., Schaeffer, L.R., 1997. Comparison of possible covariates for use in a random regression model for analysis of test day yields. J. Dairy Sci. 80, 2550–2556. Jamrozik, J., Schaeffer, L.R., Grignola, F., 1998. Genetic Parameters for Production Traits and Somatic Cell Score of Canadian Holsteins with Multiple Trait Random Regression Model. Proceedings of the 6th World Congress on Genetic Applied to Livestock Production, Armidale, Australia 23, pp. 303–306. Kettunen, A., Mantysaari, E.A., Stranden, I., Poso, J., 1997. Genetic parameters for test day milk yields of Finish Ayshires with random regression model. J. Dairy Sci. 80 (Suppl. 1), 197. Kettunen, A., Mantysaari, E.A., Stranden, I., Poso, J., 1998. Estimation of Genetic Parameters for First Lactation Test Day Production Using Random Regression Models. Proceedings of the 6th World Congress on Genetic Applied to Livestock Production, Armidale, Australia 23, pp. 307–400. Kirkpatrick, M., Lofsvold, D., Bulmer, M., 1990. Analysis of inheritance, selection and evolution of growth trajectories. Genetics 124, 979–993. Kirkpatrick, M., Hill, W.G., Thompson, R., 1994. Estimating covariance structure of traits during growth and aging, illustrated with lactation in dairy cattle. Genetics Res. Cambridge 64, 57–69. McDaniel, B.T., 1969. Accuracy of sampling procedure for estimating lactation yield: a review. J. Dairy Sci. 52, 1742– 1761. Meyer, K., 1997. An ‘average information’ restricted maximum likelihood algorithm for estimating reduced rank genetic covariance matrices of covariance functions for animal models with equal design matrices. Genet. Selection Evolution 29, 97–116. Meyer, K., 1998a. Estimating covariance functions for longitudinal data using a random regression model. Genet. Selection Evolution 30, 221–240. Meyer, K., 1998b. Modelling Repeated Records Covariance Functions and Random Regression Models to Analyse Animal Breeding Data. Proceedings of the 6th World Congress on Genetic Applied to Livestock Production, Armidale, Australia 25, pp. 517–527. Meyer, K., 1998c. ‘DXMRR’ A Program to Estimate Covariance Functions for Longitudinal Data by Restricted Maximum Likelihood. Proceedings of the 6th World Congress on Genetic Applied to Livestock Production, Armidale, Australia 27, pp. 465–466. Meyer, K., Hill, W.G., 1997. Estimation of genetic and phenotypic covariance functions for longitudinal ‘repeated’ records by restricted maximum likelihood. Livest. Prod. Sci. 47, 185–200. Meyer, K., Graser, H.U., Hammond, K., 1989. Estimates of genetic parameters for first lactation test day production of Australian black and white cows. Livest. Prod. Sci. 21, 177– 199. Olori, V.E., Brotherstone, S., Hill, W.G., McGuirk, B.J., 1997. Effect of gestation stage on milk yield and composition in Holstein Friesian dairy cattle. Livest. Prod. Sci. 52, 167–176. Olori, V.E., Brotherstone, S., McGuirk, B.J., Hill, W.G., 1999. Fit of standard models of the lactation curve to weekly records of

V.E. Olori et al. / Livestock Production Science 61 (1999) 53 – 63 milk production of cows in a single herd. Livest. Prod. Sci. 58, 55–63. Pander, B.L., Hill, W.G., 1993. Genetic evaluation of lactation yield from test day records on incomplete lactations. Livest. Prod. Sci. 37, 23–36. Pander, B.L., Hill, W.G., Thompson, R., 1992. Genetic parameters for test day records of British Holstein–Freisian heifers. Anim. Prod. 55, 11–21. Pander, B.L., Thompson, R., Hill, W.G., 1993. Phenotypic correlation among daily records of milk yields. Indian J. Anim. Sci. 63, 1282–1286. Ptak, E., Schaeffer, L.R., 1993. Use of test day yields for genetic evaluation of dairy sires and cows. Livest. Prod. Sci. 34, 23–34. Shanks, R.D., Berger, P.J., Freeman, A.E., Dickinson, F.N., 1981. Genetic aspects of lactation curves. J. Dairy Sci. 64, 1852– 1860. Strathie, R.J., McGuirk, B.J., 1995. Developments with the Genus MOET dairy breeding scheme. Brit. Cattle Breed. Club Dig. 50, 9–15. Swalve, H.H., 1995a. The effect of test day models on estimation of genetic parameters and breeding values for dairy yield traits. J. Dairy Sci. 78, 929–938. Swalve, H.H., 1995b. Test day models in the analysis of dairy

63

production data: a review. Arch. Tierzucht Dummerstorf 38, 591–612. Swalve, H.H., 1998. Use of Test Day Records for Genetic Evaluation. Proceedings of the 6th World Congress on Genetic Applied to Livestock Production, Armidale, Australia 23, pp. 295–301. Thompson, R., Hill, W.G., 1990. Univariate REML Analysis for Multivariate Data with the Animal Model. Proceedings of the 4th World Congress on Genetic Applied to Livestock Production, Edinburgh, UK 12, pp. 484–487. VanRaden, P.M., Wiggans, G.R., Ernst, C.A., 1991. Expansion of projected lactation yields to stabilise genetic variance. J. Dairy Sci. 74, 4344–4349. Wade, K.M., Quaas, R.L., Van Vleck, L.D. 1990. Mixed Linear Models with Autoregressive Error Structure. Proceedings of the 6th World Congress on Genetic Applied to Livestock Production, Edinburgh, UK 13, pp. 508–511. Wiggans, G.R., Goddard, M.E., 1996. A test day model for genetic evaluation of yield traits: Possible benefits and an approach for implementation. J. Dairy Sci. 79 (Suppl. 1), 144. Wiggans, G.R., Van Vleck, L.D., 1979. Extending partial lactation milk and fat records with a function of last-sample production. J. Dairy Sci. 62, 316–325.