Introduction

Longevity in dairy cows can be measured as the time from first parturition to animal cull or death, an interval called productive life (PL). In a producing dairy cow, PL is influenced by production level, fertility, health, etc. Of particular interest is involuntary culling, that is, for reasons beyond the control of the farmer (e.g. disease or reproductive problems). This process is known as functional productive life (FPL) when it includes cow production level in the model in an effort to correct for voluntary culling^1,2.

The study of FPL using survival analysis techniques allows incorporation of time-dependent covariables, use of incomplete (censured) records, analysis of large data bases and application of mixed models that permit national level genetic evaluations. However, the way these evaluations are done and the models used to do so, can vary depending on circumstances^3,4,5.

Previous studies using data from Mexico have analyzed the effect of milk production level on PL either based on a Kaplan-Meier estimator and a Weibull regression model⁶, or using mixed linear models and treating longevity as the ability to remain in the herd at 48 mo of age with a PL extending to the third lactation⁷. No studies have been done, however, addressing the use of survival analysis in predicting genetic values in a Holstein population in Mexico under milk recording. Genetic values are needed to develop effective genetic improvement systems for this trait on artificial insemination bulls used in Mexico. This trait is vital to genetic improvement of dairy cows for reasons of economy (e.g. production system sustainability) and animal health (improved PL reduces risk of disease)^6,7.

The present study objectives were to estimate heritabilities and predict genetic values for FPL by applying a proportional risk model with a Weibull distribution in a Holstein population in Mexico.

Material and methods

Data was for Holstein cows with first parturition between 1995 and 2008, and was provided by the Mexican Holstein Association (Asociación Holstein de México) from official milk recording and genealogical records. Productive life (PL) was calculated as the number of days between date of first parturition and date of culling or censure. Maximum credit for lactation was 305 d and censured animals were defined as records of cows sold to other farms for production, of live cows in a herd when its production stopped being recorded, or of cows with non-consecutive lactations. For greater accuracy, and to ensure model convergence, lactations were eliminated for any one of four reasons: 1) Abnormal standardized milk production, i.e. <2,500 kg (1^st percentile) or >17,000 kg (99^th percentile); 2) No data for first lactation; 3) First lactation at <17 mo of age; and 4) Daughters of sires with less than five daughters. The final database consisted of 36,507 PL records.

Parameter estimations and genetic value predictions were calculated with the Survival Kit ver. 3.12 program^8,9. A sire-maternal grandsire proportional risk model was used in which a Weibull distribution was assumed for the basic risk function. The FPL analysis model was:

In this model, h(t) represents a cow elimination risk in time t. The basic risk function is λρ(λt)^ρ-1, where ρ is the distribution shape parameter and λ corresponds to the distribution scale parameter. HY_j (t) is the random effect of the j-th herd-year of first parturition. This had a log-gamma function associated with parameter γ, included to acknowledge that animal elimination is a decision made by comparing animals within similar or contemporaneous groups, the effect of which is not explained in this study⁸. AP_k (t) is the k-th effect of age at first parturition in months, divided into four categories (1, ≤23 mo; 2, 24-25 mo; 3, 26-27 mo; and 4, ≥28 mo). LP_lm (t) is the time-dependent effect of the l-th lactation phase divided into three times (1-29, 30-249 and 250-305 d after parturition), within the m-th lactation (i.e. 1, 2, 3 and ≥4), assuming that risk does not change within each segment. MPL_n (t) is the effect of the n-th milk production level within herd-year of first parturition; production levels were calculated by organizing production records in ascending order to generate the mean and standard deviation. Finally, deciles were calculated based on normal milk production to represent production level per lactation, divided into ten categories, 1 being the lowest and 10 the highest. Sire (S _q ) and maternal grandsire (S _mg ) were random ancestor effects and were grouped into one vector (S) for analysis, assuming that they followed a normal multivariate distribution with Aσs2 as the variance - covariance matrix.

Heritabilities were calculated on a logarithmic scale (hlog2)^4,9, that is, the transformation of time (T) in a logarithmic scale, and in the original scale (ho2), following Ducrocq¹⁰ as described in Chirinos et al⁴:

Where: Ψ¹ γ_h, is the trigamma function evaluated in γ_h ; var(s) is sire’s variance; and π26 is the extreme value distribution’s variance.

Where: υ is the negative Euler number; ρ is the shape of the basic risk function; and hlog2 is the heritability on a logarithmic scale.

To confirm that the original scale heritability expression depended on the Weibull distribution estimators as reported elsewhere^4,10, effective and equivalent heritability were calculated as follows:

Where ρ is the proportion of complete records, which is applied in the present study.

The reliability of the predicted genetic values for sires was calculated based on records of their daughters¹⁰:

Where: R is the genetic value reliability of sires; N is the number of daughters with recorded FPL; and he2 is the equivalent heritability.

Model effects were tested with a partial likelihood ratio test, and estimation was done using the maximum likelihood method⁸.

Results

In the total sample of 36,507 FPL records, censoring rate was 25.5 %, cows had an average of 2.23 lactations, average time to cull/death was 615 d and average time to censoring was 731 d. The data included the daughters of 1,116 sires and 1,684 maternal grandsires.

In terms of relative risk (average risk= 1) under the influence of environmental and genetic effects, the youngest cows (categ. 1, ≤23 mo) had a greater culling risk than the cows in lactation categories 2, 3 and ≥4. Age at first parturition apparently had no effect on FPL (Table 1, Figure 1).

Table 1 Proportion (%) and number (N) of observations by phase within lactation categories and phases 

Figure 1 Relative culling risk (RR) in four categories of age (in months) at first parturition for Holstein cattle in Mexico. 

Standardized production level (SPL) data showed the cows in the lowest levels as more likely to be culled than a cow with an average SPL (categ. 6); this was 49 times more probable in SPL 1, 5 times greater in SPL 2 and twice as likely in SPL 3 (Table 2, Figure 2). Risk diminished in levels 7, 8 and 9, and then increased in level 10, the most productive animals.

Table 2 Proportion (%) and number (N) of observations in each standardized production level (SPL). 

Figure 2 Relative culling risk (RR) in ten standardized production levels (SPL); 1 is lowest production level and 10 is highest. 

In cows with consecutive lactations, relative risk of discard tended to increase concurrently with milking days and lactation number, reaching a maximum of 305 d after parturition (Figure 3). The culling rate was particularly intense during the third lactation phase (LP; 250-305 d).

Figure 3 Relative culling risk (RR) in cows with consecutive lactations (1, 2, 3 and ≤4); each lactation is divided into three phases (0-29, 30-249 and 250-305 d). 

Table 3 Estimated genetic parameters values for survival of Holstein cattle in Mexico 

The estimators for the predicted genetic values varied from -0.68 to 1.62, with a mean of zero. Relative risk of culling for FPL in the daughters of the studied sires ranged from 0.51 to 5.03. Negative genetic values indicated a low cull risk for the daughters of a sire and consequently higher FPL. Among the 2,800 analyzed sires, genetic value reliability varied from 0 to 0.95, with 179 sires having a level ≥0.50 (Figure 4).

Figure 4 Genetic value (GV) reliabilities for studied sires based on number of daughters with complete records. 

Genetic tendencies for FPL of the studied Mexican Holstein sires by birth year showed increasing risk from 1997 to 2003, followed by a slight decrease in 2003 (Figure 5).

Figure 5 Genetic trends of estimated genetic values (EGV) for length of functional productive life in Mexican Holstein sires born between 1995 and 2003. 

Discussion

Lactations per PL was low on average in the present data, suggesting that the cows had not yet expressed their maximum productive potential^11,12,13. This implies that they were discarded shortly after completing two lactations, which coincides with previous reports^14-17.

The 25.5 % censoring rate in the records due to incomplete data is within the 9.4 to 73.4 % range reported in other studies using the same methodology in Holstein populations^3,4,5.

All the effects analyzed in the model were significant (P<0.0001), although SPL and LP made the greatest contribution to calculating culling risk. Age at first parturition made a lesser contribution. In other studies, this effect has not demonstrated a clear tendency to explain relative culling risk^4,11. However, in the present data, it did show a slight increase in cows with first parturition at ≤23 mo. This may be due to the fact that this group includes cows that gave birth at young ages because they were pregnant before 14 mo of age, as well as those that did not complete the full 9 mo of their first gestation, which can cause problems in later parturitions. In older age categories, age at first parturition had no apparent effect on FPL, as has been reported for other Holstein populations^4,8.

Relative culling risk associated with standardized production levels (SPL) showed the lowest levels (1, 2 and 3) to have the highest risk. Voluntary culling in response to low SPL is clearly significant in this population. Previous studies done with data from the same population and other Holstein populations also documented greater culling risk in less productive animals^4,11. Risk also increased slightly beginning at SPL 8, suggesting that high producing cows experience more intense physiological stress than lower producers, which can affect their FPL^8,9.

Culling risk also increased in the final days of a lactation in cows with consecutive lactations. This reflects elimination of non-gestating cows at the end of lactation, possibly due to reproductive or health problems^16,17. These tendencies coincide with previous studies reporting increased culling rates at the end of each lactation. The way LP is defined can influence these results since the first LP level will always be more homogeneous than the last ones. At later LP levels, age effects begin to accumulate but this are not included in calculation of the PL by days in lactation beyond 305 d, nor in inter-parturition periods, which are not part of this methodology.

The estimated parameter was 2.37, a value within the 0.36 to 5.0 reported for other Holstein populations^4,9,10. This value supports the recorded data since values greater than 1 indicate increased culling risk as animals age. The gamma parameter represents the variance within herd-year. Values greater than 1 indicate heterogeneity in culling practices⁸. The 4.28 gamma value in the present study is similar to those reported in other Holstein populations^3,4,13. Variation in culling practices may be related to varying levels of mechanization throughout Mexico, an aspect requiring further research.

Heritabilities in this study were 0.08 in the logarithmic scale, 0.13 in the original scale, 0.12 in the effective scale and 0.10 for equivalent heritability. All these are within published ranges: 0.02 to 0.11 for the logarithmic scale and 0.038 to 0.22 for the original scale^3,4,8; and 0.048 to 0.108 for the effective and equivalent scales^4,10,18. However, the estimators were higher than reported in a previous study in the same population using linear models⁷. These discrepancies in heritabilities occur because the survival model used in the present study is superior to previous models for explaining environmental effects. It can consequently isolate the additive genetic effects associated with the trait of interest, and is better for estimating genetic parameters and predicting FPL genetic values. Calculations done using these procedures produced heritability values similar to those reported elsewhere^8,10,18. Recent research has demonstrated extensions of the formulas to other models for the genetic effects of other animal species¹⁹. These do not depend on ñ or the extreme value distribution, which was not used in the study.

Reliability of the sires’ genetic values ranged from 0.0 to 0.96. This variability is caused by the fact that the reliability of genetic predictions made using survival models is calculated based on the number of daughters with complete (observed) records and not a sire’s total number of daughters. Because of this, younger sires may be at a disadvantage since most of their daughters will still be alive at the time of evaluation, substantially lowering prediction reliability^8,9. The same reliability behavior has been reported in countries such as Spain, France and Switzerland^4,10,11, where survival models help produce better predictors of sire FPL genetic values.

Sire FPL genetic values calculated based on daughter longevity in Mexico exhibited an increase in relative culling risk as sire birth year advanced. This increase implies an important decrease in genetic values that may be due to increased selection of sires by breeders for predicted transmission of improved milk production in recent years. Aimed at producing high milk production cows, has affected daughter FPL since high milk production is associated with shorter FPL, probably due to health and fertility problems, as observed in the present study. These genetic tendencies differ from those reported for countries such as France, Germany and Canada where culling risk has declined over time for the daughters of Holstein sires. Most likely this is in response to inclusion of FPL in genetic improvement programs for many years^9-12.

Conclusions and implications

The estimated heritabilities observed in the present study (0.08 for logarithmic scale; 0.13 for original scale; 0.12 for effective scale and 0.10 for equivalent scale) indicate that this trait can be effectively integrated into genetic improvement programs, as has occurred in other Holstein populations. The Weibull distribution and associated regression model adequately represented length of functional productive life for Holstein cattle in Mexico and are useful in calculating heritability and predicting genetic values for this trait.