Estimation of co-variance components, genetic parameters, and genetic trends of reproductive traits in community-based breeding program of Bonga sheep in Ethiopia

Objective The objectives of the study were to evaluate reproductive performance and selection response through genetic trend of community-based breeding programs (CBBPs) of Bonga sheep. Methods Reproduction traits data were collected between 2012 and 2018 from Bonga sheep CBBPs. Phenotypic performance was analyzed using the general linear model procedures of Statistical Analysis System. Genetic parameters were estimated by univariate animal model for age at first lambing (AFL) and repeatability models for lambing interval (LI), litter size (LS), and annual reproductive rate (ARR) traits using restricted maximum likelihood method of WOMBAT. For correlations bivariate animal model was used. Best model was chosen based on likelihood ratio test. The genetic trends were estimated by the weighted regression of the average breeding value of the animals on the year of birth/lambing. Results The overall least squares mean±standard error of AFL, LI, LS, and ARR were 375± 12.5, 284±9.9, 1.45±0.010, and 2.31±0.050, respectively. Direct heritability estimates for AFL, LI, LS, and ARR were 0.07±0.190, 0.06±0.120, 0.18±0.070, and 0.25±0.203, respectively. The low heritability for both AFL and LI showed that these traits respond little to selection programs but rather highly depend on animal management options. The annual genetic gains were −0.0281 days, −0.016 days, −0.0002 lambs and 0.0003 lambs for AFL, LI, LS, and ARR, respectively. Conclusion Implications of the result to future improvement programs were improving management of animals, conservation of prolific flocks and out scaling the CBBP to get better results.


INTRODUCTION
The total sheep population in Ethiopia was estimated at 31.30 million with 9 breeds [1]. Sheep production is one of the major livestock production systems in the South Western part of Ethiopia [2].The breeds make immense contributions which are both tangible and intangible in nature. Some of the tangible benefits of sheep are immediate cash income, meat, milk, skin, and manure while the intangible benefit includes social prestige among the community members. Moreover, sheep play great role in the economy of the country by being exportable items and thus, are sources of much needed foreign currency [3].
Bonga sheep, one of the wellknown and largest breed of Ethiopia, is characterized by long and wide fat tail with tapering and twisted end, both male and female are polled, short and smooth hair, mainly convex facial profile of male and predominantly light red coat color. The breed is known by its docile temperament, good fattening potential, fast growth and prolificacy [4].
To improve productivity of local breeds, crossbreeding based on imported sires has been followed for many years though with little success [5]. Therefore, improving locally adapted and diverse breeds through selective breeding has been considered as an option in developing countries to satisfy the growing demand for animal products [6]. Con sequently, communitybased breeding programs (CBBPs) were initiated in Ethiopia in 2009 by the International Centre for Agricultural Research in Dry Areas (ICARDA), the In ternational Livestock Research Institute (ILRI) and the Austrian University of Natural Resources and Life Sciences (BOKU) in collaboration with the National and Regional Agricultural Research Systems in Ethiopia involving four breeds of sheep (Bonga, Afar, Horro, and Menz) [6]. The breeding programs for all breeds except Afar are still being implemented with active participation of the communities keeping the animals. For Afar; flock mobility, very high tem perature, frequent droughts and poor infrastructure in the pastoral system; so far limits designing and implementation of community based breeding programs in pastoral areas [7]. Bonga CBBP being the most successful program in Ethiopia has been upscaled and additional CBBPs have been established [6].
Reproductive performances of sheep together with, survival and growth traits are important determinants of productivity [8]. The breeding objective traits identified by the community in Bonga CBBPs were growth rate, tail type, presence or ab sence of horn, twining rate, mothering ability and coat color [9]. For accurate genetic evaluation and selection, estimates of genetic parameters for traits of importance should be known [10]. In order to achieve the largest possible gains, a thorough evaluation of the program is needed [11]. Besides, it is a rel evant tool to show that the promised benefits for farmers can be achieved and that livestock breeding is a sustainable intervention strategy under CBBP.
As has been indicated above, Bonga CBBP started in 2009 and through the upscale program an additional 14 CBBP were established since 2012. Frequent updating the genetic parameter estimate through evaluation is important and an integral component of breeding program. Therefore, the objectives of this study were to evaluate the genetic trend for reproductive traits and to estimate genetic parameters which help for optimization of the programs.

Description of study area
The study was conducted in four districts, namely Adiyo, Gesha, ShishoEnde, and Tello, of Kaffa zone of Southern Nation Nationalities People Region, Ethiopia. The area is characterized by mixed croplivestock production system. It has one major rainy season that extends from May to Octo ber and a dry season that lasts from October to April [9]. The altitude range of study area is 1,600 to 3,348 meters above sea level and the minimum and maximum temperature was 14°C and 32°C with an average of 24°C. Similarly, the re corded minimum and maximum rain fall within the period of data 2012 to 2018 collections was 1,079 and 2,032 mm per year with an average of 1,617 mm/yr.

Breeding program and animal management
Among the 14 CBBPs, five (Abeta, Buta, Dacha, Dirbedo, and Shosha) were established in 2012 while the remaining were during 2014. Selection of male lambs is being carried at two stages: screening of heavy weaners at weaning (3month) followed by selection at six months of age (post weaning) by using their estimated breeding values (EBV). When the Bonga CBBP started, selection was carried out at six months. How ever, because of fast growth potential of the breed, it was noted that many lambs are sold before they reach the selection age. Therefore, the twostage selection was implemented to keep the best ram lambs within the flock. Candidate lambs which had horn, short tail or black coat color were culled regardless of their EBVs.
Flocks being indoors at night in pens made up of bamboo walls and by any locally available corrugated roofing materials whereas, some farmers kept their flock around homestead at night. Flocks were tethered especially adult male and females during crop cultivation period. Therefore, main feed resource was pasture but additionally crop residue and kitchen left overs were used. Feed availability and abundance vary with rainfall patterns. Comparatively huge amount of feed was available in the rain season whereas feed was less in both quality and quantity during the dry season.

Data collection
Performance and pedigree data for this study were collected from the breeding database of outscaled Bonga sheep CBBPs. The reproductive performance data used for this study in cluded age at first lambing (AFL), lambing interval (LI), annual reproductive rate (ARR), and litter size (LS). Annual reproductive rate refers the number of young produced per breeding ewe female per year and it was calculated as: ARR = (365×LS)/LI. The detailed data structure and number of records for studied traits was indicated in Table 1.

Data analysis
First, the data were checked for pedigree structure using pedigree viewer then Proc Univariate in Statistical Analysis System (SAS) [12] was employed prior to any analysis to check for analysis of variances assumptions. The phenotypic evaluation was done using the general linear model proce dures of the SAS fitting nongenetic factors like year of birth (7 levels: 2012 through 2018), type of birth (3 level: single, twin, triple, and above), season of birth (2 level: wet and dry), dam parity (7 levels; parity 1 to 6, and 7 and above) and CBBP cooperative (13 levels) as fixed effect for AFL. Whereas lamb ing year, lambing season, lambing type, CBBP cooperatives and lambing parity at the same level with AFL considered for LI, LS, and ARR. Significant least square means were separated using Adjusted TukeyKramer method in SAS. All significant fixed effects were included in the genetic anal ysis.
Co(variance) components, genetic parameters and breed ing values (EBVs) were estimated by restricted maximum likelihood fitting an animal model using WOMBAT software [13]. The following six univariate for AFL and repeatability for LI, LS, and ARR animal models were tested for each trait. The statistical models used were: Model 1: y = X b +Z 1a +e Model 2: y = X b +Z 1a +Z 3pe +e Model 3: y = X b +Z 1a +Z 2m +e with Cov(a,m) = 0 Model 4: y = X b +Z 1a +Z 2m +e with Cov(a,m) = Aσ am Model 5: y = X b +Z 1a +Z 2m +Z 3pe +e with Cov(a,m) = 0 Model 6: y = X b +Z 1a +Z 2m +Z 3pe +e with Cov(a,m) = Aσ am Where, y is n×1 vector of observations for each trait, b is a vector of fixed effects (year, parity, season, sex, CBBP coop erative and birth type), a, m, pe, and e are vector of random effects for direct additive genetic effects, maternal additive genetic effects, animal permanent environmental effect and residual effects, respectively. X, Z 1 , Z 2 , and Z 3 are the inci dence matrices of fixed effect, direct additive genetic effect, maternal genetic effect and animal permanent environmen tal effect for LI, LS, and ARR but permanent environmental effect of the dam for AFL. A is the numerator relationship matrix between animals, and σ am covariance between direct and maternal genetic effects. According to El Fadili et al [14] the (co)variance structure of the random effects were: , σ am , and σ e 2 are direct additive genetic variance, maternal additive genetic variance, animal permanent environmental variance, maternal permanent environmental variance, directmaternal genetic covariance, and residual variance, respectively. I d and I n are identity ma trices of an order equal to the number of dams and the number of lambs, respectively.
Estimates of additive direct (h 2 a ) and additive maternal (h 2 m ) heritability, ratio of animal permanent environmental variance with phenotypic variance (pe 2 ) and ratio of maternal permanent environmental variance with phenotypic vari ance (c 2 ) were calculated as ratios of estimates of additive direct (σ a 2 ), additive maternal (σ m 2 ), animal permanent envi ronmental (σ pe 2 ), and maternal permanent environmental (σ c 2 ) variances to the phenotypic variance (σ p 2 ), respectively. Total heritability was calculated according to the following equation [15]: The genetic correlation between direct and maternal ge netic effects (r am ) was estimated as the ratio of the estimates of the σ am to the product of the square roots of the estimates of σ 2 a and σ 2 m [16]. The genetic correlation (r g ) between traits were estimated as the ratio of the estimates of the genetic covariance between the traits 1 and 2 to the product of the square roots of the es timates of genetic variance for trait 1 and genetic variance for trait 2.
Genetic trends of the traits were estimated by regression of predicted breeding values on the birth year [17]. Genetic gain was calculated as the difference between the EBVs of last and first year of the program [18].
Repeatability (r) was estimated according to Mokhtari et al [19]: r = (σ 2 a +σ 2 pe )/σ 2 p Where, σ 2 a = additive genetic variance; σ 2 pe = animal per manent environmental variance σ 2 p = phenotypic variance. To determine the most appropriate model likelihood ratio tests (LRT) was used. The significance of model comparison was done from univariate analysis of animal models with and without including the effects as a random effect and compared the final loglikelihoods (Maximum log L) by chisquare distribution for α = 0.05 with one degree of freedom [20]. An effect was considered to have a significant influence when its inclusion caused a significant increase in log likelihood, compared with the model in which it was ignored.
The LRT was distributed as a χ 2 statistic with degrees of freedom equal to (p f -p r ). Where LRT = Log likelihood ratio test, L(x) f = maximum likelihood for full model, L(x) r = maxi mum likelihood for reduced model, P f = number of parameter for full model, and P r = number of parameter for reduced model. If the chisquare distribution value is significance at (p<0.05) the full model is best fit the data [20].

Fixed effects
Age at first lambing: The overall least squares mean±standard error (SE) and coefficient of variation of AFL for Bonga ewes were 375±12.5 days and 19.8%, respectively (Table 2). Age at first lambing was significantly affected (p<0.05) by birth year but not by CBBP cooperatives, dam parity, birth season and birth type ( Table 2). In the previous study by Edea et al [21] AFL of Bonga was 447±93 days but in the current study it was shorter. Age at first lambing was decreasing across years from 423±24.4 to 361±14.4 days. This indicates that the breeding program had shortened AFL due to selection of fastgrowing rams and using them for breeding. According to Ayele and Urge [22] year of birth of lamb influenced AFL through its effect on feed supply and quality.
Lambing interval: The overall least squares mean±SE of LI of Bonga sheep was 284±9.9 days. It was influenced by lamb ing year, CBBP cooperative, dam parity and season of lambing (p<0.05), however, the influence of lambing type was non significant ( Table 2). The LI decreased from 302±11.6 days in 2012 to 272±5.5 days in 2017 indicating that selective breed ing under CBBP yielded positive results. Kicho cooperative had a shorter LI than the others, which was 272±10.9 days and the longest was recorded from Guta cooperative, which was 295±11.4 days. The difference of LI across cooperative was mainly due to variation in the management activities like a quicker follow up of heat signs and feeding management.
The LI was shorter in wet season (278±10.0 days) compared to dry season (289±10.0 days). The LI was longer in first parity ewes (295±10.1 days) whereas it was shorter in fifth parity ewes (273±11.3 days). The possible reason for longer LI in first parity dams may be ascribed to growing the dam result ed in poor development of reproductive system. The current result was comparable with 268±63.9 days reported by Edea et al [21] for the same breed (Bonga). Both AFL and LI are highly influenced by regular supervision for heat signs be cause the animals in the study area were mainly tethered on private land and mainly controlled breeding system was prac ticed.
Litter size: The overall least squares mean±SE of Bonga sheep was 1.45±0.010 and the coefficient of variation was 36.62%. The percentages of twins and above were 40.13%. Litter size was significantly influenced by CBBP coopera tive and dam parity (p<0.001) but nonsignificant for lambing year and lambing season ( Table 2). Perusal of LS in the co operatives showed that Omashonga cooperative had highest LS (1.57±0.020) whereas the lowest (1.38±0.010) was observed in Dacha cooperative. The difference of LS across coopera tives was mainly due to variation in the management activities. In the study area multiple lambs were provided additional milk and milk products and suckling was controlled delib erately to avoid dominating either of the lambs. Ewes in the first parity had lowest (1.40±0.008) LS whereas ewes in ≥7 parity had highest (1.48±0.030) LS. The possible reason for lower LS in first parity dams may be their poorly developed reproductive systems due to their younger age.
Annual reproductive rate: The overall least squares mean± SE of ARR was 2.31±0.050 lambs/ewe/yr and the coefficient of variation was 38.84%. Annual reproductive rate was sig nificantly affected by lambing year, CBBP cooperatives, dam parity and lambing season (p<0.01) ( Table 2). The result of dam parity indicated that there was a gradual increase in ARR from 1.56±0.050 lambs/ewe/yr (First parity) to 2.67±0.080 lambs/ewe/yr (≥7 parity). This result was corresponding to the gradual increase in the LS in succeeding parities ( Table  2).

Model comparison
Importance of including both or one of each of additive ma ternal genetic or animal permanent environment effect on direct animal genetic effect was tested using LRT to deter mine the most appropriate model fitting the dataset ( Table  3). Inclusion of both maternal additive genetic effect and an imal permanent environmental effect (model 6) to direct animal genetic effect (model 1) significantly improved the log L for LS and ARR based on LRT distributed as chisquare test statistics (Table 3). However, for LI, inclusion of animal permanent environmental effect (model 2) significantly im proved the log L but only additive genetic effect was the best fit for AFL (model 1) (Table 3).

Genetic parameter estimation (Co) variance components of reproductive traits:
The estimates of (co)variance components and resulting genetic parame ters for reproductive traits along with estimated maximum likelihood values for six models for each trait are presented in Table 3. Perusal of variance components of the best fitted model of each trait indicated that 389.22, 472.0, 0.04, and 0.13 of the total variations comprised of direct additive vari ance (σ 2 a ) for AFL, LI, LS, and ARR, respectively. These values accounted for 6.57%, 6.38%, 16.66%, and 25% of the total variance indicating limited contribution of additive effect. The ratio of animal permanent environmental variance to phenotypic variance was higher for repeatable traits (LI, LS, and ARR) which were 0.51, 0.37, and 40, respectively. This indicated that improving the animals' environment would result an improvement of these traits. This permanent effect is an environmental effect that contributes permanently to individual's phenotype and is constant across repeated mea sures and which is not transmitted to across generation. For example intrauterine environment stimuli may impact foe tal development and this permanently affects phenotype performance later in life [23].
Heritability estimation: Direct heritabilities from selected models for AFL, LI, LS, and ARR were 0.07±0.190, 0.06±0.120, 0.18±0.070, and 0.25±0.203, respectively. The heritability of LI and AFL (Table 3) was low indicating that these traits were influenced by environmental effects including ability of farmers to timely detect heat (Estrus), feeding system, and introduction of new ewes to the program and breeding sire using mechanism. Similarly, Maria et al [24] explained that low heritability estimates for these traits were expected. Direct genetic selection within the breed for AFL and LI may therefore not bring about much improvement. However, re productive traits are aggregate traits and a small improvement in these traits would mean sizeable gain in terms of overall change in the other traits and is usually realized with simul taneous change in all components. Direct heritability of LS was higher than Ethiopian Horro sheep and other exotic sheep breeds [25,26]. Heritability estimate in the current study for AFL is lower than those reported for Brazilian Table 3. Estimates of (co)variance components and genetic parameter for reproductive traits from univariate and repeatability animal model analyses , σ 2 c , σ 2 pe , σ 2 p : variance of direct, maternal, residual, maternal permanent environment, animal permanent environment and phenotypic, respectively; σ am , covariance between direct and maternal; h 2 a , h 2 m , h 2 t heritability of direct, maternal and total, respectively; c 2 , ratio of maternal permanent environmental variance to phenotypic variance; pe 2 , ratio of animal permanent environment to phenotypic variance; r am , genetic correlation between direct and maternal; r, repeatability; log L, maximum log likelihood; LRT, likelihood ratio test X 2 chi-square test value; hyphen (-) indicate for equal number of parameter models. Santa Ines sheep 0.13±0. 10 [27] but comparable with 0.04± 0.017 reported for LI [27]. According to Maria et al [24] esti mated 0.04 and 0.06 heritability for AFL and LI from multi breed meat sheep population in Brazil. Average heritability for AFL and LI was low 0.07 and 0.02, respectively for Iranian LoriBakhtiari sheep [28].
Correlation estimate: Bivariate analysis revealed negative genetic correlation (-0.66±0.03) between AFL and LS but positive for phenotypic correlation (0.03±0.052). This indi cated that a fast growing ewe lamb has the ability to achieve a higher prolificacy rate. The higher the ovulation rate, the more oocytes will be available for fertilization during the es trous and therefore increases the possibility of a larger litter [29]. Correlations between LI and LS were 0.844±0.850 (ge netic) and 0.023±0.016 (phenotypic correlation). The high genetic correlation between LI and LS indicated that similar genetic factors influence these traits. This association showed that continuous reduction in LI result a reduction of LS or twining ability. Correlated response is expected between AFL and LS due to medium and negative correlation but not ex pected for LI and LS through improvement of either of the traits due their strong and positive correlation. Also, there was strong (0.959±0.116) genetic and 0.753±0.007 (pheno typic) correlation between LS and ARR (Table 4).
Repeatability (r) estimate: The estimate of repeatability of LI, LS, and ARR was 0.57, 0.54, and 0.65, respectively, for Bonga sheep. High repeatability indicates that culling of animals based on performance in single or only few ini tially available records could be done. The result was much higher than Horro sheep for LS 0.12 [25] and Pelibuey ewes of southeastern México 0.06±0.20 and 0.12±0.04 for LI and LS, respectively [30]. Higher repeatability than heritability estimates of a trait showed that the traits were influenced by nonadditive genetic effects and permanent environmen tal effects, and to improve these traits one should improve environmental effects or management of flock in first step [31].
Genetic trends: Genetic trends of different reproductive traits in Bonga sheep ( Figure 1) were estimated using regres sion of the average predicted breeding values obtained from the best fitted model for each trait on year of birth for AFL and year of lambing for LI and LS traits. The annual genetic decreasing value for each reproductive trait were -0.0281 days, -0.016 days and -0.0002 lambs for AFL, LI, and LS, respectively whereas positive 0.0003 lambs/ewe/yr for ARR. Statistically all reproductive traits annual genetic trend was not significance. The breeder aims to reduce both AFL and LI and negatively directed EBV for these traits was in the right direction. However, for LS the breeder aims to increase twining rate but negative EBV caused decrease in LS. The possible reason for this was more emphasis for improving body weight due to selection of higher body weight sires born in single for breeding service, negatively correlated between body weights with LS and reduction in LI due to them being positively correlated. Therefore, there needs to be special attention for either onfarm or onstation conser vation for prolific flocks. Positive and better annual genetic trend for LS was estimated from Horro sheep 0.0009±0.004 [32]. Related annual trend of AFL -0.012 days/yr and LS -0.0003 lambs/yr was recorded from Brazilian Santa Ines sheep [27].

IMPLICATIONS
Difference in heritability estimates obtained from the differ ent models suggests that model choice is an important aspect for obtaining reliable parameter estimates to be used in pre diction of breeding values. Heritability estimates for age at first lambing and lambing interval were low indicating difficulty of improving these traits through direct selection because these traits are highly influenced by environment. Therefore, improvement of such traits could be made through manipu lation of production management to reduce environmental influence. The possible reason for negative trend of litter size may be a greater emphasis for improving body weight and litter size being negatively correlated with body weight. There fore, this trait needs special attention for either onfarm or onstation conservation for prolific flocks.

CONFLICT OF INTEREST
We certify that there is no conflict of interest with any financial organization regarding the material discussed in the manu script.