- Email: [email protected]
Genetic Analysis of Dystocia and Calf Mortality in Israeli-Holsteins by Threshold and Linear Models
G e n e t i c A n a l y s i s of D y s t o c i a and Calf M o r t a l i t y in Israeli-Holsteins b y T h r e s h o l d and Linear M o d e l s J. I. WELLER, 1 I. MISZTAL, and D. GIANOLA Department of Animal Sciences University of Illinois Urbana 61801
AB ST RACT
INTRODUCTION
Calvings of 106,751 Israeli Holstein heifers were analyzed for dystocia and calf mortality, scored dichotomously, and a composite trait, scored trichotomously. Dystocia was also studied with 146,973 second and third parity records. Models fitted included herd-year-season, sex of calf, calving age, calving month, sire of cow, sire of calf, and groups of sire of cow and of calf. Herd-year-season, sire of cow and calf, and residuals were random with diagonal variance-covariance matrices. Herd-year-season variance was assumed to be 10% of the residual component. Other variance components were estimated b y REML for linear models and by the counterpart of REML for threshold models. Heritability estimates were two to five times greater in threshold than in linear models, but correlations between corresponding sire evaluations were all >.9. Linear model sire evaluations were skewed positively, whereas threshold model evaluations had symmetrical distributions. Heritability for dystocia was greater in first than in later parities. Correlations between first and later parity sire evaluations were <.5. Thus, the genetic control of dystocia seems to be different for heifers and cows. Correlations between sire of cow and calf evaluations were <.3. Correlations between dystocia and calf mortality evaluations were about .7.
Dystocia (DC) and calf mortality (CM) are of major economic importance to dairy farmers. These traits have low heritabilities, high genetic correlation, and significant residual correlation (1, 2, 13, 17, 18, 19, 20, 21). Many studies have found the effect of sire of calf to be larger than that of the sire of the cow (1, 13, 17, 18, 19, 21). Although both traits are scored categorically, most studies have used standard mixed model methodology, even though this assumes the data to have a continuous distribution. Theoretical problems involved in this type of analysis have been described at length (6, 7). Several studies have suggested that the threshold model (TM), which assumes the existence of an underlying normal variable, is theoretically appropriate and computationally feasible for genetic analysis of categorical traits (3, 4, 6, 7, 8, 9, 10, 11, 14, 15, 16). However, only a few studies have actually applied this model to analysis of field data (4, 11, 14). Because TM equations are nonlinear and involve normal probability functions, computational complexity and computing resources required are greater than in a linear model (LM) analysis. Furthermore, TM genetic evaluations also require accurate estimates of variance components in the underlying scale. Although methods have been developed that are similar to REML, computations are more complex than those in standard REML algorithms (9, 10). Limitations in computing resources may pose a serious problem for sire evaluation based on large field data sets. However, recent developments in computers and programming techniques have made large-scale TM analyses feasible (16). The objectives of this study were to compare TM and LM analyses of calving traits using a large field data set o f Israeli Holsteins, to estimate genetic parameters for these traits and to find the model of choice for routine genetic evaluation.
Received January 7, 1988. Accepted May 9, 1988. 1Institute of Animal Sciences, Agricultural Research Organization, The Volcani Center, Bet Dagan, Israel 50250. 1988 J Dairy Sci 71:2491-2501
2491
2492
WELLER ET AL. MATERIALS AND METHODS
Data were 347,691 first through third parity calvings of milk-recorded Israeli Holsteins, calving between January 1978 and June 1985. Both DC and CM were scored as dichotomous traits. Dystocia was scored as 0 if the calving was normal or as 1 if it was difficult; CM was scored as 0 for a live birth or as 1 if the calf died within 48 h of birth. Records were deleted from the analysis for the following reasons: 1) sire of cow or of calf was unknown, 2) multiple births, 3) either CM or DC was undefined, 4) the calf was abnormal or died due to causes not related to calving, 5) abortions, 6) calvings of daughters of sires with less than 20 daughter calvings, and 7) calvings in which the sire of the calf sired less than 20 calves. Dystocia was recorded in all but two calvings. However, 2.1% of the records were deleted because either calf fate was unknown, the calf was abnormal, it died due to causes unrelated to calving, or it was aborted. There are no natural service matings in Israel, and about 1000 inseminations are performed for each young sire tested. Therefore restrictions 6 and 7 eliminated 1274 observations, consisting of recording mistakes, progeny of a few matings from imported semen, and progeny of sires tested at the boundaries of the time period considered. In total, 27% of the records (93,967) were deleted, mostly because the sire of the calf was unknown. All analyses were by both TM and LM. Dystocia was analyzed using all parities together and first and later parities separately. Based on these results, it was decided to analyze CM using only first parity calvings. Because these traits have a high genetic correlation, a composite trait (CT) was also analyzed. In the LM analysis, CT was the sum of the DC and CM scores; thus, each calving received a score of either 0, 1, or 2. In the TM analysis, it was assumed that these tbree classes were ordered so that the model included two thresholds. Thus five analyses were run for each of the two statistical methods: 1) DC, all parities; 2) DC, first parity; 3) DC, later parities; 4) CM, first parity; and 5) CT, first parity. The following model was used for analysis:
Journal of Dairy Science Vol. 71, No. 9, 1988
Yijklmnopq = HYSi + GSj + SIREjk + GSC I + SCIm + Sn + A o + Mp + eijklmnopq [1] where: Yijklmnopq = record on cow ijklmnopq; HYS i = random effect of herd-yearseason i; GSj = fixed effect of group j of sires of cows; SIREjk = random effect of sire of cow k in group j ; GSC 1 = fixed effect of group 1 of sires of calves; SClm = random effect of sire of calf m in group l; S n = fixed effect of sex n; A o = fixed effect of calving age o; Mp = fixed effect of calving month p; and eijklrnpopq = random residual. Calving seasons were October to March, and April to September. Sires were grouped by year of birth, with all sires born prior to 1972 in group 1. Groups 2 to 7 included all sires born during each subsequent 2-yr interval. Twelve classes were defined for calving month. Because HYS was random, there was no "confounding" between HYS and calving month. The first 10 calving ages were for primiparous cows representing calving ages from 21 to 30 mo. The 1st also included cows calving prior to 21 too, and the lOth included cows that calved after 30 mo. In the all parities analysis, second and third parity cows were included in classes 11 and 12, respectively, regardless of calving age. In the latter parity analysis, only two classes were defined for this effect. Analysis was facilitated by assuming that all random effects were mutually uneorrelated with diagonal dispersion matrices. It was further assumed that the variance of the HYS effect was known and equal to 10% of the residual variance. This value is similar to estimates calculated from US data (4). Using these
GENETIC ANALYSIS OF CALVING TRAITS assumptions, HYS effects, which had many levels in all analyses, could be absorbed. The assumption of diagonal dispersion matrices was inaccurate for the SIRE and SC effects, because sires were related. The assumption of zero correlation between effects was also an approximation, because several studies have found nonnegligible correlations between SIRE and SC effects (13, 17, 18, 21). A program based on principles described by Misztal et al. (16) and adapted to the CRAY XMP-48 supercomputer, was used for analysis. The SIRE, SC, and residual variance components were estimated by REML for LM model analyses and by the counterpart of REML for TM (9). After absorption of HYS effects, coefficient matrices were inverted at each round of iteration for both models. In TM analyses, interation involved Fisher's scoring and a modification of the expectation-maximization (EM) algorithm for variance components (16). Two to three EM rounds were performed per round of Fisher scoring iteration. At least six rounds of REML iteration in LM and six rounds of Fisher scoring iteration in TM were completed for each analysis. Iteration stopped when the change in all variance components was less than 1% of the values in the previous round. Sire evaluations were computed as the sum of the sire and sire group effects. Heritability (h 2) as a trait of the cow was calculated as: h 2 = 4 VAR(SIRE)/ [VAR(SIRE) + VAR(SC) + VAR(HYS) + VAR(E)]
[2]
where VAR(SIRE), VAR(SC), VAR(HYS), and VAR(E) are the SIRE, SC, HYS, and residual variance components, respectively. Heritability as a trait of the calf was computed in a similar manner except that VAR(SC) replaced VAR(SIRE) in the numerator. Threshold model heritabilities for dichotomous traits were also obtained dividing LM estimates by: z2/[p(1 - p)]
[3]
where z is the ordinate of the standard normal density function corresponding to p, the observed incidence in the first category (3). For CT, a trait with two thresholds, a formula of Gianola (7) was used.
2493
The unit of measurement of estimates and predictions in TM is the residual standard deviation. To facilitate comparison between LM and TM estimates, the TM solutions were multiplied by the square root of the corresponding LM residual variance component so that the units of the two models would be equivalent. Results obtained with the two different methods of analysis were compared in terms of ratios of variance components (heritability), the variance and skewness of the estimates and evaluations, and correlations between the solutions or evaluations. Although all sires included in the full data set had 20 or more calving records, some sires in the first and later parity subsets had less than 20 records. Therefore, correlations between sire evaluations were also computed for sires with at least 20 records in both analyses, and for sires with repeatability based on LM greater than .6. Repeatability was calculated as 1 - PEV/VC, where PEV is the prediction error variance of the SIRE or SC effect, and VC is the corresponding variance component. The PEV was computed as the diagonal d e m e n t of the inverse of the coefficient matrix at the final round of iteration times the residual variance. R ESU LTS
Basic statistics for the edited data set are in Table 1 for all parities combined and for first and later parities separately. Incidences of DC and CM for all parities were 5.1 and 5.4%, respectively; the incidence of either DC or CM in a single calving was 7.9%. Frequencies were higher for heifers, which is consistent with previous studies (1, 4, 13, 17, 20, 21). The distribution of records by level of fixed effects is presented in Table 2 for all parities combined, and for first and later parities separately. The ratio of male to female calves was 1.114, which is consistent with the general situation in the Israeli Holstein population. Calvings were less frequent in the early summer months because fertility is lowest in late summer. Ten rounds of Fisher scoring iteration for TM and nine rounds of REML iteration for LM were performed for the analysis of DC in the complete data set. Computing times on a CRAY XMP-48 were about 12 and 4 rain for TM and LM, respectively, and less for the single Journal of Dairy Science Vol. 71, No. 9, 1988
2494
WELLER ET AL.
TABLE 1. Basic statistics of the data sets analyzed. Parity All
1
2+ 3
Incidence (%) of: Dystocia
5.08
Calf mortality Dystocia + calf mortality
5.43 2.59
Number of: Records Herd-year-seasons Sires of cow Sires of calf
253,724 8053 436 386
parity analyses. Computing times on most standard mainframe computers would be larger (16). Estimates of SIRE and SC variance components, as a fraction of the residual variance, and of heritabilities are in Table 3 for all analyses. Variance component estimates by LM were similar to those found in previous studies (1, 2, 13, 18, 19, 20, 21). Heritability estimates were higher for TM than for LM in all cases. The TM heritability estimates obtained using the counterpart of REML were in most cases slightly lower than the values obtained from the LM estimates using the equations of Dempster and Learner (3) and Gianola (7). Except for the analysis of CM by TM, the SC variance component was always larger than the SIRE component. The SIRE and SC variance component estimates for DC based on later parities were less than a third of the corresponding first parity estimates with TM and even less for LM. Variance components using first parity records were also larger than with all parities included. Heritability of CT was higher than that of CM in both LM and TM. It was higher than that of DC only as a trait of the SIRE in the LM analyses. Solutions for fixed effects for first parity DC and CM are presented in Table 4. Similar trends were evident for both TM and LM and for both traits. As in previous studies (2, 4, 13, 14, 18, 19, 20, 21), incidence of both DC and CM was higher for births of male calves. Dystocia and CM occurred in 10.9 and 10.2% of male calves and in 5.0 and 4.9% of female calves. As Journal of Dairy Science Vol. 71, No. 9, 1988
8.05 7.66 4.17 106,751 5882 377 329
2.92 3.81 1.44 146,973 7748 399 368
expected, both DC and CM tended to decrease with increasing calving age. Dystocia was lowest in June and highest in February and March. For CM May was most favorable and February least favorable. These results differed from the effect of season on production and fertility. Due to heat stress, both production and fertility are lowest in Israel in the late summer. No clear trends were evident for the group of SIRE or SC solutions. Correlations between TM and LM estimates of fixed effects are given in Table 5. Correlations are not given for the effects of sex for all analyses and calving age for the later parity analysis because there were only two levels for each of these effects. The correlations for GS ranged from .87 for later parity DC to .99 for first parity DC, CM, and CT. The correlations for GSC ranged from .75 for all parity DC to .99 for first parity CM. These results suggest that estimates of genetic trends by LM and TM may be different. This point is being investigated further. Correlations for age and calving month were all greater than .98. Correlations between estimates of fixed effects in the three DC analyses and in the three first parity analyses, separately for TM and LM, are given in Table 6. None of the correlations between first and later parity estimates of fixed effects for calving difficulty differed from zero (P>.05). Correlations between the DC and CM solutions were high for calving age and m o n t h but were not significant for group effects with the exception of GSC by LM.
GENETIC ANALYSIS OF CALVING TRAITS
2495
TABLE 2. Distribution of records by level of fixed effects. Parity All
1
2+3
Group of sires of cow 1 2 3 4 5 6
121,110 35,312 47,575 29,387 15,390 4753
35,296 14,749 24,755 19,152 8241 4558
85,814 20,563 22,820 10,235 7149 195
Group of sires of calf 1 2 3 4 5 6 7
46,041 26,876 42,960 90,317 18,470 15,792 13,261
18,641 14,271 22,282 45,705 3414 1576 862
27,400 12,605 20,678 44,612 15,056 14,216 12,399
133,805 119,919
55,361 51,390
78,444 68,529
Calving age, mo <22 22 23 24 25 26 27 28 29 >29 2nd parity 3rd parity
963 6609 27,704 31,897 18,386 10,017 5536 3126 1816 697 91,448 55,525
963 6609 27,704 31,897 18,386 10,017 5536 3126 1816 697 ... ...
... ... ... . .. ... ... ... . .. ... 911448 55,525
Calving month January February March April May June July August September October November December
27,575 24,061 22,042 17,306 15,040 11,623 13,457 16,532 24,051 26,546 26,892 28,599
11,233 9491 9233 8026 6921 4966 5242 6380 10, 335 11,159 11,470 12,295
16,342 14,570 12,809 9280 8119 6657 8215 10,152 13,716 15,387 15,422 16,304
Sex Male Female
Variances a n d s k e w n e s s values o f t h e S I R E a n d SC e v a l u a t i o n s ( s u m s o f t h e c o r r e s p o n d i n g sire a n d g r o u p s o l u t i o n s ) are in Table 7 f o r b o t h t h e T M a n d LM a n a l y s e s . Similar r e s u l t s w e r e o b t a i n e d f o r t h e S I R E a n d SC s o l u t i o n s a n d
t h e r e f o r e a r e n o t p r e s e n t e d . In all cases, varia n c e s derived f r o m T M a n a l y s e s w e r e h i g h e r t h a n t h e c o r r e s p o n d i n g LM values, e v e n t h o u g h t h e r e s i d u a l v a r i a n c e w a s set e q u a l f o r c o m p a r i son purposes. As expected from the magnitude Journal of Dairy Science Vol. 71, No. 9, 1988
2496
WELLER ET AL.
TABLE 3. "Variance component and heritability estimates obtained with the threshold and linear model analyses. Variance components a
Heritabilities 3
Trait
Parity
Model 1
SIRE
SC
SIRE
SC
Dystocia
All
TM LM TM LM TM LM
.0154 .0043 .0287 .0088
.0261 .0065 .0406 .0181
.054(.060) .015 .098(.100) .031
.0061
.0121
.022(.026)
.0011
.0022
.004
.091(.101) .023 .139(.204) .063 .042(.052) .008
First Later Calf mortality
First
TM LM
.0241 .0076
.0202 .0088
.084(.094) .027
.071(.108) .031
Compositetrait
First
TM LM
.0257 .0101
.0310 .0177
.089(.065) .036
.107(.110) .061
1TH = Threshold model; LM = linear model. ~As a fraction of residual variance. 3Herd_year_season component was assumed to be .1 of the residual variance. Heritability estimates derived from equations of Dempster and Lerner (3) and of Gianola (7) are in parentheses.
of the variance c o m p o n e n t estimates, the variances of the SC evaluations were higher in m o s t cases than the corresponding S I R E values. Coefficients of skewness were positive for all LM analyses and higher than the corresponding TM values. Skewness in the LM analyses probably results f r o m fitting a LM to discontinuous data with most of the observations in an e x t r e m e category. Correlations b e t w e e n TM and LM S I R E and SC evaluations are given in Table 8 for sires with at least 20 records in each analysis. Results were similar for correlations b e t w e e n sire solutions so these are not presented. Correlations ranged from .92 (SC, all parities DC) to .99 (SIRE and SC in CM evaluations). Even t h o u g h first parity correlations for DC were greater than .97, the median rank changes were 8 for S I R E and 4 for SC. M a x i m u m rank changes were 57 and 20, respectively. A l t h o u g h not shown, correlations were marginally higher for high repeatability sires. Correlations b e t w e e n evaluations obtained in different analyses for sires with at least 20 records in each analysis are s h o w n in Table 9. Results of the TM and LM analyses were similar in all cases. Correlations a m o n g high repeatability sires were marginally higher (data not shown). Correlations b e t w e e n first and later Journal of Dairy Science Vol. 71, No. 9, 1988
parity DC evaluations were less than .5 for b o t h SIRE and SC. This agrees with m o s t studies (1, 13, 17, 21) and strongly suggests that the genetic control of DC is different for first and later parities. Correlations b e t w e e n first parity DC and CM were about .7 for b o t h S I R E and SC evaluations. These results suggest that although the genetic correlation b e t w e e n these traits is high, there is still genetic variation in CM that is i n d e p e n d e n t of variation in DC. Correlations b e t w e e n CT and DC or CM evaluations were all larger t h a n .9. Correlations b e t w e e n S I R E and SC evaluations are given in Table 10 for b o t h TM and LM analyses. Sires w i t h less than 20 records for either evaluation were n o t included in the computations. Correlations were less than .3, but differed (P<.05) f r o m zero. This agrees with previous results obtained with LM analyses (17, 18, 19, 21). Even though the assumption of a zero covariance b e t w e e n these effects is incorrect, it m a y provide a reasonable a p p r o x i m a t i o n to the true variance structure.
DISCUSSION
Sire evaluation by TM in large data sets is feasible with currently available programs and c o m p u t i n g facilities. The program used here
GENETIC ANALYSIS OF CALVING TRAITS
2497
TABLE 4. Solutions for fixed effects for the first parity analyses of dystocia and calf mortality. 1 Dystocia
Calf mortality
Threshold model
Linear model
Threshold model
Linear model
Group o f s i r e s o f c o w 1 2 3 4 5 6
-1.82 -.07 -3.93 .36 .31 0
-.85 -.50 --1.75 .17 .05 0
.83 1.68 -1.39 1.14 .51 0
.52 1.04 -.76 ,60 .23 0
Group o f s i r e s o f c a l f 1 2 3 4 5 6 7
--1.63 -2.86 .18 --.52 --1.69 --.68 0
--.82 --.98 .31 --.24 --1.17 --.44 0
-1.40 -2.52 .02 --1.50 --4.90 --2.75 0
--.61 --1.30 .37 --.59 --2.67 --1.60 0
Sex Male Female
11.48 0
5.95 0
10.11 0
5.27 0
Calving, age, mo <22 22 23 24 25 26 27 28 29 >29
3.58 1.55 .31 --.20 --.92 --1.50 --1.22 --1.26 --1.53 0
2.07 .84 .36 .07 --.34 --.72 --.54 --.63 --.99 0
3.29 1.31 -.05 -.38 --.44 --.83 --.25 --.81 --1.13 0
2.07 ,84 .35 .06 --.34 --.72 --.54 --.63 --.99 0
Calving m o n t h January February March April May June July August September October November December
2.03 3.14 3.22 .19 --1.97 -3.39 -1.12 .68 .27 - .75 -.76 0
.98 1.64 1.61 --.13 --1.03 -1.69 -.61 .28 .32 - . 33 --. 30 0
.95 1.12 .74 -2.80 --3.51 -2.47 -.43 .08 -.71 - 1.06 --.95 0
.52 .64 .37 -1.59 --1.80 -1.36 -.30 --.02 -.29 -.47 --.45 0
(%)
Threshold model solutions were multiplied by the square root of the residual variance c o m p o n e n t from the linear model analyses.
Journal of Dairy Science Vol. 71, No. 9, 1988
2498
WELLER ET AL.
TABLE 5. Correlations between threslaold and linear model estimates of fixed effects. Effects Group of sires of
Calving
Trait
Parity
Cow
Calf
Age
Month
Dystocia
All First Later
.95 .99 .87
.75 .90 .96
.99 .99 . .
.98 .99 99
First First
.99 .99
.99 .98
.99 .99
Calf mortality Composite trait
.
.
.99 .99
•9 in all cases. R a n k changes were small in m o s t cases, a l t h o u g h 5% o f t h e sires h a d r a n k changes greater t h a n 30 f o r S I R E e v a l u a t i o n s a n d larger t h a n 14 for SC evaluations. T h e LM evaluations h a d higher skewness t h a n t h e TM evaluations. This implies t h a t p r o b a b i l i t y s t a t e m e n t s f o r t h e s e traits based o n sire e v a l u a t i o n s o b t a i n e d w i t h LM assuming n o r m a l i t y are s o m e w h a t d i s t o r t e d . A significant d i s a d v a n t a g e o f TM is t h a t at p r e s e n t t h e r e is n o k n o w n simple m e t h o d f o r assessing p r e d i c t i o n e r r o r variances
(16) was able t o e s t i m a t e v a r i a n c e c o m p o n e n t s a n d p r e d i c t b r e e d i n g values. O n c e variance comp o n e n t s are e s t i m a t e d , p r e d i c t i o n o f b r e e d i n g values is m u c h less t i m e c o n s u m i n g , especially if s o l u t i o n s are o b t a i n e d b y i t e r a t i o n o n d a t a w i t h o u t f o r m i n g t h e c o e f f i c i e n t m a t r i x (16). However, t h e advantages of TM over LM appeared t o b e slight in t h e d a t a set studied. T h e variances o f t h e T M e v a l u a t i o n s were larger t h a n t h o s e o f t h e c o r r e s p o n d i n g LM evaluations, b u t t h e i r c o r r e l a t i o n s w e r e g r e a t e r t h a n
TABLE 6. Correlations between estimated fixed effects in different analyses. Effects Group of sires of Analysis Dystocia: All vs. first parity All vs. later parity First vs. later parity First parity: DC vs. CM DC vs. CT CM vs. CT
Model
Cow
Calving Calf
Age
.99 .99
.83 .89 .37 .29 -.20 -.15
.99
.79 .81 .95 .96 .94 .94
Threshold Linear Threshold Linear Threshold Linear
.83 .83 .47 .49 -.07 .06
.88 .84 .85 .53 .68 .74
Threshold Linear Threshold Linear Threshold Linear
.67 .56 .96 .90 .85 .87
.62 .84 .85 .92 .93 .98
1DC = Dystocia, CM = calf mortality, CT = composite trait• Journal of Dairy Science Vol. 71, No. 9, 1988
.99
.99 .99 .98 .99
Month
GENETIC ANALYSIS OF CALVING TRAITS
2499
TABLE 7. Variance and skewness of sire evaluations by threshold and linear model analyses. Evaluations Sire of cow Variance (× 104)
Trait
Parity
Model I
Dystocia
All
TM LM TM LM TM LM
3.006 .962 9.720 2.813 .298 .061
.119 .523 --.004 .449 --.119 .121
6.137 1.801 8.842 4.581 1.863 .327
--.017 1.006 .052 .792 .227 .579
First Later
Skewness
Sire of calf Variance (× 104)
Skewness
Calf mortality
First
TM LM
6.523 2.157
.141 .417
6.134 2.629
.121 .508
Composite trait
First
TM LM
24.432 9.769
.074 .472
22.715 14.736
.282 .733
1TM = Threshold model, LM = linear model.
in large models. F u r t h e r study is required to d e t e r m i n e w h e t h e r m e t h o d s suitable for LM can be adapted to TM. Estimates o f variance c o m p o n e n t s obtained by R E M L were similar to previous estimates for the Israeli p o p u l a t i o n by Henderson's M e t h o d 3 using a slightly different m o d e l (21). These results agree w i t h o t h e r studies p e r f o r m e d on large data sets (13). A l t h o u g h conflicting opinions have been presented (1, 2, 13, 17, 18, 19, 20, 21), the results of this study clearly indicate that first and later parity CD should n o t be viewed as the same trait in genetic evaluation or parameter estimation. Heritabilities are lower for later parities, even in the threshold m o d e l analysis, which
should account for d e p e n d e n c e b e t w e e n m e a n incidence and variance. F u r t h e r m o r e , correlations b e t w e e n first and later parity evaluations w e r e low, and the estimates o f t h e same fixed effects were dissimilar. A multitrait analysis has been suggested (20), b u t little w o u l d be gained if later parity heritability and t h e genetic correlation b e t w e e n first and later parities are as low as the results of this study suggest. Also, a multitrait threshold analysis is considerably m o r e c o m p l e x than a single trait analysis (5). An advantage of including the relationship m a t r i x or of conducting a m u l t i p a r i t y analysis is the possibility o f obtaining evaluations for animals with no records. F o r example, sires m a t e d o n l y to cows but not to heifers w o u l d
TABLE 8. Correlations between sire evaluations obtained with threshold and linear models. 1 Evaluation Trait
Parity
Sire of cow
Sire of calf
Dystocia
All First Later
.96 (436) .98 (313) .96 (399)
.92 (386) .97 (166) .97 (368)
Calf mortality Composite trait
First First
.99 ( 313) .98 (313)
.99 (166) .98 (166)
1Number of sires in each correlation is listed in parentheses. Only sires with at least 20 records are included. Journal of Dairy Science Vol. 71, No. 9, 1988
2500
WELLER ET AL.
TABLE 9. Correlations between sire evaluations in different analyses. ~ Evaluations Analysis Dystocia: All vs. first parity All vs. later parity First vs. later parity First parity: DC vs. CM DC vs. CT CM vs. CT
Model
Sire of cow
Sire of calf
Threshold Linear Threshold Linear Threshold Linear
.90 .88 .65 .59 .29 .25
.88 .88 .80 .72 .49 .47
Threshold Linear Threshold Linear Threshold Linear
.69 (313) .69 ( 313) .91 (313) .92 (313) .91 (313) .91 (313)
( 313) (313) (362) (362) (244) (244)
(166) (166) (359) (359) (145) (145)
.73 (166) .75 (166) .94 (166) .94 (166) .90 (166) .92 (166)
Number of sires in each correlation is listed in parentheses. Only sires with at least 20 records are included.
have e v a l u a t i o n s for first p a r i t y d y s t o c i a based o n later p a r i t y records o f t h e i r m a t e s and o n t h e first and l a t e r p a r i t y p e r f o r m a n c e o f t h e i r relatives. However, t h e s e e v a l u a t i o n s will have low r e p e a t a b i l i t y even using b o t h sources o f i n f o r m a t i o n . T h u s , a n analysis using o n l y first p a r i t y r e c o r d s seems m o s t practical f o r r o u t i n e sire evaluation. Most studies have f o u n d a negative genetic c o r r e l a t i o n b e t w e e n direct a n d m a t e r n a l g e n e t i c effects o n DC (13, 18, 19). T h e small positive c o r r e l a t i o n s f o u n d b e t w e e n S I R E a n d SC evaluations a n d s o l u t i o n s do n o t c o n t r a d i c t t h e s e results. T h e SC variable includes 1/2 o f t h e direct genetic effect, w h e r e a s t h e SIRE
variable includes 1/2 o f t h e m a t e r n a l a n d 1 / 4 o f t h e d i r e c t genetic effects. T h u s , t h e c o v a r i a n c e b e t w e e n t h e SC a n d S I R E variables includes 1 / 4 of t h e c o v a r i a n c e b e t w e e n t h e d i r e c t a n d m a t e r nal genetic effects, a n d 1/8 of t h e direct g e n e t i c variance (18). T h e r e f o r e , even if t h e c o v a r i a n c e b e t w e e n t h e m a t e r n a l a n d d i r e c t genetic comp o n e n t s is negative, t h e c o r r e l a t i o n b e t w e e n t h e S I R E a n d SC s o l u t i o n s c o u l d still b e positive. Heritabilities were similar f o r DC a n d CM, b u t c o r r e l a t i o n s b e t w e e n first p a r i t y sire evaluations for t h e s e traits were less t h a n u n i t y , even for high r e p e a t a b i l i t y sires. This suggests t h a t t h e genetic c o r r e l a t i o n b e t w e e n t h e s e traits is n o t perfect. F u r t h e r m o r e , if t h e genetic corre-
TABLE 10. Correlations between sire of cow and sire of calf evaluations by threshold and linear models. 1 Models Trait
Parity
Threshold
Linear
Dystocia
All First Later
.21 (305) .29 (144) .26 (235)
.14 (305) .25 (144) .19 (235)
Calf mortality Composite trait
First First
.23 (144) .27 (144)
.22 (144) .26 (144)
1 Number of sires in each correlation is listed in parentheses. Journal of Dairy Science Vol. 71, No. 9, 1988
GENETIC ANALYSIS OF CALVING TRAITS lation b e t w e e n t h e s e traits is p e r f e c t , but t h e residual c o r r e l a t i o n is less t h a n 1, c o m b i n a t i o n o f DC and CM into CT should have increased heritability, b u t this was n o t t h e case. It w o u l d s e e m t h a t a m u l t i t r a i t evaluation o f first parity DC and CM w o u l d be a p p r o p r i a t e for these traits. T h e c o m p u t i n g strategy for such an evaluation w o u l d be facilitated in p o p u l a t i o n s w i t h g o o d r e c o r d i n g s y s t e m s , because virtually all calvings w o u l d have a record f o r b o t h traits. Multitrait analysis is greatly simplified if all animals have records o n all traits (12), b u t similar s i m p l i f i c a t i o n s are n o t k n o w n to exist in TM. If t h e residual c o r r e l a t i o n is in fact l o w e r t h a n t h e genetic correlation, a m u l t i t r a i t evaluation w o u l d yield m o r e accurate evaluations f o r b o t h traits t h a n t h o s e o b t a i n e d w i t h t h e t w o traits a n a l y z e d separately. ACKNOWLEDGMENTS This research was s u p p o r t e d b y t h e USIsrael Binational Agricultural R e s e a r c h and D e v e l o p m e n t F u n d , Project N u m b e r US-80584, and b y t h e National C e n t e r for Superc o m p u t i n g A p p l i c a t i o n s , C h a m p a i g n , IL. REFERENCES 1 Bar-Anan, R., M. Soller, and J. C. Bowman. 1976. Genetic and environmental factors affecting the incidence of difficult calving and perinatal calf mortality in Israeli Friesian dairy herds. Anim. Prod. 22:299. 2 Berger, P. J., and A. E. Freeman. 1978. Prediction of sire merit for calving difficulty. J. Dairy Sci. 61:1146. 3 Dempster, E. R., and 1. M. Lerner. 1950. Heritability of threshold characters (with Appendix by A. Robertson). Genetics 35:212. 4 Djemali, M., P. J. Berger, and A. E. Freeman. 1987. Ordered categorical sire evaluations for dystocia in Holsteins. J. Dairy Sci. 70:2374. 5 Foulley, J. L., S. lm, D. Gianola, and lna Hoschele. 1987. Empirical Bayes estimation of parameters for n polygenic binary traits. Genet. Set. Evot. 19: 197. 6 Gianola, D. 1980. A method of sire evaluation for
2501
dichotomies. J. Anirn. Sci. 51:1266. 7 Gianola, D. 1982. Theory and analysis of threshold characters. J. Anim. Sci. 54:1079. 8 Gianola, D., and J. L. Foulley. 1983. Sire evaluation for ordered categorical data with a threshold model. Genet. Sel. Evol. 15:201. 9 Gilmour, A. R., R. D. Anderson, and A. L. Rae. 1987. Variance components on an underlying scale for ordered multiple threshold categorical data using a generalized linear mixed model. Z. Tierz. Zuchtungsbiol. 104:149. 10 Harville, D. A., and R. W. Mee. 1984. A mixedmodel procedure for analyzing ordered categorical data. Biometrics 40:393. 11 Jensen, J. 1986. Sire evaluation for type traits with linear and nonlinear procedures. Livest. Prod. Sci. 15:165. 12 Lee, A. J. 1979. Mixed model, multiple trait evaluation of related sires when all traits are recorded. J. Anita. Sci. 48:1079. 13 Meijering, A. 1984. Dystocia and stillbirth in cattle-a review of causes, relations and implications. Livest. Prod. Sci. 11:143. 14 Meijering, A. 1985. Sire evaluation for calving traits by best linear unbiased prediction and nonlinear methodology. Z. Tierz. Zuchtungsbiol. 102:95. 15 Meijering, A., and D. Gianola. 1985. Linear vs. nonlinear methods of sire evaluation for categorical traits: a simulation study. Genet. Sel. Evol. 17: 115. 16 Misztal, I., D. Gianola, and J. L. Foulley. 1988. Computing aspects of a nonlinear method of sire evaluation for categorical data. J. Dairy Sci. 71: (in press). 17 Ron, M., R. Bar-Anan, and J. I. Welter. 1986. Sire and maternal grandsire effects on calving difficulty and calf mortality in Israeli Holsteins. J. Dairy Sci. 69:243. 18 Thompson, J. R., A. E. Freeman, and P. J. Berger. 1981. Age of dam and maternal effects for dystocia in Holsteins. J. Dairy Sci. 64:1603. 19 Thompson, J. R., and J.E.O. Rege. 1984. Influences of dam on calving difficulty and early calf mortality. J. Dairy Sci. 67:847. 20 Van Vleck, L. D., and K. M. Edlin. 1984. Multiple trait evaluation of bulls for calving ease. J. Dairy Sci. 67:3025. 21 Weller, J. I., E. Ezra, and R. Bar-Anan. 1986. Studies on the model of choice for genetic analysis of calving traits. J. Dairy Sci. 69(Suppl. 1):124. (Abstr.)
Journal of Dairy Science Vol. 71, No. 9, 1988
AB ST RACT
INTRODUCTION
Calvings of 106,751 Israeli Holstein heifers were analyzed for dystocia and calf mortality, scored dichotomously, and a composite trait, scored trichotomously. Dystocia was also studied with 146,973 second and third parity records. Models fitted included herd-year-season, sex of calf, calving age, calving month, sire of cow, sire of calf, and groups of sire of cow and of calf. Herd-year-season, sire of cow and calf, and residuals were random with diagonal variance-covariance matrices. Herd-year-season variance was assumed to be 10% of the residual component. Other variance components were estimated b y REML for linear models and by the counterpart of REML for threshold models. Heritability estimates were two to five times greater in threshold than in linear models, but correlations between corresponding sire evaluations were all >.9. Linear model sire evaluations were skewed positively, whereas threshold model evaluations had symmetrical distributions. Heritability for dystocia was greater in first than in later parities. Correlations between first and later parity sire evaluations were <.5. Thus, the genetic control of dystocia seems to be different for heifers and cows. Correlations between sire of cow and calf evaluations were <.3. Correlations between dystocia and calf mortality evaluations were about .7.
Dystocia (DC) and calf mortality (CM) are of major economic importance to dairy farmers. These traits have low heritabilities, high genetic correlation, and significant residual correlation (1, 2, 13, 17, 18, 19, 20, 21). Many studies have found the effect of sire of calf to be larger than that of the sire of the cow (1, 13, 17, 18, 19, 21). Although both traits are scored categorically, most studies have used standard mixed model methodology, even though this assumes the data to have a continuous distribution. Theoretical problems involved in this type of analysis have been described at length (6, 7). Several studies have suggested that the threshold model (TM), which assumes the existence of an underlying normal variable, is theoretically appropriate and computationally feasible for genetic analysis of categorical traits (3, 4, 6, 7, 8, 9, 10, 11, 14, 15, 16). However, only a few studies have actually applied this model to analysis of field data (4, 11, 14). Because TM equations are nonlinear and involve normal probability functions, computational complexity and computing resources required are greater than in a linear model (LM) analysis. Furthermore, TM genetic evaluations also require accurate estimates of variance components in the underlying scale. Although methods have been developed that are similar to REML, computations are more complex than those in standard REML algorithms (9, 10). Limitations in computing resources may pose a serious problem for sire evaluation based on large field data sets. However, recent developments in computers and programming techniques have made large-scale TM analyses feasible (16). The objectives of this study were to compare TM and LM analyses of calving traits using a large field data set o f Israeli Holsteins, to estimate genetic parameters for these traits and to find the model of choice for routine genetic evaluation.
Received January 7, 1988. Accepted May 9, 1988. 1Institute of Animal Sciences, Agricultural Research Organization, The Volcani Center, Bet Dagan, Israel 50250. 1988 J Dairy Sci 71:2491-2501
2491
2492
WELLER ET AL. MATERIALS AND METHODS
Data were 347,691 first through third parity calvings of milk-recorded Israeli Holsteins, calving between January 1978 and June 1985. Both DC and CM were scored as dichotomous traits. Dystocia was scored as 0 if the calving was normal or as 1 if it was difficult; CM was scored as 0 for a live birth or as 1 if the calf died within 48 h of birth. Records were deleted from the analysis for the following reasons: 1) sire of cow or of calf was unknown, 2) multiple births, 3) either CM or DC was undefined, 4) the calf was abnormal or died due to causes not related to calving, 5) abortions, 6) calvings of daughters of sires with less than 20 daughter calvings, and 7) calvings in which the sire of the calf sired less than 20 calves. Dystocia was recorded in all but two calvings. However, 2.1% of the records were deleted because either calf fate was unknown, the calf was abnormal, it died due to causes unrelated to calving, or it was aborted. There are no natural service matings in Israel, and about 1000 inseminations are performed for each young sire tested. Therefore restrictions 6 and 7 eliminated 1274 observations, consisting of recording mistakes, progeny of a few matings from imported semen, and progeny of sires tested at the boundaries of the time period considered. In total, 27% of the records (93,967) were deleted, mostly because the sire of the calf was unknown. All analyses were by both TM and LM. Dystocia was analyzed using all parities together and first and later parities separately. Based on these results, it was decided to analyze CM using only first parity calvings. Because these traits have a high genetic correlation, a composite trait (CT) was also analyzed. In the LM analysis, CT was the sum of the DC and CM scores; thus, each calving received a score of either 0, 1, or 2. In the TM analysis, it was assumed that these tbree classes were ordered so that the model included two thresholds. Thus five analyses were run for each of the two statistical methods: 1) DC, all parities; 2) DC, first parity; 3) DC, later parities; 4) CM, first parity; and 5) CT, first parity. The following model was used for analysis:
Journal of Dairy Science Vol. 71, No. 9, 1988
Yijklmnopq = HYSi + GSj + SIREjk + GSC I + SCIm + Sn + A o + Mp + eijklmnopq [1] where: Yijklmnopq = record on cow ijklmnopq; HYS i = random effect of herd-yearseason i; GSj = fixed effect of group j of sires of cows; SIREjk = random effect of sire of cow k in group j ; GSC 1 = fixed effect of group 1 of sires of calves; SClm = random effect of sire of calf m in group l; S n = fixed effect of sex n; A o = fixed effect of calving age o; Mp = fixed effect of calving month p; and eijklrnpopq = random residual. Calving seasons were October to March, and April to September. Sires were grouped by year of birth, with all sires born prior to 1972 in group 1. Groups 2 to 7 included all sires born during each subsequent 2-yr interval. Twelve classes were defined for calving month. Because HYS was random, there was no "confounding" between HYS and calving month. The first 10 calving ages were for primiparous cows representing calving ages from 21 to 30 mo. The 1st also included cows calving prior to 21 too, and the lOth included cows that calved after 30 mo. In the all parities analysis, second and third parity cows were included in classes 11 and 12, respectively, regardless of calving age. In the latter parity analysis, only two classes were defined for this effect. Analysis was facilitated by assuming that all random effects were mutually uneorrelated with diagonal dispersion matrices. It was further assumed that the variance of the HYS effect was known and equal to 10% of the residual variance. This value is similar to estimates calculated from US data (4). Using these
GENETIC ANALYSIS OF CALVING TRAITS assumptions, HYS effects, which had many levels in all analyses, could be absorbed. The assumption of diagonal dispersion matrices was inaccurate for the SIRE and SC effects, because sires were related. The assumption of zero correlation between effects was also an approximation, because several studies have found nonnegligible correlations between SIRE and SC effects (13, 17, 18, 21). A program based on principles described by Misztal et al. (16) and adapted to the CRAY XMP-48 supercomputer, was used for analysis. The SIRE, SC, and residual variance components were estimated by REML for LM model analyses and by the counterpart of REML for TM (9). After absorption of HYS effects, coefficient matrices were inverted at each round of iteration for both models. In TM analyses, interation involved Fisher's scoring and a modification of the expectation-maximization (EM) algorithm for variance components (16). Two to three EM rounds were performed per round of Fisher scoring iteration. At least six rounds of REML iteration in LM and six rounds of Fisher scoring iteration in TM were completed for each analysis. Iteration stopped when the change in all variance components was less than 1% of the values in the previous round. Sire evaluations were computed as the sum of the sire and sire group effects. Heritability (h 2) as a trait of the cow was calculated as: h 2 = 4 VAR(SIRE)/ [VAR(SIRE) + VAR(SC) + VAR(HYS) + VAR(E)]
[2]
where VAR(SIRE), VAR(SC), VAR(HYS), and VAR(E) are the SIRE, SC, HYS, and residual variance components, respectively. Heritability as a trait of the calf was computed in a similar manner except that VAR(SC) replaced VAR(SIRE) in the numerator. Threshold model heritabilities for dichotomous traits were also obtained dividing LM estimates by: z2/[p(1 - p)]
[3]
where z is the ordinate of the standard normal density function corresponding to p, the observed incidence in the first category (3). For CT, a trait with two thresholds, a formula of Gianola (7) was used.
2493
The unit of measurement of estimates and predictions in TM is the residual standard deviation. To facilitate comparison between LM and TM estimates, the TM solutions were multiplied by the square root of the corresponding LM residual variance component so that the units of the two models would be equivalent. Results obtained with the two different methods of analysis were compared in terms of ratios of variance components (heritability), the variance and skewness of the estimates and evaluations, and correlations between the solutions or evaluations. Although all sires included in the full data set had 20 or more calving records, some sires in the first and later parity subsets had less than 20 records. Therefore, correlations between sire evaluations were also computed for sires with at least 20 records in both analyses, and for sires with repeatability based on LM greater than .6. Repeatability was calculated as 1 - PEV/VC, where PEV is the prediction error variance of the SIRE or SC effect, and VC is the corresponding variance component. The PEV was computed as the diagonal d e m e n t of the inverse of the coefficient matrix at the final round of iteration times the residual variance. R ESU LTS
Basic statistics for the edited data set are in Table 1 for all parities combined and for first and later parities separately. Incidences of DC and CM for all parities were 5.1 and 5.4%, respectively; the incidence of either DC or CM in a single calving was 7.9%. Frequencies were higher for heifers, which is consistent with previous studies (1, 4, 13, 17, 20, 21). The distribution of records by level of fixed effects is presented in Table 2 for all parities combined, and for first and later parities separately. The ratio of male to female calves was 1.114, which is consistent with the general situation in the Israeli Holstein population. Calvings were less frequent in the early summer months because fertility is lowest in late summer. Ten rounds of Fisher scoring iteration for TM and nine rounds of REML iteration for LM were performed for the analysis of DC in the complete data set. Computing times on a CRAY XMP-48 were about 12 and 4 rain for TM and LM, respectively, and less for the single Journal of Dairy Science Vol. 71, No. 9, 1988
2494
WELLER ET AL.
TABLE 1. Basic statistics of the data sets analyzed. Parity All
1
2+ 3
Incidence (%) of: Dystocia
5.08
Calf mortality Dystocia + calf mortality
5.43 2.59
Number of: Records Herd-year-seasons Sires of cow Sires of calf
253,724 8053 436 386
parity analyses. Computing times on most standard mainframe computers would be larger (16). Estimates of SIRE and SC variance components, as a fraction of the residual variance, and of heritabilities are in Table 3 for all analyses. Variance component estimates by LM were similar to those found in previous studies (1, 2, 13, 18, 19, 20, 21). Heritability estimates were higher for TM than for LM in all cases. The TM heritability estimates obtained using the counterpart of REML were in most cases slightly lower than the values obtained from the LM estimates using the equations of Dempster and Learner (3) and Gianola (7). Except for the analysis of CM by TM, the SC variance component was always larger than the SIRE component. The SIRE and SC variance component estimates for DC based on later parities were less than a third of the corresponding first parity estimates with TM and even less for LM. Variance components using first parity records were also larger than with all parities included. Heritability of CT was higher than that of CM in both LM and TM. It was higher than that of DC only as a trait of the SIRE in the LM analyses. Solutions for fixed effects for first parity DC and CM are presented in Table 4. Similar trends were evident for both TM and LM and for both traits. As in previous studies (2, 4, 13, 14, 18, 19, 20, 21), incidence of both DC and CM was higher for births of male calves. Dystocia and CM occurred in 10.9 and 10.2% of male calves and in 5.0 and 4.9% of female calves. As Journal of Dairy Science Vol. 71, No. 9, 1988
8.05 7.66 4.17 106,751 5882 377 329
2.92 3.81 1.44 146,973 7748 399 368
expected, both DC and CM tended to decrease with increasing calving age. Dystocia was lowest in June and highest in February and March. For CM May was most favorable and February least favorable. These results differed from the effect of season on production and fertility. Due to heat stress, both production and fertility are lowest in Israel in the late summer. No clear trends were evident for the group of SIRE or SC solutions. Correlations between TM and LM estimates of fixed effects are given in Table 5. Correlations are not given for the effects of sex for all analyses and calving age for the later parity analysis because there were only two levels for each of these effects. The correlations for GS ranged from .87 for later parity DC to .99 for first parity DC, CM, and CT. The correlations for GSC ranged from .75 for all parity DC to .99 for first parity CM. These results suggest that estimates of genetic trends by LM and TM may be different. This point is being investigated further. Correlations for age and calving month were all greater than .98. Correlations between estimates of fixed effects in the three DC analyses and in the three first parity analyses, separately for TM and LM, are given in Table 6. None of the correlations between first and later parity estimates of fixed effects for calving difficulty differed from zero (P>.05). Correlations between the DC and CM solutions were high for calving age and m o n t h but were not significant for group effects with the exception of GSC by LM.
GENETIC ANALYSIS OF CALVING TRAITS
2495
TABLE 2. Distribution of records by level of fixed effects. Parity All
1
2+3
Group of sires of cow 1 2 3 4 5 6
121,110 35,312 47,575 29,387 15,390 4753
35,296 14,749 24,755 19,152 8241 4558
85,814 20,563 22,820 10,235 7149 195
Group of sires of calf 1 2 3 4 5 6 7
46,041 26,876 42,960 90,317 18,470 15,792 13,261
18,641 14,271 22,282 45,705 3414 1576 862
27,400 12,605 20,678 44,612 15,056 14,216 12,399
133,805 119,919
55,361 51,390
78,444 68,529
Calving age, mo <22 22 23 24 25 26 27 28 29 >29 2nd parity 3rd parity
963 6609 27,704 31,897 18,386 10,017 5536 3126 1816 697 91,448 55,525
963 6609 27,704 31,897 18,386 10,017 5536 3126 1816 697 ... ...
... ... ... . .. ... ... ... . .. ... 911448 55,525
Calving month January February March April May June July August September October November December
27,575 24,061 22,042 17,306 15,040 11,623 13,457 16,532 24,051 26,546 26,892 28,599
11,233 9491 9233 8026 6921 4966 5242 6380 10, 335 11,159 11,470 12,295
16,342 14,570 12,809 9280 8119 6657 8215 10,152 13,716 15,387 15,422 16,304
Sex Male Female
Variances a n d s k e w n e s s values o f t h e S I R E a n d SC e v a l u a t i o n s ( s u m s o f t h e c o r r e s p o n d i n g sire a n d g r o u p s o l u t i o n s ) are in Table 7 f o r b o t h t h e T M a n d LM a n a l y s e s . Similar r e s u l t s w e r e o b t a i n e d f o r t h e S I R E a n d SC s o l u t i o n s a n d
t h e r e f o r e a r e n o t p r e s e n t e d . In all cases, varia n c e s derived f r o m T M a n a l y s e s w e r e h i g h e r t h a n t h e c o r r e s p o n d i n g LM values, e v e n t h o u g h t h e r e s i d u a l v a r i a n c e w a s set e q u a l f o r c o m p a r i son purposes. As expected from the magnitude Journal of Dairy Science Vol. 71, No. 9, 1988
2496
WELLER ET AL.
TABLE 3. "Variance component and heritability estimates obtained with the threshold and linear model analyses. Variance components a
Heritabilities 3
Trait
Parity
Model 1
SIRE
SC
SIRE
SC
Dystocia
All
TM LM TM LM TM LM
.0154 .0043 .0287 .0088
.0261 .0065 .0406 .0181
.054(.060) .015 .098(.100) .031
.0061
.0121
.022(.026)
.0011
.0022
.004
.091(.101) .023 .139(.204) .063 .042(.052) .008
First Later Calf mortality
First
TM LM
.0241 .0076
.0202 .0088
.084(.094) .027
.071(.108) .031
Compositetrait
First
TM LM
.0257 .0101
.0310 .0177
.089(.065) .036
.107(.110) .061
1TH = Threshold model; LM = linear model. ~As a fraction of residual variance. 3Herd_year_season component was assumed to be .1 of the residual variance. Heritability estimates derived from equations of Dempster and Lerner (3) and of Gianola (7) are in parentheses.
of the variance c o m p o n e n t estimates, the variances of the SC evaluations were higher in m o s t cases than the corresponding S I R E values. Coefficients of skewness were positive for all LM analyses and higher than the corresponding TM values. Skewness in the LM analyses probably results f r o m fitting a LM to discontinuous data with most of the observations in an e x t r e m e category. Correlations b e t w e e n TM and LM S I R E and SC evaluations are given in Table 8 for sires with at least 20 records in each analysis. Results were similar for correlations b e t w e e n sire solutions so these are not presented. Correlations ranged from .92 (SC, all parities DC) to .99 (SIRE and SC in CM evaluations). Even t h o u g h first parity correlations for DC were greater than .97, the median rank changes were 8 for S I R E and 4 for SC. M a x i m u m rank changes were 57 and 20, respectively. A l t h o u g h not shown, correlations were marginally higher for high repeatability sires. Correlations b e t w e e n evaluations obtained in different analyses for sires with at least 20 records in each analysis are s h o w n in Table 9. Results of the TM and LM analyses were similar in all cases. Correlations a m o n g high repeatability sires were marginally higher (data not shown). Correlations b e t w e e n first and later Journal of Dairy Science Vol. 71, No. 9, 1988
parity DC evaluations were less than .5 for b o t h SIRE and SC. This agrees with m o s t studies (1, 13, 17, 21) and strongly suggests that the genetic control of DC is different for first and later parities. Correlations b e t w e e n first parity DC and CM were about .7 for b o t h S I R E and SC evaluations. These results suggest that although the genetic correlation b e t w e e n these traits is high, there is still genetic variation in CM that is i n d e p e n d e n t of variation in DC. Correlations b e t w e e n CT and DC or CM evaluations were all larger t h a n .9. Correlations b e t w e e n S I R E and SC evaluations are given in Table 10 for b o t h TM and LM analyses. Sires w i t h less than 20 records for either evaluation were n o t included in the computations. Correlations were less than .3, but differed (P<.05) f r o m zero. This agrees with previous results obtained with LM analyses (17, 18, 19, 21). Even though the assumption of a zero covariance b e t w e e n these effects is incorrect, it m a y provide a reasonable a p p r o x i m a t i o n to the true variance structure.
DISCUSSION
Sire evaluation by TM in large data sets is feasible with currently available programs and c o m p u t i n g facilities. The program used here
GENETIC ANALYSIS OF CALVING TRAITS
2497
TABLE 4. Solutions for fixed effects for the first parity analyses of dystocia and calf mortality. 1 Dystocia
Calf mortality
Threshold model
Linear model
Threshold model
Linear model
Group o f s i r e s o f c o w 1 2 3 4 5 6
-1.82 -.07 -3.93 .36 .31 0
-.85 -.50 --1.75 .17 .05 0
.83 1.68 -1.39 1.14 .51 0
.52 1.04 -.76 ,60 .23 0
Group o f s i r e s o f c a l f 1 2 3 4 5 6 7
--1.63 -2.86 .18 --.52 --1.69 --.68 0
--.82 --.98 .31 --.24 --1.17 --.44 0
-1.40 -2.52 .02 --1.50 --4.90 --2.75 0
--.61 --1.30 .37 --.59 --2.67 --1.60 0
Sex Male Female
11.48 0
5.95 0
10.11 0
5.27 0
Calving, age, mo <22 22 23 24 25 26 27 28 29 >29
3.58 1.55 .31 --.20 --.92 --1.50 --1.22 --1.26 --1.53 0
2.07 .84 .36 .07 --.34 --.72 --.54 --.63 --.99 0
3.29 1.31 -.05 -.38 --.44 --.83 --.25 --.81 --1.13 0
2.07 ,84 .35 .06 --.34 --.72 --.54 --.63 --.99 0
Calving m o n t h January February March April May June July August September October November December
2.03 3.14 3.22 .19 --1.97 -3.39 -1.12 .68 .27 - .75 -.76 0
.98 1.64 1.61 --.13 --1.03 -1.69 -.61 .28 .32 - . 33 --. 30 0
.95 1.12 .74 -2.80 --3.51 -2.47 -.43 .08 -.71 - 1.06 --.95 0
.52 .64 .37 -1.59 --1.80 -1.36 -.30 --.02 -.29 -.47 --.45 0
(%)
Threshold model solutions were multiplied by the square root of the residual variance c o m p o n e n t from the linear model analyses.
Journal of Dairy Science Vol. 71, No. 9, 1988
2498
WELLER ET AL.
TABLE 5. Correlations between threslaold and linear model estimates of fixed effects. Effects Group of sires of
Calving
Trait
Parity
Cow
Calf
Age
Month
Dystocia
All First Later
.95 .99 .87
.75 .90 .96
.99 .99 . .
.98 .99 99
First First
.99 .99
.99 .98
.99 .99
Calf mortality Composite trait
.
.
.99 .99
•9 in all cases. R a n k changes were small in m o s t cases, a l t h o u g h 5% o f t h e sires h a d r a n k changes greater t h a n 30 f o r S I R E e v a l u a t i o n s a n d larger t h a n 14 for SC evaluations. T h e LM evaluations h a d higher skewness t h a n t h e TM evaluations. This implies t h a t p r o b a b i l i t y s t a t e m e n t s f o r t h e s e traits based o n sire e v a l u a t i o n s o b t a i n e d w i t h LM assuming n o r m a l i t y are s o m e w h a t d i s t o r t e d . A significant d i s a d v a n t a g e o f TM is t h a t at p r e s e n t t h e r e is n o k n o w n simple m e t h o d f o r assessing p r e d i c t i o n e r r o r variances
(16) was able t o e s t i m a t e v a r i a n c e c o m p o n e n t s a n d p r e d i c t b r e e d i n g values. O n c e variance comp o n e n t s are e s t i m a t e d , p r e d i c t i o n o f b r e e d i n g values is m u c h less t i m e c o n s u m i n g , especially if s o l u t i o n s are o b t a i n e d b y i t e r a t i o n o n d a t a w i t h o u t f o r m i n g t h e c o e f f i c i e n t m a t r i x (16). However, t h e advantages of TM over LM appeared t o b e slight in t h e d a t a set studied. T h e variances o f t h e T M e v a l u a t i o n s were larger t h a n t h o s e o f t h e c o r r e s p o n d i n g LM evaluations, b u t t h e i r c o r r e l a t i o n s w e r e g r e a t e r t h a n
TABLE 6. Correlations between estimated fixed effects in different analyses. Effects Group of sires of Analysis Dystocia: All vs. first parity All vs. later parity First vs. later parity First parity: DC vs. CM DC vs. CT CM vs. CT
Model
Cow
Calving Calf
Age
.99 .99
.83 .89 .37 .29 -.20 -.15
.99
.79 .81 .95 .96 .94 .94
Threshold Linear Threshold Linear Threshold Linear
.83 .83 .47 .49 -.07 .06
.88 .84 .85 .53 .68 .74
Threshold Linear Threshold Linear Threshold Linear
.67 .56 .96 .90 .85 .87
.62 .84 .85 .92 .93 .98
1DC = Dystocia, CM = calf mortality, CT = composite trait• Journal of Dairy Science Vol. 71, No. 9, 1988
.99
.99 .99 .98 .99
Month
GENETIC ANALYSIS OF CALVING TRAITS
2499
TABLE 7. Variance and skewness of sire evaluations by threshold and linear model analyses. Evaluations Sire of cow Variance (× 104)
Trait
Parity
Model I
Dystocia
All
TM LM TM LM TM LM
3.006 .962 9.720 2.813 .298 .061
.119 .523 --.004 .449 --.119 .121
6.137 1.801 8.842 4.581 1.863 .327
--.017 1.006 .052 .792 .227 .579
First Later
Skewness
Sire of calf Variance (× 104)
Skewness
Calf mortality
First
TM LM
6.523 2.157
.141 .417
6.134 2.629
.121 .508
Composite trait
First
TM LM
24.432 9.769
.074 .472
22.715 14.736
.282 .733
1TM = Threshold model, LM = linear model.
in large models. F u r t h e r study is required to d e t e r m i n e w h e t h e r m e t h o d s suitable for LM can be adapted to TM. Estimates o f variance c o m p o n e n t s obtained by R E M L were similar to previous estimates for the Israeli p o p u l a t i o n by Henderson's M e t h o d 3 using a slightly different m o d e l (21). These results agree w i t h o t h e r studies p e r f o r m e d on large data sets (13). A l t h o u g h conflicting opinions have been presented (1, 2, 13, 17, 18, 19, 20, 21), the results of this study clearly indicate that first and later parity CD should n o t be viewed as the same trait in genetic evaluation or parameter estimation. Heritabilities are lower for later parities, even in the threshold m o d e l analysis, which
should account for d e p e n d e n c e b e t w e e n m e a n incidence and variance. F u r t h e r m o r e , correlations b e t w e e n first and later parity evaluations w e r e low, and the estimates o f t h e same fixed effects were dissimilar. A multitrait analysis has been suggested (20), b u t little w o u l d be gained if later parity heritability and t h e genetic correlation b e t w e e n first and later parities are as low as the results of this study suggest. Also, a multitrait threshold analysis is considerably m o r e c o m p l e x than a single trait analysis (5). An advantage of including the relationship m a t r i x or of conducting a m u l t i p a r i t y analysis is the possibility o f obtaining evaluations for animals with no records. F o r example, sires m a t e d o n l y to cows but not to heifers w o u l d
TABLE 8. Correlations between sire evaluations obtained with threshold and linear models. 1 Evaluation Trait
Parity
Sire of cow
Sire of calf
Dystocia
All First Later
.96 (436) .98 (313) .96 (399)
.92 (386) .97 (166) .97 (368)
Calf mortality Composite trait
First First
.99 ( 313) .98 (313)
.99 (166) .98 (166)
1Number of sires in each correlation is listed in parentheses. Only sires with at least 20 records are included. Journal of Dairy Science Vol. 71, No. 9, 1988
2500
WELLER ET AL.
TABLE 9. Correlations between sire evaluations in different analyses. ~ Evaluations Analysis Dystocia: All vs. first parity All vs. later parity First vs. later parity First parity: DC vs. CM DC vs. CT CM vs. CT
Model
Sire of cow
Sire of calf
Threshold Linear Threshold Linear Threshold Linear
.90 .88 .65 .59 .29 .25
.88 .88 .80 .72 .49 .47
Threshold Linear Threshold Linear Threshold Linear
.69 (313) .69 ( 313) .91 (313) .92 (313) .91 (313) .91 (313)
( 313) (313) (362) (362) (244) (244)
(166) (166) (359) (359) (145) (145)
.73 (166) .75 (166) .94 (166) .94 (166) .90 (166) .92 (166)
Number of sires in each correlation is listed in parentheses. Only sires with at least 20 records are included.
have e v a l u a t i o n s for first p a r i t y d y s t o c i a based o n later p a r i t y records o f t h e i r m a t e s and o n t h e first and l a t e r p a r i t y p e r f o r m a n c e o f t h e i r relatives. However, t h e s e e v a l u a t i o n s will have low r e p e a t a b i l i t y even using b o t h sources o f i n f o r m a t i o n . T h u s , a n analysis using o n l y first p a r i t y r e c o r d s seems m o s t practical f o r r o u t i n e sire evaluation. Most studies have f o u n d a negative genetic c o r r e l a t i o n b e t w e e n direct a n d m a t e r n a l g e n e t i c effects o n DC (13, 18, 19). T h e small positive c o r r e l a t i o n s f o u n d b e t w e e n S I R E a n d SC evaluations a n d s o l u t i o n s do n o t c o n t r a d i c t t h e s e results. T h e SC variable includes 1/2 o f t h e direct genetic effect, w h e r e a s t h e SIRE
variable includes 1/2 o f t h e m a t e r n a l a n d 1 / 4 o f t h e d i r e c t genetic effects. T h u s , t h e c o v a r i a n c e b e t w e e n t h e SC a n d S I R E variables includes 1 / 4 of t h e c o v a r i a n c e b e t w e e n t h e d i r e c t a n d m a t e r nal genetic effects, a n d 1/8 of t h e direct g e n e t i c variance (18). T h e r e f o r e , even if t h e c o v a r i a n c e b e t w e e n t h e m a t e r n a l a n d d i r e c t genetic comp o n e n t s is negative, t h e c o r r e l a t i o n b e t w e e n t h e S I R E a n d SC s o l u t i o n s c o u l d still b e positive. Heritabilities were similar f o r DC a n d CM, b u t c o r r e l a t i o n s b e t w e e n first p a r i t y sire evaluations for t h e s e traits were less t h a n u n i t y , even for high r e p e a t a b i l i t y sires. This suggests t h a t t h e genetic c o r r e l a t i o n b e t w e e n t h e s e traits is n o t perfect. F u r t h e r m o r e , if t h e genetic corre-
TABLE 10. Correlations between sire of cow and sire of calf evaluations by threshold and linear models. 1 Models Trait
Parity
Threshold
Linear
Dystocia
All First Later
.21 (305) .29 (144) .26 (235)
.14 (305) .25 (144) .19 (235)
Calf mortality Composite trait
First First
.23 (144) .27 (144)
.22 (144) .26 (144)
1 Number of sires in each correlation is listed in parentheses. Journal of Dairy Science Vol. 71, No. 9, 1988
GENETIC ANALYSIS OF CALVING TRAITS lation b e t w e e n t h e s e traits is p e r f e c t , but t h e residual c o r r e l a t i o n is less t h a n 1, c o m b i n a t i o n o f DC and CM into CT should have increased heritability, b u t this was n o t t h e case. It w o u l d s e e m t h a t a m u l t i t r a i t evaluation o f first parity DC and CM w o u l d be a p p r o p r i a t e for these traits. T h e c o m p u t i n g strategy for such an evaluation w o u l d be facilitated in p o p u l a t i o n s w i t h g o o d r e c o r d i n g s y s t e m s , because virtually all calvings w o u l d have a record f o r b o t h traits. Multitrait analysis is greatly simplified if all animals have records o n all traits (12), b u t similar s i m p l i f i c a t i o n s are n o t k n o w n to exist in TM. If t h e residual c o r r e l a t i o n is in fact l o w e r t h a n t h e genetic correlation, a m u l t i t r a i t evaluation w o u l d yield m o r e accurate evaluations f o r b o t h traits t h a n t h o s e o b t a i n e d w i t h t h e t w o traits a n a l y z e d separately. ACKNOWLEDGMENTS This research was s u p p o r t e d b y t h e USIsrael Binational Agricultural R e s e a r c h and D e v e l o p m e n t F u n d , Project N u m b e r US-80584, and b y t h e National C e n t e r for Superc o m p u t i n g A p p l i c a t i o n s , C h a m p a i g n , IL. REFERENCES 1 Bar-Anan, R., M. Soller, and J. C. Bowman. 1976. Genetic and environmental factors affecting the incidence of difficult calving and perinatal calf mortality in Israeli Friesian dairy herds. Anim. Prod. 22:299. 2 Berger, P. J., and A. E. Freeman. 1978. Prediction of sire merit for calving difficulty. J. Dairy Sci. 61:1146. 3 Dempster, E. R., and 1. M. Lerner. 1950. Heritability of threshold characters (with Appendix by A. Robertson). Genetics 35:212. 4 Djemali, M., P. J. Berger, and A. E. Freeman. 1987. Ordered categorical sire evaluations for dystocia in Holsteins. J. Dairy Sci. 70:2374. 5 Foulley, J. L., S. lm, D. Gianola, and lna Hoschele. 1987. Empirical Bayes estimation of parameters for n polygenic binary traits. Genet. Set. Evot. 19: 197. 6 Gianola, D. 1980. A method of sire evaluation for
2501
dichotomies. J. Anirn. Sci. 51:1266. 7 Gianola, D. 1982. Theory and analysis of threshold characters. J. Anim. Sci. 54:1079. 8 Gianola, D., and J. L. Foulley. 1983. Sire evaluation for ordered categorical data with a threshold model. Genet. Sel. Evol. 15:201. 9 Gilmour, A. R., R. D. Anderson, and A. L. Rae. 1987. Variance components on an underlying scale for ordered multiple threshold categorical data using a generalized linear mixed model. Z. Tierz. Zuchtungsbiol. 104:149. 10 Harville, D. A., and R. W. Mee. 1984. A mixedmodel procedure for analyzing ordered categorical data. Biometrics 40:393. 11 Jensen, J. 1986. Sire evaluation for type traits with linear and nonlinear procedures. Livest. Prod. Sci. 15:165. 12 Lee, A. J. 1979. Mixed model, multiple trait evaluation of related sires when all traits are recorded. J. Anita. Sci. 48:1079. 13 Meijering, A. 1984. Dystocia and stillbirth in cattle-a review of causes, relations and implications. Livest. Prod. Sci. 11:143. 14 Meijering, A. 1985. Sire evaluation for calving traits by best linear unbiased prediction and nonlinear methodology. Z. Tierz. Zuchtungsbiol. 102:95. 15 Meijering, A., and D. Gianola. 1985. Linear vs. nonlinear methods of sire evaluation for categorical traits: a simulation study. Genet. Sel. Evol. 17: 115. 16 Misztal, I., D. Gianola, and J. L. Foulley. 1988. Computing aspects of a nonlinear method of sire evaluation for categorical data. J. Dairy Sci. 71: (in press). 17 Ron, M., R. Bar-Anan, and J. I. Welter. 1986. Sire and maternal grandsire effects on calving difficulty and calf mortality in Israeli Holsteins. J. Dairy Sci. 69:243. 18 Thompson, J. R., A. E. Freeman, and P. J. Berger. 1981. Age of dam and maternal effects for dystocia in Holsteins. J. Dairy Sci. 64:1603. 19 Thompson, J. R., and J.E.O. Rege. 1984. Influences of dam on calving difficulty and early calf mortality. J. Dairy Sci. 67:847. 20 Van Vleck, L. D., and K. M. Edlin. 1984. Multiple trait evaluation of bulls for calving ease. J. Dairy Sci. 67:3025. 21 Weller, J. I., E. Ezra, and R. Bar-Anan. 1986. Studies on the model of choice for genetic analysis of calving traits. J. Dairy Sci. 69(Suppl. 1):124. (Abstr.)
Journal of Dairy Science Vol. 71, No. 9, 1988