Breeding beef cattle for an extended productive life: Evaluation of selection criteria in the Retinta breed

We characterized genetically the longevity and the productive life of beef cows considering three different traits: length of true life (LTL), the length of productive life (LPL) and a number of calvings (NC), in the Retinta breed. Data considered in the analysis were collected in 30 Spanish herds and consisted of 3187 cows born between 1993 and 2013: 1802 cows had complete information while 1385 cows had censored information. Pedigree information accounted for 7359 individuals and the genetic evaluation was performed under an animal model using the Weibull proportional hazard model for survival analysis. Results indicated that cow effects (age at first calving, herd, season-year of cow birth), as well as some calf effects (breed group, season-year of birth), were significant (P < 0.05) for all traits; however, sex of calves did not affect the traits. The low heritabilities obtained for LTL and LPL (0.14 ± 0.01 and 0.14 ± 0.01 respectively) and modest heritability for NC (0.30 ± 0.01), suggest that a higher response to selection can be expected for NC than for LTL and LPL.


Introduction
In general, beef breeding programs are focused on increasing live weight or short-term weight gain in order to produce more kilograms of meat per individual. However, it has hardly been studied how those criteria affect cow longevity. Cow longevity is an economically important trait directly related to the productive efficiency and meat yield in beef production systems (Van Melis et al., 2007). A shorter productive life of the cows forces to keep more heifers in the herd for replacement (Lopez de Maturana et al., 2007) and decreases the possibility of obtaining a higher number of calves during the cow productive life. As a result, the costly period from birth to first parity will be best amortized in cows with an increased longevity because the cost per weaned calf will decrease proportionally (Dákay et al., 2006). The inclusion of this trait in beef cattle breeding programs would allow selecting heifers with an increased chance of remaining in the herd for a longer period of time and, at the same time, farmers would be able to determine the age of culling depending on their breeding objectives (Caetano et al., 2013).
The cow longevity is a highly significant trait in Retinta, an autochthonous breed widely utilized in the Southwest of the Iberian peninsula. Nowadays, its census is slightly over 200,000 breeding cows, which are reared under the extensive regime of the Dehesa ecosystem, characterized by marginal pasturelands and dry and hot climate (Serrano et al., 1992). Although pure breeding is largely preserved in Retinta breed, the breed is also widely used in crossbreeding as a maternal line with other continental breeds, such as Charolais and Limousin, yielding highly profitable crossbred individuals. To maximize the herd profit under the extensive and low-input production systems, cows should have a long productive life to produce a large number of calves. It would be also necessary to know the economic optimum age to cull cows. From simulation studies, maximum herd lives ranging from 8 to 11 years (Rogers, 1972;Melton, 1980). It is generally assumed that the use of linear models in the genetic analysis of longevity traits is inadequate due to a violation of assumptions of normality (Lagakos, 1979). In the same way, random regression linear models, which are frequently used to analyze the length of productive life (LPL), cannot be fully considered since censored data and the non-linearity of the factors cannot be treated appropriately (Caraviello et al., 2004). Survival analysis is considered as the relevant statistical approach for the genetic evaluation of cow longevity (Olechnowicz et al., 2016). It has been widely performed in most livestock species, such as dairy cattle (Famula, 1981;Caraviello et al., 2004;Jenko et al., 2013), beef cattle (Rogers et al., 2004), pigs (Meszaros et al., 2010) and more recently in horses (Sole et al., 2017). Venot et al. (2013) carried out a study to incorporate productive life, through longevity and productive efficiency assessed by the number of calvings at a target age, in a new genetic evaluations of French beef cows. However, this analysis is still frequently overlooked in beef cattle breeding programs and despite their economic importance the number of applied studies performed to date is scarce (Rogers et al., 2004;Dákay et al., 2006;Szabó and Dákay, 2009). Similarly, there are few studies estimating genetic parameters for longevity traits in beef cattle (Tanida et al., 1988;Rogers et al., 2004;Dákay et al., 2006;Forabosco et al., 2006;Phocas and Ducrocq, 2006;Szabó and Dákay, 2009).
One of the reasons for this lack of information is the fact that a large number of animals is still alive during the genetic evaluation and therefore, only the lower bound of their eventual productive life is known (censored information). The exclusion of these records from the evaluation or their incorrect use as definitive values would lead to biased results (Raguž, 2014). This problem could be corrected using survival analysis (Ducrocq, 1997). This methodology combines the information on dead (uncensored) and alive (censored) individuals and enables a proper statistical treatment of censored records by accounting for nonlinear characteristics of longevity data (Vukasinovic, 1999).
The more accurate and widely used model in the genetic evaluations for LPL is the Weibull proportional hazard models (Sewalem et al., 2005;Mészáros et al., 2008;Bonetti et al., 2009;Strapáková et al., 2013). Due to a lack of exhaustive recording of reproduction records and cow presence inventory in beef cattle herds, the number of longevity studies of native breeds is very limited and integration of longevity as selection criteria into beef cattle breeding program is scarce.
Herein, we performed the genetic characterization of the productive life of Retinta breed cattle analyzing three different selection criteria as measures of longevity-related traits using a survival analysis based on the Weibull proportional hazard model.

Material and methods
We analyzed the reproductive data recorded in the selection nucleus of the Retinta breed. Data were collected between 1995 and 2016 by the National Retinta Breeders Association. A total of 31,639 calving records from 6212 cows belonging to 30 genetically linked herds were available. Genetic links were due to the use of artificial insemination sires that connected the 30 herds. Data were initially filtered following two standardization criteria: 1) only cows whose age at first calving ranged between 20 and 42 months were considered in the analysis, and 2) cows whose calving records fell outside the ± 3 SD range from the average calving interval were excluded. Cows with an incomplete reproductive history (age at the first birth unknown or some intermediate calving not controlled) were eliminated. After editing, 16,100 calving records from 3187 dams and 398 sires were considered in the analysis. Pedigree of the cows was extended to include all the available information in the stud book breed database (7359 animals). The cows were classified as presenting uncensored (dead) or censored (alive) records. Right-censored data (43.46%) included animals alive or with an unknown status at the time of analysis and less than 3 years since the last calving. Our dataset did not present left-truncated records because only cows born in or after 1993 with calvings started after 1994 (at the age of 2 or more years) were included. Beef cattle longevity was evaluated using three different productive traits: length of true life (LTL), estimated as the number of months between the birth of the cow and the birth of the last calf on record or censoring, LPL, estimated as the number of months between the birth of the first and the birth of the last calf on record or censoring (Meyer, 2009); and number of calving (NC), estimated as the number of calvings for each cow until its last calving or censoring time.
A preliminary analysis was performed using a Weibull hazard survival model to determine the significance of the effects included in the model using a likelihood random test, in which comparisons were performed between the full model and the same model but excluding one effect at a time. The effects evaluated were the age at first calving, herd, season-year of cow birth as time-independent effects, and age at calving, calf breed group, season-year of calving and the combination of sex and type of calving (simple or twin births) as time-dependent effects. Time-dependent effects account for effect changes at every time point (i.e., in each calving), while time-independent effects remain equal in each calving. Calving difficulty was not tested in this study since Retinta is a maternal breed in which calving problems are negligible: less than 0.2% of calvings exhibited some difficulties (MAPAMA, 2016).
Factors potentially associated with cow longevity were analyzed with the following model: where h0(t) is the baseline hazard function, agefirst i is the time-independent effect of the age at first calving (monthly intervals from 20 to 26.5; > 26.5 to 32.5; > 32.5 to 38.5; > 38.5); Herd is the time-independent effect of the herd (30 herds); cow*year*season k is the timeindependent effect of the season combined with the year of cow birth (79 classes); breed group l is the time-dependent effect of the calf breed group (1= Retinta, 2= Limousin x Retinta, 3= Charolais x Retinta, 4= cross with others breeds); year*season m is the time-dependent effect of the season combined with the year of calf birth (83 classes); age n is the time-independent effect of the age at calving (monthly intervals from 20 to 45; > 45 to 70; > 70 to 95; > 95 to 120; > 120 to 145; > 145 to 170; > 170 m.); g o is the random animal additive genetic value. Genetic parameters and expected breeding value (EBV) predictions were based on the Weibull proportional hazards method (Ducrocq et al., 1988). The three dependent variables (LTL, LPL, and NC) were examined with the same model described above. Such models are equivalent to a fully parametric model, in which the baseline is estimated at every discrete time point and the definition of hazard is modified to consider a discrete time scale (Ducrocq, 1999). A special case of proportional hazard models are the grouped data models, where failure times are grouped into intervals A i = [a i-1 , a i ), i = 1…r, with a 0 = 0, a r = +∞, and λ as the hazard function: where λ 0 (t) is the baseline hazard function and z i is the design vector for explanatory variables β for cow i. The last calving of the cow or censoring during the time interval A i is denoted as t i .
Rates of last calving risk were expressed via risk ratios, which denoted the relative risk of a cow in a certain class to end its productive life compared with cows belonging to a reference class with a risk ratio equal to 1.
The heritabilities (h 2 ) of the traits analyzed were estimated using the methodology described by Yazdi et al. (2002): where σ a 2 = genetic variance and P = proportion of uncensored records.
The reliability R 2 of each estimated breeding value (EBV) was based on the prediction error variance for each estimate as R where pev is the prediction error variance calculated as the square of the standard error for each EBV.
All the analyses were carried out using the Survival Kit v6 software (Meszaros et al., 2013) 3. Results Table 1 shows the distribution of censored and uncensored records in each trait studied (LTL, LPL and NC). Average censoring ages after first calving 87.62 were and 54.95 months for LTL and LPL respectively; an average number of calves produced per cow was 4.71 for NC. Average failure ages after first calving were 97.23 and 65.02 months for LTL and LPL, and 5.31 calvings for NC.
Likelihood ratio test results for each effect included in the model are presented in Table 2. Only the combination of sex and type of calving was not significant (P > 0.05) and therefore, excluded from the analysis. All other effects were highly significant (P < 0.001) except for age at first calving in LTL and NC. It is noteworthy that higher hazard ratios for the end of productive life were observed in cows which were older at first calving ( Fig. 1-a). Moreover, the highest risk (1.51) was for cows whose first calving age was above 38 months, while the lowest risk (0.42) was for cows starting to calve at early ages (> 20-26.5). The relative risks showed large differences among herds, which reflected the importance of this effect; these were from 0.16 to 2.78 for LTL, from 0.15 to 2.62 for LPL and from 0.19 to 2.31 for NC ( Fig. 1-b). There was also a wide variation in the relative risks from one to another cow birth seasonyear, from 0.01 to 9.65 for LTL, from 0.01 to 14.41 for LPL and from 0.01to 5.28 for CN ( Fig. 1-c). The highest risk for age at calving was for > 170 months for the three traits and for LTL also for 20-45 months ( Fig. 2-a). Regarding breed of calf, the Limousine-Retinta crossbreeding showed the highest risk of failure, and the Retinta pure breeding the lowest (Fig. 2-b). Relative risks based on season-year of calf birth ranged from 0.001 to 43.23 for LTL, from 0.001 to 365.83 for LPL and from 0.001to 1.82 for CN (Fig. 2-c).
Heritability, EBV averages, and relative risk averages were estimated for all the traits studied (Table 3). The EBV reliability frequency distribution of the cows, for the three traits follows a normal distribution with a slight right bias for CN (Fig. 3). Correlations between EBV of all the pedigree animals of the three traits analyzed were very high, 0.99 for LTL and LPL, 0.97 for LTL and NC and 0.96 for LPL and NC. The standardized average cow EBV trend over time is shown in Fig. 4. The standardization of genetic values was carried out using the average genetic value of the year 1994 (beginning of the breeding program) and dividing by the standard deviation of the global VGs. The result was multiplied by a coefficient so that the limits of this standardized VG oscillated between 80 and 120 with an average of 100. A positive genetic trend was found over the 20 years of the study. Table 4 showed positive but moderate correlations between EBV of the current growth selection (weight at 120 and 180 days) and EBV for the three longevity traits (LTL, LPL and NC) in Retinta breed.

Discussion
The lifespan of livestock females is divided into two well-differentiated parts: (i) the period from birth to first parity, in which their economic return is zero, and (ii) the productive period, from their first calving to their disposal. This second period, in which the individual effectively produces and generates an economic income, is commonly called longevity (Essl, 1998). In dairy cattle, LPL has been defined as the number of days from first calving to death or culling, having a substantial impact on the net returns of dairy production (Reinhardt, 1996). Ducrocq (1994) defined two types of LPL: true LPL is the ability of a cow to avoid culling no matter the reasons, being estimated as the period of elapsed time between first calving and disposal of the cow; functional LPL is the ability of a cow to avoid involuntary culling and corresponds to LPL corrected for the level of production. True LPL is strongly related to the production performance of the cow. On the contrary, functional LPL is more related to the ability of the cow to avoid non-productive causes of culling, such as infertility, mastitis, lameness or any other disease.
In the present study, we analyzed three different traits to determine their adjustment as longevity selection criteria to be included in a breeding program.
The accuracy of genetic evaluation by survival analysis is highly influenced by the ratio of censored and uncensored records available, being increased as the proportion of censored records decreases. In our case, 43% of the data employed was censored, fitting into the limits proposed by Vukasinovic et al. (1997), who recommended not much more than 40% of censored records to obtain accurate results. Estimations of genetic variance in survival models have been usually based on sire or sire-maternal grand sire models (Tarres et al., 2006). In our case, we adopted an animal model despite it was believed it could not fit with the Survival Kit because of this alleged poor performance of the Laplace approximation (Ducrocq, 2006). However, with an appropriate data structure, it has been shown that this concern was not justified (Ducrocq, 2006). In this sense, Meszaros et al. (2010) suggested that the use of these models would produce a lower estimation of the genetic variance and heritability. However, the same authors also suggested that the results obtained would be more accurate due to a better correction of environmental effects and be less contaminated by other genetic effects.  According to our results, longevity traits were heavily influenced by several factors associated with the productive cycle. Among them, the age of first calving showed that cows delivering for the first time at older ages (> 38 months) had higher hazard ratios of ending their productive life at early stages of their LPL. This observation agrees with Rogers et al. (2004), who demonstrated that beef cows younger than 24 months at first calving had fewer chances of being culled than cows delivering at older ages. Phocas and Ducrocq (2006) reported that Charolais cows having their first calving below 28 months of age had 10% less risk to be culled than heifers calving between 28 and 39 months. Conversely, Ducrocq (1994) and Vukasinovic et al. (2001), among others, suggested that the percentage of variation in longevity explained by the age at first calving for dairy cows was minimal, and therefore, removed this effect from the analysis. In our case, the age at first calving was non-significant for LTL and CN, however, it was significant for LPL. Therefore, the results obtained in this study should be taken into account in the modeling of these traits for including them in future breeding programs.
In our study, culling risk was highly influenced by the herd, the year-season of cow birth, and the year-season of calf birth. However, our data do not allow us to disaggregate the causes of the variability observed. In any case, there are large differences among herds that should be included in breeding programs. Additionally, a more detailed analysis, including nutritional and health-related factors, should be performed to determine the origin of these differences in order to adjust the management practices employed to improve the longevity of the cows.
Retinta cows were widely used as a maternal line in two-breed crosses with Charolais and Limousin bulls by farmers due to the high yield observed in the offspring. In this system, the average weaning weight (determined at 180 days) was 199.8 kg for pure Retinta, 206.4 kg for LimousinxRetinta (LR), 209.6 kg for CharolaisxRetinta (CR) crosses, and 185.1 kg in Retinta crosses with other breeds. However, LR and CR crosses increased the risk of premature culling by 32.1% and 15.9%, respectively. These results suggest that despite their decreased weaning weight, the breeding of pure Retinta individuals could increase the whole productivity of the system since the total number of calvings was increased. Similar results were previously observed by Venot et al. (2013), who found differences among the productive life among French beef cattle breeds. However, to our knowledge, this is the first study that compares pure and crossbreed individuals that have the same maternal breed.
Previous studies determined that heritability of longevity traits in beef cattle ranged between 0.04 and 0.21 (Sewalem et al., 2005;  and Ducrocq, 2006;Mészáros et al., 2008;Bonetti et al., 2009;Strapáková et al., 2013). In our case, heritability estimates were low for LTL and LPL (0.14). These results suggest that the genetic improvement of these traits will be difficult. This hypothesis was also suggested by Rogers et al. (2004) in crossbred cows and by Snelling et al. (1995) in Angus cattle, supporting our findings. Another interesting point is that EBV correlations among the three traits were higher than 0.96, showing that the traits are biologically similar. However, the h 2 determined in NC was moderate (0.30) and, therefore, this trait could be used as a selection criterion to obtain a faster genetic improvement of longevity. In addition, the trend of NC over cow age (relative productivity) could allow a close monitoring of the development of each individual throughout the productive life.
Finally, we also found a positive genetic trend in the traits analyzed over the last twenty years (Fig. 4). Cows born in recent years had a higher probability of achieving a longer life and, as a consequence, a higher probability of getting an increased number of calvings. This improvement can be due to the current selection for weaning weight  due to the existence of moderate positive correlations between weaning weight and longevity traits (Table 4). However, Martinez et al. (2004) found that selection for weights at weaning and yearling ages had a little genetic effect on any measure of lifetime production in Hereford cows.
One possible way to increase the reliability of breeding values is to combine longevity with other positively correlated traits measured early in life, such as conformation traits. So far, little is known about the relationship between longevity and conformation traits in beef cattle. Forabosco et al. (2004) found that muscularity could be used as an early predictor of longevity in Chianina breed. Van Melis et al. (2010) found that heifer pregnancy rate can be used as an early predictor of longevity as well as post-weaning weight and muscle score in the Nellore breed.

Conclusion
According to the results obtained herein, age at first calving, herd, and season-year of cow birth, as well as calf breed and year-season of birth were important factors affecting longevity in Retinta beef cattle. Producers can increase productive life by modifying these parameters. Since differences among herds are very large, it is important to take care of the reproductive, sanitary and feeding management and to choose the more appropriate breed for crossbreeding purposes.
The high heritability for NC indicates that the number of calvings is the more relevant selection criterion for a beneficial response to direct genetic selection for longevity. The high correlations among the breeding values of the three traits indicate that the same kind of biological improvement would occur by selecting for one or another. Nevertheless, we would recommend NC because of its higher heritability and reliability, and because the collection of data on the number of calvings is easier and more accurate.