Go to Top Go to Bottom
Gebreyohannes, Koonawootrittriron, Elzo, and Suwanasopee: Variance Components and Genetic Parameters for Milk Production and Lactation Pattern in an Ethiopian Multibreed Dairy Cattle Population


The objective of this study was to estimate variance components and genetic parameters for lactation milk yield (LY), lactation length (LL), average milk yield per day (YD), initial milk yield (IY), peak milk yield (PY), days to peak (DP) and parameters (ln(a) and c) of the modified incomplete gamma function (MIG) in an Ethiopian multibreed dairy cattle population. The dataset was composed of 5,507 lactation records collected from 1,639 cows in three locations (Bako, Debre Zeit and Holetta) in Ethiopia from 1977 to 2010. Parameters for MIG were obtained from regression analysis of monthly test-day milk data on days in milk. The cows were purebred (Bos indicus) Boran (B) and Horro (H) and their crosses with different fractions of Friesian (F), Jersey (J) and Simmental (S). There were 23 breed groups (B, H, and their crossbreds with F, J, and S) in the population. Fixed and mixed models were used to analyse the data. The fixed model considered herd-year-season, parity and breed group as fixed effects, and residual as random. The single and two-traits mixed animal repeatability models, considered the fixed effects of herd-year-season and parity subclasses, breed as a function of cow H, F, J, and S breed fractions and general heterosis as a function of heterozygosity, and the random additive animal, permanent environment, and residual effects. For the analysis of LY, LL was added as a fixed covariate to all models. Variance components and genetic parameters were estimated using average information restricted maximum likelihood procedures. The results indicated that all traits were affected (p<0.001) by the considered fixed effects. High grade B×F cows (3/16B 13/16F) had the highest least squares means (LSM) for LY (2,490±178.9 kg), IY (10.5±0.8 kg), PY (12.7±0.9 kg), YD (7.6±0.55 kg) and LL (361.4±31.2 d), while B cows had the lowest LSM values for these traits. The LSM of LY, IY, YD, and PY tended to increase from the first to the fifth parity. Single-trait analyses yielded low heritability (0.03±0.03 and 0.08±0.02) and repeatability (0.14±0.01 to 0.24±0.02) estimates for LL, DP and parameter c. Medium heritability (0.21±0.03 to 0.33±0.04) and repeatability (0.27±0.02 to 0.53±0.01) estimates were obtained for LY, IY, PY, YD and ln(a). Genetic correlations between LY, IY, PY, YD, ln(a), and LL ranged from 0.59 to 0.99. Spearman’s rank correlations between sire estimated breeding values for LY, LL, IY, PY, YD, ln(a) and c were positive (0.67 to 0.99, p<0.001). These results suggested that selection for IY, PY, YD, or LY would genetically improve lactation milk yield in this Ethiopian dairy cattle population.


The majority of Ethiopian cattle (e.g., Boran and Horro) are low producers of milk with high individual variability (Yohannes et al., 2002) which is attributed to natural selection to the tropical environment and management. Milk yield could be increased through within and between breed selection and crossbreeding with cattle from improved dairy breeds. The strategy currently employed in Ethiopia to increase milk yield is by crossbreeding locally adapted Bos indicus breeds capable of withstanding heat stress, diseases, low feed supply and low management level with Bos taurus dairy breeds of high genetic ability for milk yield. This mating scheme has resulted in a variety of crossbred groups with different levels of additive and non-additive genetic abilities.
Variance components and genetic parameters are needed for genetic improvement programs to predict the breeding values of candidates for genetic selection, to choose among mating plans and to predict selection response (Montaldo et al., 2012). Knowledge of the type and amount of genetic variation and distribution of animals for traits considered for selection in the population can help design optimum breeding programs (Willham and Pollak, 1985). Additive genetic variation has been instrumental for the genetic improvement of traits of economic importance in dairy cattle populations (Willham and Pollak, 1985).
Numerous estimates of genetic parameters for production traits have been reported (Kahi et al., 2000; Demeke et al., 2004; Espinoza et al., 2007). These parameter estimates are specific to particular populations. Estimates of genetic parameters could change due to selection, migration of genes from one population to another, or changing environmental conditions (Van Der Werf and De Boer, 1989). Parameter estimates for one population or breed may not be appropriate for use in another population. Limited research has been done to estimate variance components and genetic parameters for lactation pattern and milk production traits and parameters of a lactation curve function in Ethiopian dairy cattle populations. Thus, the objective of this study was to estimate variance components and genetic parameters (i.e., heritability, repeatability, genetic and phenotypic correlations) for lactation pattern and milk production traits and parameters of a lactation curve function in a multibreed Ethiopian dairy cattle population.


Description of the study area, animals and data

The study was based on dairy cattle milk data from Bako, Debre Zeit and Holetta Research Centers, Ethiopia. The Bako Agricultural Research Center is located 250 km West of Addis Ababa at an altitude of 1,650 m above sea level. The area receives a mean annual rainfall of 1,200 mm in a bimodal distribution, 80% of which falls from May to September. In this area, the average relative humidity is 59%, and the average minimum and maximum temperatures are 13.5°C and 27°C, respectively (Gebreyohannes et al., 2013).
The Debre Zeit Research Station of the International Livestock Research Institute (ILRI) is located approximately 50 km south East of Addis Ababa, in the Ethiopian highlands (9°N and 39°E), at an altitude of approximately 1,850 metres above sea level. The average annual rainfall in the Debre Zeit area is approximately 866 mm. The annual average temperature is 18.7°C and the average monthly relative humidity is 52.4% (Haile et al., 2011).
The Holetta Agricultural Research Centre is located 45 km west of Addis Ababa at 38.5°E longitude and 9.8°N latitude, and elevation of 2,400 metres above sea level. It is situated in the central highlands of Ethiopia. Average annual rainfall is approximately 1,200 mm. The annual average temperature is 18°C and the average monthly relative humidity is 60% (Demeke et al., 2004; Haile et al., 2011).
The herds in these three centers were genetically linked by the use of common sires supplied from the National Artificial Insemination Centre (NAIC). All of them used at least one common indigenous breed as dams of the replacement cows. Mating took place throughout the year using both artificial insemination and natural service. The bulls for natural services were those selected from the available males in the herd by visual appraisal. The Boran (B) and Horro (H) cattle (Bos indicus) used as a foundation stock for crossbreeding were bought from farmers in the Southern and Western Ethiopia, respectively. The herd from Debre Zeit research center was B and their crosses with Friesian (F).
Animals in this population were purebred B and H (Alberro and Hailemariam, 1982), and their crosses with F, Jersey (J) and Simmental (S) recorded by the Bako and Holetta research centers, and B and F×B crosses in Debre Zeit research station. This mating strategy generated 23 straightbred and crossbred groups (Table 1).

Lactation curve function and parameters

The modified incomplete gamma function (MIG; Papajcsik and Bodero, 1988) was chosen to be used in this study for its goodness of fit to the monthly test-day milk data and predictive ability of lactation milk yields (Gebreyohannes et al., 2013). The incomplete gamma (IG) is represented by yt = atbe−ct where t is the length of time since calving. The MIG assumes the parameter b of the IG function to be equal to one, i.e., yt = ate−ct. Here, the log-transformed linear form of MIG ((ln(yt/t) = ln(a)+(−ct)) was used to fit to the monthly test-day milk data using the regression procedure (SAS, 2003). The yt is milk yield (kg) at time t (d) after calving, and a, b and c are parameters of the functions. The parameter a is a scaling factor associated with the average yield, b is related to pre-peak curvature and c is related to post-peak curvature. Parameters of the MIG function were obtained as an output from the regression analysis fitting the MIG function to each lactation of each cow. Cows with lactation lengths shorter than 90 d (three monthly test-day records) were excluded from the analysis and longer lactations were truncated to 305-d lactation length. A total of 59,413 monthly test-day milk records were used to estimate the parameters of the MIG function.

Data analysis

The data covered a period from 1977 to 2010 for Bako and Holetta, and 1989 to 2006 for Debre Zeit research center. The analysis used 5,507 lactation milk records of 1,639 cows from 283 sires. The traits considered were lactation milk yield (LY, kg), initial milk yield (IY, kg), peak milk yield (PY, kg), average daily milk yield (YD, milk yield per day of lactation, kg), days to peak (DP, d), lactation length (LL, d), and parameters of the MIG function (i.e., ln(a) and c). Phenotypic means by breed group and parity were compared using least squares means (LSM). The model used to compute LSM included the effects of herd-year-season of calving, parity and breed group subclasses for all traits (LY, IY, PY, YD, DP, LL, ln(a) and c). In addition, the model for LY included LL as a covariate to account for LL differences among cows. Computations were done with the GLM procedure of the Statistical Analysis System (SAS, 2003). The model for IY, PY, YD, DP, LL, ln(a) and c was:
  • Yijkl = Observation on lth cow that calved in ith herd-year-season of calving, jth parity, kth breed group subclasses,

  • μ = Overall mean,

  • HYSi = ith calving herd-year-season subclasses (i = 1 to 325),

  • Pj = jth parity (j = 1, 2, 3, 4, 5, 6 and ≥7),

  • BGk= kth breed group (k = 1 to 23; Table 1), and

  • eijkl = residual associated with yijkl.

In addition to all the effects in the model above, the model for LY included the covariate term b LLijkl, where b = regression coefficient of LY on LL, and LLijkl = lactation length of the lth cow that calved in ith herd-yr-season of calving, jth parity, and from kth breed group. The LSM for breed groups and parities were compared using Bonferroni t-tests (SAS, 2003).

Estimation of variance components and genetic parameters

Variance components were estimated using an Average Information Restricted Maximum Likelihood (AI-REML) procedure with the ASREML software (Gilmour et al., 2009). Single and two-trait animal repeatability models composed of fixed (i.e., herd-year-season, parity and regressions on cow H, F, J, and S breed fraction and general heterosis) and random (i.e., animal additive genetic, permanent environment and residual) effects were used to estimate variance components and to predict breeding values. Lactation length was considered as a covariate only for LY analysis. The general heterosis of cow involved interbreed interactions between alleles of any two different breeds among the 5 breeds (H, B, F, J, and S) present in the population. The H, F, J, and S breed effects were computed as a deviation from the Boran. Eight single-trait (LY, IY, PY, YD, DP, LL, ln(a) and c) analyses were performed to estimate variances, heritability and repeatability, and predict breeding values. The single trait animal repeatability model in matrix notation could be described as:
  • y = the vector of observations (i.e., LY, IY, PY, YD, DP, LL, ln(a) and c),

  • b = vector of fixed effects of herd-year-season and parity subclasses, and LL covariate for LY only,

  • g = vector of fixed breed effects,

  • a = vector of random animal additive genetic effects,

  • pe = vector of random permanent environmental effects,

  • e = vector of random residual effects,

  • X = incidence matrices relating records to fixed herd-year-season and parity, and LL covariate for LY only,

  • Q = matrix relating observations to breed effect (through H, F, J, and S breed fractions) and general heterosis (through heterozygosity fraction),

  • Z = incidence matrices relating records to animal additive genetic effects, and

  • W = incidence matrices relating records to permanent environmental effects.

The model assumed an expected value of y to be Xb+ Qg. The vector of direct animal additive genetic effects was assumed to have a normal distribution with mean zero and variance V(a)=Aσa2, where A is the additive numerator relationship matrix among animals in the population, and σa2 is the additive genetic variance. The vectors of permanent environmental and random residual effects were assumed to have independent normal distributions with mean zero and variances v(pe)=Iσpe2 and v(e)=Iσe2 where I is the identity matrix, and σpe2 and σe2 are the permanent environment and residual variances, respectively. The phenotypic variance was σy2=σa2+σpe2+σe2. Heritabilities ( σa2/σy2) and repeatabilities (( σa2+σpe2)/ σy2) were calculated from the animal additive genetic, permanent environmental and residual variances, which were estimated in the single-trait analyses. Genetic and phenotypic correlations between traits (LY, IY, PY, YD, DP, LL, ln(a) and c) were computed using two-trait animal repeatability models. The two-trait animal repeatability model in matrix notation could be described as:
where yi is vector of observation for the ith trait; bi is vector of fixed effects (i.e., HYS, parity, and LL covariate for LY only), gi is vector of cow breed and general heterosis effects for the ith trait; ai is vector of animal direct genetic effects for the ith trait; Xi, Zi, and Wi are incidence matrices that relate fixed (herd-year-season and parity subclasses and covariate), animal direct genetic and permanent environmental effects to observations for the ith trait, respectively and Qi is matrix relating observations to breed effect (through H, F, J, and S breed fractions) and general heterosis (through heterozygosity fraction). The model assumed the expected value for trait i (yi) to be Xibi+Qigi. The vectors of direct animal additive genetic effects, permanent environment and residual effects for each trait were assumed to be normally distributed with mean zero. The variance-covariance structure of the random effects for a bivariate animal repeatability model could be described as:
where σa12 and σa22 are variances due to additive genetic effects for traits 1 and 2, σpe12 and σpe22 are variances due to permanent environmental effects for traits 1 and 2, σe12 and σe22 are residual variances for traits 1 and 2, σa1a2, σpe1pe2 and σe1e2 are additive, permanent environment and residual covariances between traits 1 and 2, A is the numerator relationship matrix, and I is an identity matrix.

Estimation of breeding values and sire rank correlations

Sire estimated breeding values (EBV) were computed from the Best Linear Unbiased Prediction (BLUP) solutions obtained from ASREML (Gilmour et al., 2009). The EBV for an animal was calculated as the sum of the generalized least squares (GLS) estimate of the difference between its breed group and B plus its predicted random additive genetic value. In matrix notation the EBV for each animal could be computed as:
Where is the EBV of an individual for a particular trait, is a vector of GLS estimates of differences between breeds H, F, J, and S and breed B, p' is the transpose of the vector of fractions of H, F, J, and S breeds in an individual, and is the predicted random animal additive genetic value (i.e., from p' ; Elzo and Famula, 1985; Arnold et al., 1992; Koonawootrittriron et al., 2002). The ranks of the sires EBV for each trait were tested using Spearman’s Rank correlations for significant association between the ranks of the sires for pairs of traits using the CORR procedure of the SAS (SAS, 2003). Genetic trends across years were evaluated using a regression analysis of mean yearly EBV for each trait calculated from solutions for each herd-year-season contemporary group and regressed on year. The year data collection started, i.e., 1977, was considered as the base year.


Factors affecting the lactation pattern and milk production traits

Herd-year-season, parity and breed group subclasses were important factors (p<0.001) for all traits (LY, IY, PY, YD, DP, LL, ln(a) and c). Crossbred cows had higher least squares means (LSM) for LY, IY, PY, and YD, and longer LSM for LL than B and H breeds (Table 2). The 13/16F 3/16B crossbred cows had the highest LSM for LY (2,490±178.9 kg), IY (10.5±0.8 kg), PY (12.7±0.9 kg) and YD (7.6±0.55 kg) and the longest LL (361.4±31.2 days). Upgrading crossbred F×B or F×H from 50% to 87.5% F showed no significant increase in LSM for LY, IY, PY, and YD. However, the H cows yielded higher LSM for LY (1,210±37.9 kg vs 947±42.3 kg), IY (3.9±0.2 kg vs 1.0±0.2 kg), PY (5.5±0.2 kg vs 3.5±0.2 kg) and YD (2.8±0.12 kg vs 1.8±0.13 kg) and milked longer (234.9±6.5 kg vs 211.1±7.2 kg) than B cows. These differences among breed groups agreed with previous reports (Kahi et al., 2000; Demeke et al., 2004; Wolf et al., 2005) suggesting that biological differences among indigenous (B and H) and exotic parental populations (F, J, and S) may help explain differences in additive and non-additive genetic effects among breed groups in this study.
Upgrading from 50% to higher F fractions for both F×B and F×H crosses showed no significant differences for milk yields (Table 2). This outcome was different from Koonawootrittriron et al. (2012), who found that the LSM for milk yield tended to increase as Holstein fraction increased from less than 0.625 (4,075±109 kg) to 0.875 (4,285±47 kg) and then it decreased towards purebred Holstein (4,120±76 kg). This difference may be attributed to limitations in the provision of the necessary management that high percentage F cows require, and differences in adaptation to the tropical environment and feed availability. In addition, the preferential treatment given to improved breed groups relative to low milk yield indigenous cows may also have contributed to the higher performance observed in crossbred cows than that of B and H cows. The milk yield and lactation length for B and H cows reported in this study were higher than values reported for the same breeds by Demeke et al. (2004). This could be due to inclusion of cows with lactation lengths above 90 d in the present study compared to the 60 d considered in the previous study (Demeke et al., 2004). Among the three exotic sire breeds (F, J, and S), the F were superior followed by S and J for the considered traits.
First parity cows had longer DP (39.6±1.1 d) and LL (322.6±4.2 d) compared to the second and later (≥3) parity cows (p<0.001; Table 3) while LSM for LY, IY, PY, and YD were lower for the first parity cows followed by cows in the second parity compared to cows in later parity (≥3). Other authors have observed increased daily milk yields with parity number (Epaphras et al., 2004). The difference among parities could be associated with the maturity of the cows, high feed intake of mature cows than heifers, differences in partitioning of more feed to maintenance and growth at the expense of milk production in heifers than mature cows and early culling of cow with lower milk yield (Gradiz et al., 2009). Longer DP and LL observed in first parity cows compared to later parities could be associated to the better persistency of first lactation cows (Stanton et al., 1992).

Covariance component and genetic parameter estimates

Variance component estimates (Table 4) from single-trait analysis for LY, IY, PY, YD, DP, LL, ln(a), and c produced low heritability (0.03±0.03 to 0.08±0.02) estimates for LL, DP, and parameter c of MIG. Heritability estimates for LY, IY, PY, YD, LL, and ln(a) were medium (0.21±0.03 to 0.33±0.04). Repeatability estimates ranged from 0.21±0.02 to 0.53±0.01 for LY, IY, PY, YD, LL, ln(a), and c, and 0.14±0.01 for DP (Table 4).
The heritability estimate for LY found here was comparable to values reported by Albuquerque et al. (1996) and Mashhadi et al. (2008) for Holstein cows, but higher than estimates (0.24±0.04) from Demeke et al. (2004) for B and F×B crosses, Yilmaz et al. (2011) for Brown Swiss cattle (0.25) and Ayied et al. (2011) for Holstein cattle (range 0.05 to 0.21). However, it was lower than values obtained in a Thai multibreed dairy population (0.43±0.24; Seangjun et al., 2009) and in a Dutch Black and White dairy cattle population (0.48; Hoekstra et al., 1994). Further, the heritabilities estimated in this study were higher for IY and PY, but lower for LL and DP than the values reported for a multibreed Thai dairy population (0.23±0.18 for IY; 0.20± 0.18 for PY and 0.52±0.28 for DP; Seangjun et al., 2009). Heritability estimates are specific to a particular population (Van Der Werf and De Boer, 1989). Heritability estimates vary with lactation (Albuquerque et al., 1996; Cilek and Sahin, 2009). Cilek and Sahin (2009) attributed the lower heritability for all lactation milk yields than the heritabilities of first and second lactations to culling and selection that may have reduced the phenotypic variance of older cows compared to first-lactation cows. Hoekstra et al. (1994), for instance, attributed the high heritability estimate in their study to the high milk production. Ayied et al. (2011) compared heritability estimates for milk yield and lactation length using four methods of estimation (MINQUE, ANOVA, ML, and REML) and reported heritability estimates ranging from 0.07 (REML) to 0.21 (MINQUE) for milk yield and 0.02 to 0.09 for lactation length. Espinoza et al. (2007) compared heritability estimated for the first, second, third and fourth parities lactation milk yield using univariate analysis with animal repeatability model and reported lower heritability values (ranged 0.10 to 0.15) from the univariate analysis than from the animal repeatability model (0.16).
The heritability estimate for YD (0.33±0.04, Table 4) was comparable to the heritability estimate of 0.28 reported by Lopez-Villalobos and Spelman (2010). Atili and Khattab (2005) reported lower heritability estimates for lactation milk yield and lactation length for Holstein Friesian cows using animal repeatability model compared to heritability estimated from sire model. In general, moderate heritability estimates indicate sufficient additive genetic variance that would allow selection within the population. Thus, considering heritability estimates, LY, IY, PY, YD, and ln(a) could be improved through selection of either one or a combination of these traits in this dairy cattle population. The low heritability and repeatability estimates for LL, DP, and c indicate that these traits were largely affected by environmental factors, perhaps largely related to the management and feeding of the animals. The medium repeatability obtained for LY, IY, PY, and YD (0.39 to 0.53) in this study suggested that information from the first parity could be used for early prediction of EBV and cow selection. These early EBV would later be updated as data from later lactations became available.
The genetic correlations between LY, IY, PY, YD, LL, and ln(a) were high (range 0.59±0.15 to 0.99±0.03) and positive (Table 5). Conversely, the genetic correlations coefficients between DP and LY (−0.27±0.15), IY (−0.57± 0.12), PY (−0.41±0.15), YD (−0.31±0.16), ln(a)(−0.42±0.16) were negative. Cows with higher IY and PY tended to reach their peak yield earlier than cows with lower IY and PY. The high and positive genetic correlations obtained here between LY, IY, PY, YD, ln(a), and LL indicated that cows with high genetic potential for milk production traits (LY, IY, PY, YD, and ln(a)) tended to show longer LL. Thus, selection for longer LL would tend to improve lactation milk yield in this Ethiopian population. To improve LY in this Ethiopian population, selection should consider sires and dams with higher IY, PY, and YD and longer LL. Similar correlations were also reported in different studies (Batra et al., 1987; Varona et al., 1998; Tekerli et al., 2000; Gradiz et al., 2009). Cows producing at higher levels in the beginning of lactation tended to have the highest peak yield and consequently the highest 305-d yield (Ali and Schaeffer, 1987; Tekerli et al., 2000).

Estimated breeding values and genetic trends

Estimates for general heterosis were 489±44.9 kg for LY, 2.6±0.2 kg for IY, 3.0±0.2 kg for PY, 1.7±0.1 kg for YD, −4.7±1.7 d for DP, 50.9±6.8 d for LL, 0.33±0.03 kg for ln(a), and 0.002±0.0003 kg/d for c. Table 6 show estimates of breed differences between H, F, J, and S relative to Boran. These values were used to estimate breeding values. The estimates of breed effects for F, S, and J relative to B were large and influenced the ranking of sires across breed groups for all traits. Differences between F, S, and J breeds and B were positive for all traits, except for DP that were negative. Sire EBV ranged from −392.5 to 1,962 kg for LY, −2.9 to 10.9 kg for IY, −2.5 to 10.6 kg for PY, −1.8 to 7.4 kg for YD, −26.3 to 139.4 d for LL, −28.7 to 10.1 d for DP, −0.495 to 1.301 kg for ln(a), and −0.002 to 0.01 kg/d for c, respectively (Table 7). Comparison of EBV across sire breed groups showed that the highest EBV were observed for straightbred exotic sires (Friesian, Jersey and Simmental) for LY, IY, PY, YD, LL, and ln(a). Spearman’s rank correlations between sire EBV were significant for all pairs of traits (p<0.001) except for the correlation between DP and parameter c (Table 8). The rank correlation for pairs of traits involving LY, LL, IY, PY, YD, ln(a), and c, were positive and significant which indicated that sires with higher EBV for LY would likely also have higher EBV for the other positively correlated traits. Conversely, the rank correlation between DP with LY, IY, PY, and YD was negative which could be explained by the low negative genetic correlation (range −0.27±0.15 to −0.57±0.12) observed in this study (Table 5).
The only trait with a significant regression of EBV on calving year was LL. The EBV for LL declined across calving years (b = −0.91; p = 0.02). In contrast, the EBV for LY (corrected for LL) across years was highly variable (Figure 1) perhaps preventing the upward trend of EBV yearly means becoming significant. The high degree of variation of EBV for LY across years may have been influenced by genetic variability among sires used across years and by management variability across years within and across the three locations involved in this study.


Crossbreeding and upgrading could help increase milk production in this Ethiopian dairy population. The high additive genetic variance, heritability and repeatability estimates obtained here suggested the possibility of improving LY, IY, PY, and YD by genetic selection. High genetic and phenotypic correlation estimates among the traits (LY, IY, PY, and YD) also suggested that genetic improvement in one of these traits would result in improvement of other genetic correlated traits. The medium to high repeatability estimates suggested that first lactations could be used for early prediction of animal genetic abilities in this population. These early predictions will later be updated as information from later lactations becomes available.


The authors would like to acknowledge the Rural Capacity Building Project of the Ethiopian Ministry of Agriculture and Rural Development for the financial support and the Tigray Agricultural Research Institute for the study leave offered to the first author. The study is based on data collected in the Holetta Research Center of the Ethiopian Institute of Agriculture Research, Bako research center of the Oromia Agricultural Research Institute and Debre Zeit research station of the International Livestock Research Institute (ILRI). The authors are grateful to all institutions.

Figure 1
Yearly mean estimated breeding value (EBV) of sires for lactation milk yield (LY, kg).
Table 1
Number of observations in each breed group
Breed group number Breed group1 Number of records
1 Boran (B) 222
2 1/4F 3/4B 38
3 1/2F 1/2B 1,857
4 9/16F 7/16B 15
5 5/8F 3/8B 183
6 3/4F 1/4B 375
7 13/16F 3/16B 10
8 7/8F 1/8B 43
9 1/2J 1/2B 642
10 9/16J 7/16B 16
11 5/8J 3/8B 37
12 3/4J 1/4B 115
13 1/2S 1/2B 310
14 Horro (H) 362
15 1/4F 3/4H 84
16 3/8F 5/8H 17
17 1/2F 1/2H 409
18 3/4F 1/4H 64
19 1/4J 3/4H 11
20 1/2J 1/2H 372
21 3/4J 1/4H 64
22 1/2S 1/2H 240
23 3/4S 1/4H 21
All 5,507

1 F = Friesian, J = Jersey, and S = Simmental.

Table 2
Least squares means (±standard errors) for lactation pattern and milk production traits and lactation curve parameters1
Breed group2 LY (kg) IY (kg) PY (kg) YD (kg) DP (d) LL (d) ln(a) c
Pure breed
  Boran (B) 947±42.3e 1.0±0.2d 3.5±0.2e 1.8±0.13e 53.9±1.9a 211.1±7.2c −2.25±0.03d 0.016±0.0003c
  Horro (H) 1,210±37.9de 3.9±0.2c 5.5±0.2d 2.8±0.12d 29.1±1.7b 234.9±6.5bc −1.95±0.03cd 0.017±0.0003c
Two breed cross
  1/4 F 3/4B 1,476±95.2cd 5.9±0.4b 8.2±0.5bc 4.3±0.29bc 33.1±4.4b 279.0±16.6abc −1.51±0.07ab 0.015±0.0007bc
  1/2F 1/2B 2,031±20.9a 8.4±0.1a 11.0±0.1a 6.4±0.06a 32.6±1.0b 337.2±3.6a −1.32±0.02a 0.012±0.0001a
  9/16F 7/16B 1,682±147.4bcd 7.1±0.7b 9.8±0.7ab 5.6±0.46ab 26.0±6.8b 323.3±25.7a −1.33±0.11a 0.013±0.0011abc
  5/8F 3/8B 2,134±52.1a 8.8±0.2a 11.4±0.2a 6.6±0.16a 31.5±2.4b 329.1±9.1a −1.28±0.04a 0.012±0.0004ab
  3/4F 1/4B 2,240±35.9a 8.9±0.2a 11.6±0.2a 7.0±0.11a 32.3±1.6b 343.2±6.3a −1.25±0.03a 0.011±0.0003a
  13/16F 3/16B 2,490±178.9a 10.5±0.8a 12.7±0.9a 7.6±0.55a 34.8±8.2b 361.4±31.2a −1.10±0.13a 0.012±0.0013ab
  7/8F 1/8B 2,410±92.4a 9.8±0.4a 12.1±0.4a 7.6±0.29a 24.8±4.2b 299.0±16.1ab −1.20±0.07a 0.012±0.0007ab
  1/2J 1/2B 1,788±26.5bc 7.0±0.1b 9.4±0.1b 5.6±0.08ab 33.3±1.2b 315.3±.6a −1.51±0.02ab 0.011±0.0002a
  9/16J 7/16B 1,496±141.7cd 5.4±0.7bc 7.8±0.7bcd 4.7±0.44bc 46.8±6.5ab 312.4±24.7a −1.69±0.10abc 0.011±0.0010a
  5/8J 3/8B 1,707±96.7bc 6.8±0.5b 9.0±0.5bc 5.4±0.30ab 26.3±4.4b 320.6±16.9a −1.54±0.07ab 0.012±0.0007ab
  3/4J 1/4B 1,832±56.0abc 6.8±0.3b 9.2±0.3b 5.7±0.17ab 37.5±2.6b 302.8±9.8a −1.57±0.04ab 0.011±0.0004a
  1/2S 1/2B 1,900±35.0ab 7.3±0.2ab 10.4±0.2a 6.0±0.11a 38.6±1.6b 339.1±6.1a −1.48±0.03a 0.011±0.0003a
  1/4F 3/4H 1,346±68.2 cde 5.2±0.3bc 7.0±0.3cd 3.7±0.21cd 30.4±3.1b 271.0±11.9abc −1.70±0.05bc 0.015±0.0005bc
  3/8F 5/8H 1,809±152.9abc 7.9±0.7ab 10.2±0.7ab 5.8±0.47ab 27.0±7.0b 324.5±26.7a −1.33±0.11a 0.013±0.0011abc
  1/2F 1/2H 1,836±31.6abc 7.1±0.1b 9.8±0.2ab 5.7±0.10ab 36.1±1.5b 321.0±5.5a −1.46±0.02a 0.012±0.0002ab
  3/4F 1/4H 2,184±72.8a 8.7±0.3a 11.5±0.3a 6.8±0.23a 38.2±3.3b 360.7±12.7a −1.31±0.05a 0.011±0.0005a
  1/4J 3/4H 1,404±174.3cde 5.4±0.8bc 7.8±0.8bcd 4.1±0.54bcd 34.6±8.0b 285.0±30.4abc −1.78±0.13bcd 0.012±0.0012ab
  1/2J 1/2H 1,621±33.1cd 6.5±0.2b 8.5±0.2bc 4.9±0.10b 32.5±1.5b 303.8±5.8a −1.60±0.02ab 0.012±0.0002ab
  3/4J 1/4H 1,724±73.9bc 7.4±0.3ab 9.7±0.4ab 5.5±0.23ab 29.4±3.4b 329.0±12.9a −1.53±0.05ab 0.011±0.0005a
  1/2S 1/2H 1,662±40.0cd 6.5±0.2b 8.8±0.2bc 5.2±0.12b 34.3±1.8b 299.2±7.0a −1.57±0.03ab 0.012±0.0003ab
  3/4S 1/4H 1,876±125.9abc 6.6±0.6b 10.1±0.6ab 5.9±0.39ab 47.0±5.8ab 309.3±22.0a −1.42±0.09a 0.012±0.0009ab

1 LY = Lactation milk yield; IY = Initial milk yield; PY = Peak milk yield; YD = Average yield per d; DP = Days to peak; LL = Lactation length; ln(a) and c = Parameters of the MIG function.

2 F = Friesian, J = Jersey and S = Simmental.

a,b,c,d Least squares means within a column group with different letters differ significantly (p<0.001).

Table 3
Least squares means (±standard errors) for lactation pattern and milk production traits and lactation curve parameters1
Parity N LY (kg) IY (kg) PY (kg) YD (kg) DP (d) LL (d) ln(a) c
1 1,361 1,443±24.2c 4.9±0.1d 7.5±0.1d 4.4±0.07d 39.6±1.1a 322.6±4.2a −1.79±0.02d 0.012±0.0002a
2 1,138 1,680±25.1b 6.4±0.1c 8.9±0.1c 5.1±0.08c 35.4±1.2b 317.2±4.4ab −1.59±0.02c 0.012±0.0002ab
3 917 1,805±26.9a 7.0±0.1b 9.5±0.1b 5.5±0.08b 35.5±1.2b 311.2±4.7ab −1.49±0.02b 0.013±0.0002bc
4 739 1,878±28.7a 7.5±0.1a 9.9±0.1ab 5.7±0.09ab 33.6±1.3bc 314.6±5.0ab −1.45±0.02ab 0.013±0.0002bc
5 540 1,892±31.8a 7.7±0.1a 10.1±0.2a 5.8±0.10a 31.0±1.5c 305.9±5.5bc −1.40±0.02a 0.013±0.0002cd
6 380 1,837±35.8a 7.3±0.2ab 9.7±0.2ab 5.6±0.11ab 34.9±1.6bc 302.3±6.2bc −1.45±0.03ab 0.013±0.0003cd
≥7 432 1,884±35.6a 7.6±0.2a 9.8±0.2ab 5.8±0.11ab 30.4±1.6c 290.8±6.2c −1.38±0.03a 0.014±0.0003d

1 LY = Lactation milk yield; IY = Initial milk yield; PY = Peak milk yield; YD = Average yield per d; DP = Days to peak; LL = Lactation length; ln(a) and c = parameters of the MIG function.

a,b,c,d Least squares means within a column group with different letters differ significantly (p<0.001).

Table 4
Variance component and genetic parameter estimates for lactation pattern and milk production traits and lactation curve parameters
Trait1 Variance components2
Genetic parameters3
σa2 σpe2 σe2 σp2 h2 r
LY (kg) 84,166 71,259 139,090 294,514 0.29±0.04 0.53±0.01
IY (kg) 2.2 0.5 4.1 6.8 0.33±0.04 0.39±0.02
PY (kg) 2.1 1.2 3.4 6.8 0.32±0.04 0.49±0.02
YD (kg) 0.9 0.6 1.4 2.9 0.30±0.04 0.51±0.02
DP (d) 52 38 557 647 0.08±0.02 0.14±0.01
LL (d) 547 1,395 7,358 9,300 0.06±0.03 0.21±0.02
ln(a) (kg) 0.04 0.01 0.12 0.16 0.21±0.03 0.27±0.02
c (kg/d) 0.0000004 0.0000034 0.0000120 0.0000157 0.03±0.03 0.24±0.02

1 LY = Lactation milk yield; IY = Initial milk yield; PY = Peak milk yield; YD = Average yield per d; DP = Days to peak; LL = Lactation length; ln(a) and c = parameters of the MIG function.

2 σa2 = Additive genetic variance, σe2 = Residual variance, σpe2 = Permanent environmental variance, σp2 = Phenotypic variance.

3 h2 = Heritability, r = Repeatability.

Table 5
Genetic correlations (above diagonal) and phenotypic correlations (below diagonal) for lactation pattern and milk production traits
Traits Traits1
LY IY PY YD DP LL ln(a) c
LY 0.91±0.03 0.96±0.02 0.98±0.01 −0.27±0.15 0.73±0.12 0.93±0.05 −0.02±0.29
IY 0.63±0.01 0.96±0.02 0.92±0.03 −0.57±0.12 0.59±0.15 0.90±0.05 −0.09±0.28
PY 0.78±0.01 0.80±0.01 0.97±0.02 −0.41±0.15 0.61±0.16 0.99±0.03 −0.04±0.28
YD 0.88±0.01 0.64±0.01 0.81±0.01 −0.31±0.16 0.67±0.16 0.96±0.04 0.13±0.25
DP 0.04±0.02 −0.30±0.01 −0.10±0.02 0.04±0.02 ** −0.42±0.16 **
LL 0.50±0.01 0.19±0.02 0.21±0.02 0.08±0.02 0.04±0.01 0.67±0.24 −0.09±0.77
ln(a) 0.57±0.01 0.53±0.01 0.65±0.01 0.63±0.01 0.09±0.01 −0.13±0.02 0.001±0.27
c 0.15±0.02 0.08±0.02 0.11±0.02 0.16±0.02 0.15±0.01 0.54±0.01 −0.43±0.01

1 LY = Lactation milk yield; IY = Initial milk yield; PY = Peak milk yield; YD = Average yield per d; DP = Days to peak; LL = Lactation length; ln(a) and c = parameters of the MIG function.

** Genetic correlation between DP and LL and DP and c could not be estimated.

Table 6
Estimates (±standard errors) for breed differences from Boran for lactation pattern and milk production traits
Trait Breed
Horro Friesian Jersey Simmental
LY (kg) 10.9±69.5 1,134±89.3 616±98.3 707±128.9
IY (kg) 0.6±0.3 6.3±0.4 3.9±0.5 3.8±0.6
PY (kg) −0.1±0.3 6.5±0.4 3.5±0.5 4.3±0.6
YD (kg) 0.1±0.2 4.5±0.3 2.9±0.3 3.4±0.4
DP (d) −10.0±2.5 −11.5±3.2 −9.8±3.5 −0.9±4.5
LL (d) −13.7±9.8 115.7±12.6 71.4±13.9 91.4±17.8
ln(a) (kg) 0.01±0.05 0.81±0.06 0.42±0.07 0.48±0.08
c (kg/d) −0.001±0.0004 0.01±0.001 0.01±0.001 0.01±0.001

1 LY = Lactation milk yield; IY = Initial milk yield; PY = Peak milk yield; YD = Average yield per d; DP = Days to peak; LL = Lactation length; ln(a) and c = Parameters of the MIG function.

Table 7
Number of sires (N) and ranges (minimum and maximum) of estimated breeding values of sires for lactation pattern and milk production traits across different sire breed groups
Sire breed group N Traits
Min Max Min Max Min Max Min Max Min Max Min Max Min Max Min Max
B 21 −393 196 −2.9 0.9 −2.5 1.0 −1.8 0.6 −4.5 10.1 −17.4 16.2 −0.50 0.17 −0.0003 0.001
H 28 −43 277 0.0 2.9 −0.4 1.6 −0.3 1.1 −22.0 −3.5 −26.3 5.5 −0.10 0.25 −0.0017 −0.001
F 77 729 1,962 3.5 10.9 3.7 10.6 3.0 7.4 −28.7 0.2 94.5 139.4 0.52 1.30 0.0046 0.005
J 32 505 1,069 2.7 5.4 2.5 6.2 2.3 4.3 −15.1 −4.5 47.4 114.8 0.32 0.74 0.0047 0.005
S 37 503 1,287 2.8 5.1 3.2 6.5 2.8 4.9 −6.9 2.8 62.4 108.5 0.34 0.67 0.0056 0.006
1/2F 1/2B 31 150 709 1.3 4.5 1.5 4.1 0.9 2.9 −13.4 −0.7 36.9 77.2 0.26 0.53 0.0022 0.003
3/4F 1/4B 4 718 1,006 4.3 5.6 4.4 5.7 2.9 4.1 −11.2 −6.6 72.0 90.5 0.56 0.70 0.0036 0.004
1/2J 1/2B 17 144 333 0.6 2.6 1.0 2.5 1.1 1.7 −8.8 −2.3 26.1 61.7 0.07 0.34 0.0023 0.003
3/4J 1/4B 2 432 564 3.1 3.1 2.6 3.1 2.2 2.5 −8.6 −4.4 52.5 55.6 0.29 0.40 0.0037 0.004
1/2S 1/2B 3 236 385 1.7 3.3 1.6 3.4 1.3 2.3 −1.7 0.9 38.3 57.1 0.14 0.46 0.0029 0.003
1/2F 1/2H 11 213 694 2.4 4.1 1.2 4.1 0.9 2.9 −16.8 −8.1 38.4 61.3 0.14 0.52 0.0016 0.002
3/4F 1/4H 2 899 967 5.3 5.3 5.1 5.4 3.6 3.8 −12.4 −9.8 77.9 94.5 0.68 0.69 0.0034 0.004
1/2J 1/2H 11 113 567 1.3 3.2 0.3 2.7 0.7 2.0 −12.5 −4.9 18.8 41.6 −0.01 0.36 0.0019 0.002
3/4J 1/4H 1 496 496 3.2 3.2 2.7 2.7 2.5 2.5 −8.4 −8.4 59.8 59.8 0.32 0.32 0.0036 0.004
1/2S 1/2H 6 37 463 1.1 2.3 0.3 2.6 1.1 2.2 −13.0 −0.8 25.8 50.7 0.06 0.33 0.0022 0.003

1 LY = Lactation milk yield; IY = Initial milk yield; PY = Peak milk yield; YD = Average yield per d; DP = Days to peak; LL = Lactation length; ln(a) and c = Parameters of the MIG function.

Table 8
Spearman’s Rank Correlation coefficients of sires estimated breeding values for lactation pattern and milk production traits and lactation curve function parameters
Traits Trait1
IY PY YD DP LL ln(a) c
LY 0.97* 0.99* 0.99* −0.33* 0.95* 0.96* 0.77*
IY 0.98* 0.96* −0.44* 0.93* 0.94* 0.73*
PY 0.99* −0.34* 0.95* 0.97* 0.76*
YD −0.30* 0.96* 0.95* 0.81*
DP −0.26* −0.37* 0.07ns
LL 0.92* 0.82*
ln(a) 0.67*

1 LY = Lactation milk yield; IY = Initial milk yield; PY = Peak milk yield; YD = Average yield per d; DP = Days to peak; LL = Lactation length; ln(a) and c = parameters of the MIG function.

* p<0.001, ns = Not significant.


Albuquerque LG, Keown JF, Van Vleck LD. 1996. Genetic parameters of milk yield, fat and protein yields in the first three lactations using an animal model and restricted maximum likelihood. Brazilian J Genet 19:79–86.

Ali TE, Schaeffer LR. 1987. Accounting for covariance among test-day milk yields in dairy cows. Can J Anim Sci 67:637–644.
Arnold JW, Bertrand JK, Benyshek LL. 1992. Animal model for genetic evaluation of multibreed data. J Anim Sci 70:3322–3332.
crossref pmid
Atil H, Khattab AS. 2005. Estimation of genetic trends for productive and reproductive traits of Holstein Friesian cows in Turkey. Pakistan J Bio Sci 8:202–205.
Batra TR. 1986. Comparison of two mathematical models in fitting lactation curves for pure line and cross line dairy cows. Can J Anim Sci 66:405–414.
Cilek S, Sahin E. 2009. Estimation of some genetic parameters (heritability and repeatability) for milk yield in the Anatolian population of Holstein cows. Archiva Zootechnica 12:57–64.

Demeke S, Neser FWC, Schoeman SJ. 2004. Estimates of genetic parameters for Boran, Friesian, and crosses of Friesian and Jersey with the Boran cattle in the tropical highlands of Ethiopia: milk production traits and cow weight. J Anim Breed Genet 121:163–175.
Elzo MA, Famula TR. 1985. Multibreed sire evaluation within a country. J Anim Sci 60:942–952.
Epaphras A, Karimuribo ED, Msellem SN. 2004. Effect of season and parity on lactation of crossbred Ayrshire cows reared under coastal tropical climate in Tanzania. Livest Res Rural Dev. 6:16http://www.cipav.org.co/lrrd/lrrd16/6/epap16042.htm consulted on 10/06/2012)

Espinoza AP, Villavicencio JLE, Gonzalez-Pena D, Iglesias DG, Pena RLDL, Almeida FR. 2007. Estimation of covariance components for the first four lactations in Holstein cattle according to different models. Zootecnia Trop 25:9–18.

Gebreyohannes G, Koonawootrittriron S, Elzo MA, Suwanasopee T. 2013. Fitness of lactation curve functions to daily and monthly test-day milk data in an Ethiopian dairy cattle population. Kasetsart J Nat Sci 47:60–73.

Gilmour AR, Cullis BR, Welham SJ, Thompson R. 2009. ASREML Discovery Reference Manual. NSW Agric.; Australia:

Gradiz L, Alvarado L, Kahi AK, Hirooka H. 2009. Fit of Wood's function to daily milk records and estimation of environmental and additive and non-additive genetic effects on lactation curve and lactation parameters of crossbred dual purpose cattle. Livest Sci 124:321–329.
Haile A, Joshi BK, Ayalew W, Tegegne A, Singh A. 2011. Genetic evaluation of Ethiopian Boran cattle and their crosses with Holstein Friesian for growth performance in central Ethiopia. J Anim Breed Genet 128:133–140.
crossref pmid
Hoekstra J, Van Der Lugt AW, Van Der Werf JHJ, Ouweltjes W. 1994. Genetic and phenotypic parameters for milk production and fertility traits in upgraded dairy cattle. Livest Prod Sci 40:225–232.
Kahi AK, Thorpe W, Nitter G, Baker RL. 2000. Crossbreeding for dairy production in the lowland tropics of Kenya I. Estimation of individual crossbreeding effects on milk production and reproductive traits and on cow live weight. Livest Prod Sci 63:39–54.
Koonawootrittriron S, Elzo MA, Tumwasorn S, Nithichai K. 2002. Estimation of covariance components and prediction of additive genetic effects for first lactation 305-d milk and fat yields in a Thai multibreed dairy population. Thai J Agric Sci 35:245–258.

Koonawootrittriron S, Yodklaew P, Elzo MA, Suwanasopee T. 2012. Association between milk production and Holstein fraction of upgraded dairy cattle in the Thai tropics. In : ADSA-AMPA-ASAS-CSAS-WSASAA Joint Annual Meeting; Phoenix, AZ. July 15–19, 2012.

Lopez-Villalobos N, Spelman RJ. 2010. Estimation of genetic and crossbreeding parameters for clinical mastitis, somatic cell score and daily yields of milk, fat and protein in new zealand dairy cattle. http://www.kongressband.de/wcgalp2010/assets/pdf/0534.pdf

Mashhadi MH, Kashan NEJ, Nassiry MR, Torshizi RV. 2008. Prediction breeding value and genetic parameter in Iranian Holstein bulls for milk production traits. Pak J Biol Sci 11:108–112.
crossref pmid
Montaldo HH, Castillo-Juarez H, Valencia-Posadas M, Cienfuegos-Rivas EG, Ruiz-Lopez FJ. 2010. Genetic and environmental parameters for milk production, udder health, and fertility traits in Mexican Holstein cows. J Dairy Sci 93:2168–2175.
crossref pmid
SAS. 2003. SAS OnlineDoc 9.1.3. SAS Institute Inc; Cary, NC, USA:

Seangjun A, Koonawootrittriron S, Elzo MA. 2009. Characterization of lactation patterns and milk yield in a multibreed dairy cattle population in the central Thailand. Kasetsart J Nat Sci 43:74–82.

Stanton TL, Jones LR, Everett RW, Kachman SD. 1992. Estimating milk, fat and protein lactation curve with a test-day model. J Dairy Sci 75:1691–1770.
crossref pmid
Tekerli M, Akinci Z, Dogan I, Akcan A. 2000. Factors affecting the shape of lactation curves of Holstein cows from the Balikesir Province of Turkey. J Dairy Sci 83:1381–1386.
crossref pmid
Van Der Werf JHJ, De Boer W. 1989. Estimation of genetic parameters in a crossbred population of Black and White dairy cattle. J Dairy Sci 72:2615–2623.
Varona L, Moreno C, Cortes LAG, Altarriba J. 1998. Bayesian analysis of Wood’s lactation curve for Spanish dairy cows. J Dairy Sci 81:1469–1478.
crossref pmid
Willham RL, Pollak E. 1985. Theory of heterosis. J Dairy Sci 68:2411–2417.
Wolf J, Zavadilova L, Nemcova E. 2005. Non-additive effects on milk production in Czech dairy cows. J Anim Breed Genet 122:332–339.
crossref pmid
Yilmaz I, Eyduran E, Kaygisiz A, Javed K. 2011. Estimates of genetic parameters for lactation shape parameters with multivariate statistical technique in Brown Swiss cattle. Int J Agric Biol 13:174–178.

Yohannes G, Zelalem Y, Gizachew B, Alemu GW, Sendros D. 2002. Milk yield and reproductive performance of Boran cows and growth rate of their calves under partial suckling method. p. 367–378. In : ESAP (Ethiopian Society of Animal Production) Proceedings of 9th Annual Conference of the ESAP held in Addis Ababa; Ethiopia. August 30–31; 2001.

Editorial Office
Asian-Australasian Association of Animal Production Societies(AAAP)
Room 708 Sammo Sporex, 23, Sillim-ro 59-gil, Gwanak-gu, Seoul 08776, Korea   
TEL : +82-2-888-6558    FAX : +82-2-888-6559   
E-mail : editor@animbiosci.org               

Copyright © 2024 by Asian-Australasian Association of Animal Production Societies.

Developed in M2PI

Close layer
prev next