Estimates of genetic trends for daily milk yield and somatic cell score of primiparous Holstein cattle in Iran

ABSTRACT

Genetic and permanent environmental correlations between milk yield and somatic cell score (SCS) were estimated. The data set used in this analyse contained a total of 204,438 test-day (TD) observations of milk yield and SCS from 22,903 primiparous daughters of 554 dairy Holstein sires. A multiple-trait random regression TD model was employed in the analyses using the Gibbs sampling method. Obtained results showed that daily heritabilities were highest in the second half of lactation and varied from 0.086 to 0.250 for milk yield and 0.041 to 0.086 for SCS. Undesirable genetic relationships with milk yield were found for SCS with genetic correlations ranging from −0.093 to 0.319 that the highest was estimated between SCS days in milk (DIM) 305 and milk DIM 5. Based on these results, we can conclude that the cows with high yield potential in the initial stage of lactation had the highest genetic potential for SCS at the end of lactation. This emphasizes the need for including udder health traits in the breeding goal. However results of this investigation, based on genetic trends, indicate that more emphasis is paid to increase milk yield of Holstein cows in Iran with no or a little emphasis on SCC.

KEYWORDS:

• Genetic relationships
• Gibbs sampling
• random regression

1. Introduction

The major source of income for most Iranian dairy producers is the sale of milk. Accordingly, milk yields have been the main breeding objective in Iranian Holsteins with no or little emphasis on conformation and health. Although these traits are of primary economic importance, but many studies have shown that antagonistic phenotypic and genetic correlations of somatic cell score (SCS) trait with milk yield would lead to a genetic decline in cow health, if selection is for milk production only (Hortet et al. 1999; Juozaitiene et al. 2006). Gilmore (1977) found a positive correlation between the incidence of health disorders and milk production and reported that dairy cattle with the highest annual health costs were likely to have the highest production during first lactation. In addition to the negative impact on incidence for health problems, the increase in SCS was reported to have a detrimental effect on production, reproductive performance and decline in modern dairy cows in industry herds (Juozaitiene and Juozaitis 2005). Carlén et al. (2004) found a high genetic correlation between lactation average SCS and mastitis (mean 0.70).

In the dairy industry, mastitis is a highly undesirable disease that considerably affects overall profitability (Mrode and Swanson 1996). With incidence mastitis the mean production of the herd decreased because a greater proportion of the culling decisions are based on production and proportion of mature cows, which produce more milk than young cows, is decreased (Mrode and Swanson 2003). Fortunately improving resistance mastitis through breeding, due to the high genetic correlation with mastitis, can be based on indirect selection for low SCS (Mrode and Swanson 1996; Heringstad et al. 2000). If selection against mastitis is included in a total merit index, the genetic level of mastitis may be kept constant or even improved (Heringstad et al. 2000).

Genetic relationship between the milk yield and SCS vary with stage of lactation and with test-day (TD), so an analysis based on TD records should be more accurate (Haile-Mariam et al. 2001). This genetic relationship, if identified, would contribute to the genetic improvement of disease resistance. Therefore, the aim of this study was to evaluate the genetic relationships and annual genetic trends for daily milk yield and SCS of primiparous Holstein cattle in Iran.

2. Material and methods

2.1. Data set

Data were provided by the Animal Breeding Center of Iran. In Iran, where records are costly to obtain for farmers, somatic cell count (SCC) is recorded in some herds on a monthly basis. Data were restricted to primiparous cows, calving between 2004 and 2012, that met the following requirements: (1) age at first calving between 660 and 1000 d, (2) test-days from 5 to 305 days in milk (DIM), (3) the first TD was before 60 days after calving, (4) the number of TD records ranged from 8 to 10. All records belonging to the cows whose sires (AI sires) had less than 10 daughters were also omitted; and also herd-year of calving subclasses were required to have a minimum of 10 cows. The SCC was transformed to SCS [SCS = log₂ (SCC/10⁵) + 3] to achieve an approximate normal distribution of the TD records (Ali and Shook 1980). Finally, edited data set consisted of 204,438 records, produced by 22,903 Holstein cows with known sire (Table 1) in 159 herds. The average number of daughters per sire was 41, and about 10% of the sires had more than 100 daughters. Summary statistics of the data set are presented in Table 2. Generally, the Legendre polynomials have largely been used to fit random curves due to their ability of describing the variation along the period (Bignardi et al. 2009). In this research, third-degree (with four coefficients) Legendre polynomials (Kirkpatrick et al. 1990) were used for both fixed and random regressions (RR) on the scale from 5 to 305 DIM.

Table 1. Number of sires and daughters by class of daughters per sire.

Class of number of daughters per sire	Number of sires	Total number of daughters
10–19	235	3177
20–49	172	5121
50–99	91	6339
100–199	48	6270
>200	8	1996
Sum	554	22,903

Table 2. Number of observations, arithmetic means, standard deviations (in parenthesis) for milk yield, SCS, and DIM by test number.

Test no.	Number	DIM, d	Milk, kg	SCC	SCS
5–35	19,570	20.2 (8.9)	29.0 (6.81)	230 (655)	2.60 (1.89)
36–65	21,251	50.1 (8.5)	33.0 (6.62)	222 (1141)	2.21 (1.97)
66–95	21,151	80.2 (8.5)	33.5 (6.51)	219 (1008)	2.22 (1.97)
96–125	21,356	110.3 (8.5)	33.0 (6.54)	217 (787)	2.30 (1.95)
126–155	21,236	140.3 (8.5)	32.3 (6.67)	212 (850)	2.33 (1.93)
156–185	21,412	170.3 (8.5)	31.5 (6.71)	219 (627)	2.39 (1.97)
186–215	21,290	200.4 (8.5)	30.7 (6.67)	221 (597)	2.46 (1.96)
216–245	21,311	230.4 (8.5)	29.6 (6.76)	239 (704)	2.54 (1.97)
246–275	20,048	260.3 (8.5)	28.3 (6.83)	231 (624)	2.58 (1.94)
276–305	15,813	289.9 (8.5)	27.2 (6.90)	248 (673)	2.70 (1.93)
5–305	204,438	152.2 (84.2)	30.9 (6.99)	225 (790)	2.42 (1.95)

2.2. Model

GLM procedure of SAS software (2003) was used for fitting the fixed effects and variables which were significant at P < 0.05 were included in the model. Finally, for multiple trait analyse, the animal model fitted was as follows:$y_{ijklmpq}=y{c}_{ji}+m{t}_{ki}+ht{d}_{li}+\sum _0^{3}A{S}_{mni}{\omega }_{qn}+\sum _0^3{\alpha }_{pni}{\omega }_{qn}\sum _0^{3}p{e}_{pni}{\omega }_{qn}+e_{ijklmpq},$where y_ijklmpq is the qth TD records of ith trait; y_cji is the fixed effect of the jth year of calving of ith trait (j = 2004, … , 2012); mt_ki is the fixed effect of kth milking frequency of ith trait (k = 2 and 3); htd_li is the fixed effect of lth Herd-year-month of TD of ith trait (l = 1, … , 5734); AS_mni is the fixed regression coefficient, nested within the mth class of the age-season of calving (m = 1, … , 12); ω_qn is the nth order polynomial corresponding to DIM at test; ɑ_pni is the nth random regression coefficients (RRC) for the additive genetic (AG) of pth animal (n = 1, … , 65,794); pe_pni is the nth RRC for the permanent environmental (PE) of pth animal (n = 1, … , 22,903); e_ijklmpq is the residual random effect associated with y_ijklmpq .

The covariance structure of the model is,$\mathrm{var}\left[\begin{array}{c}a\\ p\\ e\end{array}\right]=\left[\begin{array}{ccc}G\otimes K_a& 0& 0\\ 0& I\otimes K_p& 0\\ 0& 0& I\otimes E\end{array}\right],$where K_a is the 8 × 8 covariance matrix of the AG regression coefficients; G is the AG covariance matrix among all animals; K_p is the 8 × 8 covariance matrix of the PE regression coefficients; E is a 3 × 3 diagonal matrix of residual variances; I is the identity matrix and ⊗ is the Kronecker product. Variance and covariance for curve parameters were estimated with a Bayesian approach via the Gibbs sampling algorithm as implemented by Misztal (2015). A single chain with a length of 110,000 was generated discarding the first 10,000 samples (burn-in period). Heritabilities and genetic correlations (RAG) were calculated according to the following formula:$h_i^2={\displaystyle \frac{\mathrm{diagonal}(Q_{i}Ka{Q^{\mathrm{\prime }}}_i)}{\mathrm{diagonal}(Q_{i}Ka{Q^{\mathrm{\prime }}}_i)+\mathrm{diagonal}(Q_{i}Kpe{Q^{\mathrm{\prime }}}_i)+{\sigma }_e^2}}$ $h_{(305d)}^2={\displaystyle \frac{Q_{305}Ka{Q^{\mathrm{\prime }}}_{305}}{Q_{305}Ka{Q^{\mathrm{\prime }}}_{305}+Q_{305}Kpe{Q^{\mathrm{\prime }}}_{305}+301\ast {\sigma }_e^2}}$ $r_{g(i,j)}={\displaystyle \frac{(Q_{i}K{a}_{pq}{Q^{\mathrm{\prime }}}_j)}{\sqrt{\mathrm{diagonal}(Q_{i}K{a}_{pp}{Q^{\mathrm{\prime }}}_j)\ast \mathrm{diagonal}(Q_{i}K{a}_{qq}{Q^{\mathrm{\prime }}}_j}{)}^{\mathrm{\prime }}}}$

Ka_pq indicates a matrix of AG (co)variances between traits p and q. Note that i = j if observations were on the same day, and i ≠ j if observations were on different days. Genetic trends obtained by regressing the weighted average EBVs of the animals on the birth year of them.

3. Results

3.1. Fixed effects

Figure 1 shows that fitted lactation curve for SCS resembled inverted lactation curve for milk yield. The plot shows an increase in milk yield in early lactation, and was highest in peak of lactation (i.e. 60–90 d) and subsequently decreased toward the end of lactation. Over the period considered in the current investigation, milk yield showed a consistent decline; the opposite trend was observed for SCS. For instance, the year of calving solution for milk yield and SCS, respectively, for 2004 was 31.02 and 4.00 and that for the last period (2012) was 28.25 and 4.29. About the effect of age-season of calving, as it is expected, it is observed that with increasing age of calving, though negligible, but continuously these values showed an increasing trend. For instance, predicted milk yield for different classes/categorizations of calving age, regardless of calving season, they are obtained 4.51, 5.16 and 5.84, respectively. The corresponding values for SCS were −1.52, −1.45, and −1.41, respectively.

Figure 1. Average daily milk yield (solid line) and somatic cell score (dashed line).

3.2. Variances and heritabilities

Additive genetic and PE variance components for TD records of milk yield and SCS are presented in Figure 2. The additive genetic variance was increased as lactation progressed for milk yield and were highest in second half of lactation. For milk yield, it is obvious that the beginning and the end of lactation have higher PE variance than middle part. So, this curve was typically U-shaped. As indicated in these plats, AG variance for SCS was low at the peak of lactation, rose toward the middle of the second part of lactation, and then dropped to rise again at the end. With a few exceptions, PE of test day SCS was high in the beginning of lactation, low at the peak, and rose toward the end.

Figure 2. Daily genetic (AG) and permanent environmental (PE) variances for TD milk yield (solid line) and SCS (dashed line) over the lactation (variances for SCS are multiplied by 20).

Inspection of estimated heritability curves confirmed that the values were not constant throughout the lactation (Figure 3). Heritability estimates, derived based on the (co)variances of RRC, along the lactation trajectory showed shapes similar to the genetic variation but were less extreme at the beginning and the end of lactation because of higher PE variances. Heritabilities for the first TD generally were lowest. This finding was not only attributable to high PE variance, but also to a rather small estimate of the AG variance. Heritability estimates for TD records were highest in second half of lactation and ranged from 0.086 to 0.250 for milk yield and 0.041–0.086 for SCS. Table 3 indicates that estimates of heritabilities for TD records were smaller than the corresponding complete lactation (305-d).

Figure 3. Heritability for milk yield (solid line) and SCS (dashed line) over DIM.

Table 3. Estimates of variance components^a and the parameters derived from them^b for complete lactation (5–305 days).

	Yield
	Milk	SCS
${\sigma }_a^2$	542,063.08	11,219.04
${\sigma }_{pe}^2$	1,043,872.73	67,960.95
${\sigma }_e^2$	3156.59	577.19
$h_{(305d)}^2$	0.34	0.14
PE	0.66	0.85

^a ${\sigma }_a^2$, additive genetic variance; ${\sigma }_{pe}^2$, variance of permanent environment; ${\sigma }_e^2$, residual variance; $h_{(305d)}^2$, heritability estimate for complete lactation (5–305 days).

^bPE = Relative portion of permanent environment of total variance.

3.3. Genetic and PE correlations

Within traits: Estimated genetic and environment correlations, calculated from (co)variance components, between TD milk yield or TD SCS at different stages of lactation are shown in Figure 4. As expected, the pattern between genetic and environmental correlations is similar in that the correlations decrease as the distance between test days increase. As a consequence, genetic correlations obtained in this study were all positive (>0.30), even between yields at the extremes of the lactation trajectory.

Figure 4. Additive genetic (RAG) and permanent environmental (RPE) correlations between test-day milk yield (1st column), and test-day somatic cell scores (2nd column) at different stages of lactation for the same trait.

Between traits: Genetic and PE correlations between milk yield and SCS at different stages of lactation are presented in Table 4. Generally, the genetic correlations between these traits independent on how far apart they are and on whether they are on the same day in any two traits. For instance, the genetic correlation between milk yield DIM 305 and SCS DIM 305 is −0.006, while that between milk yield DIM 5 and SCS DIM 305 is 0.319. As a consequence, genetic correlation estimates between milk yield and SCS were small to moderate and ranged from −0.093 to 0.319. In the case of environmental correlation between these traits, our estimates were negative for most of the time and slightly higher than corresponding genetic correlations (ranged from −0.357 to 0.026).

Table 4. Additive genetic and PE correlations between milk yield and SCS for different days in milk.

SCS
RAG
	DIM	5	35	65	95	125	155	185	215	245	275	305
Milk	5	−0.09	0.07	0.18	0.22	0.23	0.23	0.22	0.21	0.23	0.27	0.32
	35	−0.03	0.14	0.23	0.24	0.25	0.23	0.20	0.17	0.16	0.19	0.24
	65	0.03	0.17	0.23	0.24	0.23	0.20	0.16	0.12	0.10	0.12	0.15
	95	0.08	0.19	0.23	0.22	0.20	0.17	0.13	0.09	0.07	0.07	0.09
	125	0.11	0.20	0.22	0.21	0.19	0.16	0.11	0.07	0.05	0.04	0.05
	155	0.14	0.21	0.22	0.20	0.17	0.14	0.10	0.06	0.04	0.02	0.03
	185	0.16	0.22	0.21	0.19	0.16	0.13	0.09	0.05	0.03	0.01	0.01
	215	0.17	0.22	0.21	0.18	0.15	0.12	0.08	0.04	0.02	0.00	−0.01
	245	0.16	0.21	0.20	0.17	0.15	0.11	0.07	0.03	0.00	−0.01	−0.01
	275	0.13	0.20	0.19	0.17	0.14	0.10	0.05	0.01	−0.01	−0.02	−0.01
	305	0.09	0.17	0.18	0.16	0.13	0.08	0.03	−0.01	−0.04	−0.04	−0.01
RPE
		5	35	65	95	125	155	185	215	245	275	305
Milk	5	−0.22	−0.13	−0.04	0.01	0.03	0.02	0.01	−0.01	−0.02	−0.02	−0.01
	35	−0.21	−0.20	−0.14	−0.09	−0.07	−0.05	−0.04	−0.04	−0.04	−0.05	−0.05
	65	–0.18	−0.23	−0.21	−0.18	−0.15	−0.12	−0.10	−0.07	−0.06	−0.07	−0.09
	95	−0.15	−0.24	−0.24	−0.22	−0.20	−0.17	−0.14	−0.11	−0.10	−0.10	−0.11
	125	−0.14	−0.24	−0.26	−0.25	−0.23	−0.21	−0.18	−0.16	−0.13	−0.13	−0.12
	155	−0.14	−0.24	−0.26	−0.26	−0.25	−0.24	−0.22	−0.20	−0.18	−0.16	−0.14
	185	−0.15	−0.23	−0.25	−0.25	−0.26	−0.26	−0.25	−0.23	−0.22	−0.20	−0.16
	215	−0.15	−0.22	−0.23	−0.24	−0.25	−0.26	−0.27	−0.26	−0.25	−0.23	−0.19
	245	−0.14	−0.21	−0.22	−0.23	−0.25	−0.26	−0.27	−0.27	−0.27	−0.27	−0.24
	275	−0.11	−0.19	−0.22	−0.23	−0.24	−0.25	−0.26	−0.27	−0.28	−0.30	−0.31
	305	−0.06	−0.16	−0.21	−0.22	−0.22	−0.22	−0.22	−0.23	−0.25	−0.30	−0.36

Note: RAG, additive genetic correlation; RPE, permanent environmental correlation.

Genetic trends

Phenotypic trend and estimated breeding values (EBV) of 305-day for milk yield and SCS by year of birth are presented in Figure 5 and the corresponding linear regression coefficients are given in Table 5. Phenotypic means and EBV of 305-day milk yield increased over the years. As shown in Figure 5, in all case the genetic trend for milk yield was positive and has risen about 550 kg in 2010, so that according to Table 5 it has been improved on average 50.45 kg/year during 2002–2010. This trend demonstrated the effectiveness of selection for improving milk yield. Phenotypic trends were in opposite directions for milk yield and SCS. Phenotypic means of SCS had a small variation and there was a decreasing trend from 2002 to 2010. On the other hand, there was a decreasing EBV of 305-day trend for SCS from 2004 to 2007, but there was no obvious variation in the concentrations of SCS thereafter. However, no significant genetic trend was observed for SCS.

Figure 5. Estimated of genetic (left) and phenotypic (right) trends for milk yield (solid line) and SCS (dashed line) of primiparous Holstein cows in Iran.

Table 5. Estimation of linear regression coefficients of animal estimated breeding values for 305-day milk yield and SCS.

Trait	Intercept ± SE^a	L^b ± SE	R ^2	Mean of EBV
Milk (kg)	−100,724* ± 11,243	50.45* ± 5.61	0.92	386.13
SCS	−961 ± 1218	0.48 ± 0.61	0.08	2.43

^aStandard errors of the estimate.

^bLinear regression coefficients (genetic trends).

*P < 0.05.

4. Discussion

Although the decline in milk yield under high SCC may be compensated by an increase of overall production, the decline might be economically more important if it involves traits correlated with milk. The economic loss of SCC may be larger than that measured by milk yield only; for instance, because of increased mastitis susceptibility (Heringstad et al. 2000) and decreased cheese yield and shelf-life of consumer milk (Schutz 1994).

Figure 2 appears that AG variances for milk yield as a function of DIM increased toward the end of lactation. Some authors found the highest AG for milk yield at the beginning of lactation (Jamrozik and Schaeffer 1997; Bignardi et al. 2009; Pereira et al. 2013) but other studies (Swalve 1995; Druet et al. 2003; Kheirabadi et al. 2014) reported the lowest variance in this period. As expected, the PE variance estimates were consistently larger than AG estimates for both traits. For milk yield, decreasing PE variance during lactation suggest that PE factor had the low influence on TD milk yields during lactation except at the beginning and the end of lactation (Kheirabadi et al. 2013). A tendency for nearly resembled inverted PE variance of SCS revealed that environments or the effect of stress components suppress the efficiency of breeding programmes, if the estimate is for SCS only. Therefore, combination of associated traits can help to identify an optimal breeding scheme.

Ranges of heritabilities for milk yield across the lactation varied remarkably among previous studies. For instance, much higher and very low estimates of heritability can be found in the literature for Holstein in Canada (Jamrozik and Schaeffer 1997) and Polish Black and White cattle (Strabel and Misztal 1999), respectively. Strabel and Misztal (1999) found a slightly lower heritability for milk yield, in the range of 0.13–0.17 and explained that this was because of the low production levels. The estimates of heritability obtained for milk yield in this study lie well in the range (0.05–0.30) reported in the literature (Yazgan et al. 2010; Pereira et al. 2013). In an earlier analysis of individual TD milk yield as separate traits in a RR model, heritabilities varied from 0.10 at DIM 5 to 0.21 at DIM 305 (Kheirabadi et al. 2013). The same pattern was also reported by Kheirabadi et al. (2014) with a multiple trait RR models on the same population used in this study.

Estimates of daily heritability, genetic and PE variances for SCS were very similar to those of Kheirabadi (2018). Daily SCS heritability was in agreement with other sources, which established that SCS was the least heritable at the beginning of lactation and most heritable in mid to late lactation (Haile-Mariam et al. 2001; Mrode and Swanson 2003; Ødegard et al. 2003; Yazgan et al. 2010). Kheirabadi and Alijani (2014) also estimated small values of SCS heritabilities, especially in the first lactation (0.03–0.07).

Heritability estimates for completed lactation of milk and SCS (0.34 and 0.14, respectively) are in the range of reported estimates from other studies using linear models. Estimates obtained by Welper and Freeman (1992) with Michigan data were 0.30 and 0.16 for the first lactation of milk and SCS, respectively. However, heritabilities of these traits for Iranian Holsteins were larger than those reported in primiparous Ontario Holstein cows (Reents et al. 1995). There are several factors that might influence level of heritability estimates. Among the many reasons which cause the differences in heritability estimates different models and different levels of yield are defined (De Veer and Van Vleck 1987).

Table 4 appears that genetic correlations between daily milk and SCS independent on the interval between tests, however, it can be important. Generally, the pattern of correlations was similar to those reported in other studies (Jamrozik et al. 2010; Samoré et al. 2008), from mostly positive value in the genetic correlation to all negative correlations for PE. In general, phenotypic correlations are negative, whereas genetic correlations display a different sign, depending on lactation number (Samoré et al. 2008). The positive correlation means that genes that positively influenced milk production had direct effect on SCS. From these results, genetic selection for production should be expected to cause unfavourable increase in SCS level in Iranian Holstein. Based on the presented results, we can conclude that the cows with high yield potential in the initial stage of lactation had the highest genetic potential for SCS at the end of lactation (Table 4). A tendency for the association of lower values of the corresponding PE correlation, as found in the present study, also revealed that poor environmental do not suppresses the proportion of the genetic correlation.

Results of this investigation, based on genetic trends, indicate that more emphasis is paid to increase milk yield with no or a little emphasis on SCS. In other words, increasing trend for milk yield over time in this study maybe indicated that Iranian dairy producers were successful in choosing progressively better semen and sires from the imported and local source over the years. However, the irregular trend genetic trend in SCS, maybe provide future justification for including SCS in the Iranian national economic selection index. Similar to the current results, Khorshidie et al. (2012) and Ghavi Hossein-Zadeh (2011) reported trend was positive for milk yield of Iranian Holstein population. Those authors reported upward trends from 2000 onwards.

5. Conclusion

As a consequence, genetic and PE correlations between traits tended to be favourable level (close zero for genetic and >−0.25 for environmental) towards the end of lactation that it was clear between DIM 245 and 305. The magnitude and trend of the heritabilities and correlations may be suggested the improvement in the level of SCC (reduction) could be achieved by basing selection on yield potential of animals at the end of lactation.

Disclosure statement

No potential conflict of interest was reported by the author.

ORCID

Khabat Kheirabadi http://orcid.org/0000-0002-9979-5346