Genetic analysis of sensitivity of Bonga and Menz sheep to climate variation under tropical conditions

Climate change is an undeniable global reality, bringing about unprecedented shifts in weather patterns, temperature extremes and environmental conditions. Animals living in high ambient temperatures tend to store excess radiant energy they do not use. When high temperature is combined with high relative humidity, as is frequently observed in many regions of Sub-Saharan Africa, evaporative heat dissipation is hampered, resulting in heat stress in the animal (Bohmanova et al., 2007). According to UNDP (2011), climate projection models suggest that temperature in Ethiopia will warm by 1.4–2.9°C in all seasons by the 2050s. Consequently, animal heat stress in Sub-Saharan Africa is predicted to rise due to climate change-related global warming (Ekine-Dzivenu et al., 2022). This temperature change may reduce the intake of animals, adversely affects productive and reproductive performance, reduces the quality and availability of feed and induces frequent disease outbreaks (Baumgard and Rhoads, 2013; Gowane et al., 2017; Contrerar-Jodar et al., 2018). Small ruminants are critical for food security and livelihood, especially under extremely stressful and diverse climatic conditions. Therefore, ensuring adaptive capacity and increasing the productivity of small ruminants under such variable climatic conditions is a major challenge, and efforts are needed to understand the response of animals to environmental stressors and improve the adaptability and productivity of sheep.

The interaction of the genetic makeup of animals and the environment where they are managed affects their performance, health and adaptability. Therefore, understanding the performance response of animals to the environmental variation is quite important to design and implement an effective genetic improvement programme, which aims to improve the overall efficiency and sustainability of sheep production. In the era of climate change, the selection of animals to be minimally affected by climate variability or to rapidly return to the state pertained to before exposure to a disturbance could be an alternative means of sustaining small ruminant production (Berghof et al., 2019; Joy et al, 2020; Waters et al., 2022). Reaction norm models are one option wherein each individual’s response to climate variation is characterized by a distinct random regression curve whose trajectory is defined by a continuous environmental descriptor (Falconer, 1990; Martin et al., 2011), such as temperature-humidity index (THI). This model has been used to simulate the performance of the animals and the underlying breeding value as a function of a THI to study heat tolerance in small ruminants (Menéndez-Buxadera et al., 2014; Sánchez-Molano et al., 2019; Sánchez-Molano et al., 2020; Tsartsianidou et al., 2021), dairy cattle (Santana et al., 2015; Zhang et al., 2019; Yin et al., 2023) and beef cattle (Toral et al., 2004; Santana et al., 2013). However, the potential to survive and produce under variable climate conditions is different for different breeds. Although Bonga and Menz sheep have been the subject of numerous studies, including resilience to climate variation (Tesema et al., 2025), little is known about their performance sensitivity to climate variability and genetic architectures at a specific point of environmental descriptor.

A comprehensive investigation into the multifaceted effects of heat stress on sheep, considering both short-term and long-term impacts, is essential for informing adaptive management practices and designing breeding strategies that enhance the robustness and productivity of sheep in a changing climate. Therefore, this study aimed to evaluate the sensitivity of the productive and reproductive performance of Bonga and Menz sheep to climate variation and to estimate genetic trends related to general production capacity and specific ability to respond to climate variations.

A community-based breeding program (CBBP) was introduced in Ethiopia in 2009 in four sheep breeds (Menz, Bonga, Horro and Afar) by ICARDA, ILRI, BOKU University, the Ethiopian National Agricultural Research System and collaborators. This breeding approach aimed to increase the productivity of indigenous sheep breeds and ensure their sustainability under in-situ conditions. Accordingly, Molale, Dargegn and Sinamba villages were established in 2009 for Menz sheep. Likewise, Boqa and Shuta villages were established in the same year for Bonga sheep. Efforts have been made to enhance the productivity of small ruminants through within-breed selection, improving producers’ income and meat consumption significantly (Haile et al., 2019; Haile et al., 2020). Currently, the CBBP is the preferred technique for the genetic improvement of small ruminants, attracting global interest and engaging various stakeholders in upscaling this approach (Haile et al., 2019). The pedigree and performance data were collected as a routine practice from 2009 to 2022 and stored in DTREO databases.

Menz and Bonga sheep’s performance and pedigree data pertaining to the period 2009–2022 were obtained from the DTREO database of ICARDA. Data considered were growth (live weight at birth, three-month, six-month and yearling age) and reproductive data (number of lambs born, lambs weaned, total weight of lambs at birth, weight of lambs at weaning, annual reproductive rate (ARR) and lambing interval). Weight gain (ADG) and Kleiber ratio (KR) were computed at four growth phases: birth to three months, three months to six months, six months to yearling age and birth to six months of age. ADG (g/day) is computed as ADG = [(final weight – initial weight)/number of days between two weighing dates] × 1000. The KR is calculated as the mean daily weight gain divided by the weight at a specific growth phase raised to the power of 0.75 (Kleiber, 1947).

The ARR was computed as follows: (ARR) = (number of lambs born per ewe × 365)/lambing interval. The lambing interval is the difference between two consecutive dates of lambing. Each animal had up to four records at different growth phases for each growth trait, whereas each ewe had up to eight records for reproductive traits. Preliminary tests, such as checking the data for outliers and the pedigree of animals for identification number duplication and bisexuality, were conducted. Performance records outside the range given by the mean of the contemporary group ± 3 standard deviations were excluded and the final dataset is shown in Table 1.

The study areas lacked meteorological stations and the stations in the vicinity of the study areas were unable to provide the required daily meteorological data. Therefore, using global positioning system coordinates of the area with a spatial resolution of 0.5°, daily meteorological data (temperature, relative humidity and rainfall) for each site were extracted from the NASA POWER database. The THI is a bioclimatic index widely used to measure the degree of heat stress on dairy animals (Bohlouli et al., 2013). Besides, THI is a convenient and noninvasive way to assess heat stress in animals (Ouellet et al., 2019). Thus, THI was computed according to the NRC (1971): THI = (1.8 × T + 32) - [(0.55–0.0055 × RH) × (1.8 × T - 26)], where T is the average daily temperature (°C) and RH is the average daily relative humidity. The average rainfall and THI are shown in Figures 1 and 2.

To the best of our knowledge, the studies on the time frame for the cumulative effect of climatic variables on sheep performance are limited. Sánchez-Molano et al. (2019) used the average climatic variables during 10 days before trait measurement day and Tesema et al. (2025) used 30 days before the measurement date. Therefore, the average THI and RF for the 30 days preceding the performance measurement date were combined with the performance record of each animal to evaluate their sensitivity to the cumulative effect of climate variation. The climate (THI and RF) and animal performance data were combined and analyzed separately as follows: for reproductive traits, such as the total weight of lambs at birth per ewe and the ARR, the performance matched the mean of weather variables 30 days preceding the lambing date. For the total weight of lambs at weaning and the number of lambs weaned per ewe, the performance was matched with the average weather variables 30 days prior to the weaning date. However, the performance data for the number of lambs born per ewe and lambing interval were matched with the average climatic variables of 30 days preceding the mating date and 30 days after lambing, respectively. The mating date was not recorded in the recording book for both sheep breeds. According to Tozlu Celik et al. (2021), the average gestation length of sheep is about 150 days. Thus, the 150 days were subtracted from the lambing date to obtain the mating date.

Following a data quality assessment, data analysis was conducted using WOMBAT software (Meyer, 2007) to fit a univariate reaction norm animal model. Legendre polynomials with second order was used to model the additive genetic effect along the THI and rainfall for all the traits. The THI and RF were the climatic variables used in the reaction norm model. The residual variance was assumed to be homogenous due to sample size.

The general model for random regression model for lamb traits:

where Y_ijkmnopqst is a set of records of dependent variables; B_i is the fixed effect of birth type i; S_j is the fixed effect of sex j; X_k is the fixed effect of season k; V_l is the fixed effect of year l; P_m is the fixed effect of parity m; G_n is the fixed effect of age n (as a covariate); W_o is the fixed effect of dam post-partum weight o; b_p is a covariable modeled by a Legendre polynomial (Φ_r) of second order (r = 2), was included to represent the average response of the population along the THI and RF scale for each of the dependent variables studied; a_q is a vector of additive genetic random effects due to the q^th animal that produces the records and their antecessors without data; pe_s is an individual random permanent environmental effect due to repetitions of measurement per animal; Φ_r is a vector of second order Legendre polynomial covariables evaluated at the corresponding THI and RF value and e_ijklmopqst is a random error term common to all observations.

The general model for random regression model for reproductive traits:

where Y_ijklmnop is a set of records of dependent variables; X_i is the fixed effect of season i; V_j is the fixed effect of year j; P_k is the fixed effect of parity k; W_l is the fixed effect of dam postpartum weight l (covariate); b_m is a covariable modeled by a Legendre polynomial (Φ_r) of second order (r = 2), was included to represent the average response of the population along the THI and RF scale for each of the dependent variables studied; a_n is a vector of genetic random effects due to the n^th animal that produces the records and their antecessors without data; pe_o is an individual random permanent environmental effect due to repetitions of measurement per animal; Φ_r is a vector of second order Legendre polynomial covariables evaluated at the corresponding THI and RF value and e_ijklmnop is a random error term common to all observations. The variances – covariance function of investigated traits for any pair of THI and RF points as defined by Meyer and Hill (1997) was: G= ∅K∅′, where G was a covariance matrix of the direct additive genetic effect of order t for breeding values at t given environmental descriptor (THI and RF), ∅ a t by k matrix with orthogonal polynomial and ∅′ the transpose matrix of ∅. K was the matrix with variances and covariance of the random regression coefficients.

The correlation between the estimated breeding values of animals over THI and RF was estimated to evaluate the influence of heat stress and rainfall on animals. Correlations between the same productive and reproductive traits across different levels of environmental descriptors were plotted using R software. The average estimate of the intercept and slope of the reaction norm was regressed on year to estimate the genetic trend and the significance of the trend was evaluated using a simple linear regression model. An individual’s lamb and ewe sensitivity to climate variability was classified according to Mattar et al. (2011) based on the standard deviation of slopes of reaction norm (σ_b) as an absolute value of b_j: |b_j| < σ_b indicates robust genotypes, σ_b ≤ |b_j| < 2σ_b indicates plastic genotypes and |b_j| ≥ 2σ_b indicates extremely plastic genotypes. The ratio of variances of the slope and intercept was also calculated to infer genetic variability associated with the environmental sensitivity of the reaction norm.

The heritability estimates across different THI values are illustrated in Figures 3 and 4. A moderate heritability was observed at lower and higher THI values for the Bonga sheep breed’s WT, ADG and KR. Menz sheep’s heritability estimates for WT and KR were low to moderate across THI. Likewise, TLWW in Menz sheep had a moderate level of direct heritability over THI. The moderate level of heritability indicates a possible genetic improvement for heat tolerance by selecting for the direct genetic component of these traits under heat-stressed conditions. The observed level of genetic variation indicates the possibility of selecting animals that combine high overall performance with a stable genetic ranking across environments.

An increasing heritability pattern over THI was observed for TLBW, while a decreasing heritability pattern was observed for LSB, LI and ARR in Bonga sheep. Nevertheless, the heritability estimates for TLWW and LSW of Bonga sheep were relatively high at low and high THI values (Figure 3). For Menz sheep, although there is a difference in magnitude, LSB, LI, ARR and LSW had a lower heritability estimate under moderate THI values (57–63) than the estimates at high and low THI (Figure 4). A slight decrement and increment of heritability estimates over THI were observed for TLWW and TLBW, respectively. On the other hand, Menéndez-Buxadera et al. (2014) noted that heritability reduced as THI increases for Merino Grazalema sheep. In line with the current result, a U-shaped heritability estimate over THI was noted by Mancin et al. (2024) for the productive and fertility traits of Reggiana cattle and by Cheruiyot et al. (2020) for Australian Holstein cattle. This shape suggests that the contribution of additive genetic variation to phenotypic variation was strong at environmental extremes.

A wide range of heritability was observed for WT (0.04–0.33), ADG (0.12–0.43) and KR (0.04–0.46) of Bonga sheep. Likewise, it was wide for WT (0.08–0.23), ADG (0.005–0.13) and KR (0.0001–0.21) of Menz sheep. Regarding reproductive traits, a variability of heritability estimates was observed for the TLBW (0.004–0.119) and LSB (0.055–0.188) of Menz sheep. Similarly, a variability of heritability over THI was observed for TLWW (0.027–0.208) and ARR (0.01–0.204) of Bonga sheep. Studies done regarding this issue on sheep are limited and, thus, we compared the magnitude of estimates with other livestock species. In line with the current result, a wide range of heritability (0.10–0.35) was noted by Carabano et al. (2014) and 0.14–0.31 by Bohlouli et al. (2013) for production traits of cattle. The observed variation in direct heritability of some productive and reproductive traits across different THI values suggests that the animals within the population respond differently to selection depending on the environmental conditions in which they are evaluated for these traits. Likewise, Toral et al. (2004) state that sets of various genes influence a trait in different contexts and that the expression of these genes varies with the degree of similarity or difference between environments. Besides, Cammack et al. (2006) noted that variations in humidity and temperature can turn particular genes on and off. According to the same author, in male mice under heat stress, variations in the phenotypic expressions of reproductive traits caused by varying heat stress levels also led to variations in genetic parameter estimates. The heritability estimates of reproductive traits for both breeds were less variable over THI values compared with the heritability of production traits, corroborating with the finding of Mancin et al. (2024) for Reggiana cattle. This indicates that the genetic control of reproductive traits was relatively stable across different THI values.

The heritability estimates across different rainfall (RF) amounts for Bonga and Menz sheep are illustrated in Figures 5 and 6, respectively. A moderate heritability was observed at lower and higher rainfall amounts for the Bonga sheep breed’s WT, ADG and KR. This moderate heritability estimate suggests a possible genetic improvement for adaptation to low rainfall by selecting for the direct genetic component of these traits under low rainfall conditions. Nevertheless, the heritability estimates for reproductive traits of Bonga sheep were found to be within the low range. Despite the magnitude, the heritability estimates for LSB, ARR and LI tended to increase with the rainfall amount. The heritability estimate for TLWW was consistent across different amounts of rainfall. The heritability estimates for all traits of Menz sheep except for TLWW were low in magnitude, although some traits such as WT, KR, ARR and LSB increased with the rainfall amount. Heritability estimates for TLWW in Menz sheep decrease with increasing rainfall amounts. Consistently, a meta-analysis done by Safari et al. (2005) indicates that reproductive traits of sheep are lowly heritable. The heritability estimates for productive and reproductive traits of Bonga and Menz sheep varied across the rainfall amounts. The observed variation in direct heritability of productive and reproductive traits across different rainfall amounts suggests that the animals within the population respond differently to selection depending on the rainfall amount in which they are evaluated for these traits. Consistently, Cheruiyot et al. (2020) noted that the heterogeneity of genetic variance and heritability across environmental descriptors could lead to the re-ranking of sires in Australian Holstein cattle. The magnitude of the heritability estimate suggests that the contribution of additive genetic variation to phenotypic variation of reproductive traits at different rainfall amounts is small and the major governor of these traits is the management of the animal. Therefore, improving the management system would significantly enhance reproductive performance under a variable climate conditions. On the other hand, the moderate heritability estimate for TLWW implies that there is sufficient genetic variation to improve TLWW through within-breed selection.

The correlation of the estimated breeding value of traits between different levels of THI for the investigated traits of Bonga and Menz sheep is shown in Figures 7 and 8, respectively. The genetic correlation for TLBW, ARR, LSB and LI of Bonga sheep between THI values was greater than 0.85. Likewise, the genetic correlation of TLWW in Menz sheep between THI values was high (0.99). This high genetic correlation suggests that the genes causing genetic variance in these traits between environments behave uniformly. Therefore, the best animals in low- to medium environments (THI) also perform best regarding these traits in high environments (THI). However, the genetic correlation among different levels of THI for the other traits of Bonga sheep (WT, ADG, KR, TLWW and LSW) and Menz sheep (WT, ADG, KR, TLBW, LSB, LSW, ARR and LI), particularly at high and low THI, was found to be low and antagonistic for some traits, which indicates some level of performance variability across environmental gradients. The genetic correlation estimates were higher as the THI became more similar and lower for more distant THI values. A similar pattern of genetic correlation was reported by Menéndez-Buxadera et al. (2014) for Merino Grazalema sheep. Strandberg et al. (2000), Strandberg et al. (2009) and Chiaia et al. (2015) consistently revealed strong negative to positive associations for the reproductive and productive traits of Holstein, Nordic and Nellore cattle. These low and antagonistic genetic correlations suggest that the genes causing genetic variance in these traits between environments do not behave uniformly. Thus, selecting for improvement in the unfavorable environment (THI) might lead to a decline in performance in a favorable environment (re-ranking of animals), particularly for traits that had a significant and antagonistic correlation between environments. Although the correlation showed some level of re-ranking of animals over THI for some traits, selection breeding animals for each level of THI may not be feasible and practical within one village, however, breeding rams could be selected for each village.

The genetic correlation of the estimated breeding value of traits between different levels of rainfall for investigated traits of Bonga and Menz sheep is illustrated in Figures 9 and 10, respectively. The genetic correlations for ARR, LI and TLWW of Bonga sheep between different amounts of rainfall were high (rG ≥ 0.80). Likewise, the genetic correlation among points of rainfall for LSB of Bonga sheep was ≥ 0.60. However, the correlation of other traits of Bonga sheep (WT, ADG, KR, LSW and TLBW) between various rainfall amounts varied from 0 to 1 and the magnitude of correlation between low and high rainfall amounts was found to be lower or near zero. The genetic correlation estimates were higher as the rainfall amount became more similar and lower for more distant rainfall values. A similar pattern of genetic correlation was reported by Menéndez-Buxadera et al. (2014) for Merino Grazalema sheep and Chiaia et al. (2015) for cattle. The TLBW and TLWW of Menz sheep between different levels of rainfall were highly correlated (rG ≥ 0.85). The correlation of breeding value of WT between different rainfall levels was moderate to high (rG = 0.4–1.0). The genetic correlation for ADG, KR, ARR and LI in Menz sheep between different rainfall amounts ranged between 0 and 1. Nevertheless, the LSB and LSW of Menz sheep in low rainfall were negatively correlated with high rainfall (rG = −0.4–1.0), although the correlation was positive and high for the nearest points. The moderate-to-high and positive genetic correlation between rainfall amounts for some traits suggests that the best animals in low-to-medium rainfall also perform best regarding these traits in high rainfall amount. On the other hand, negative genetic correlation estimates for LSB and LSW of Menz sheep indicate that genes causing genetic variance in these traits between environments (between high and low rainfall) behave differently across environments.

The intercept of the regression is the relative average performance at the standard level of the environmental variable. The slope can be interpreted as the sensitivity of the performance of animals to changes in the environmental variable (Ribeiro et al., 2017). A high intercept (high overall performance) with a flat slope (robustness of performance) is preferable (Kolmodin et al., 2002). Understanding the genetic correlation between slope and intercept is important for predicting how a population may respond to changing climatic conditions. Besides, it can inform selective breeding strategies, helping to improve specific aspects of trait expression, whether in terms of performance in a mean environment or environmental responsiveness. The genetic correlation of the intercept and slope of different productive and reproductive traits of Bonga and Menz sheep using THI as an environmental descriptor is illustrated in Table 2. A low and positive genetic correlation between the intercept and slope was obtained for ADG, KR and TLWW for Bonga sheep, implying that animals with higher average breeding values in the standard environment may not respond better to improvement in environmental conditions (THI). The correlation between the intercept and slope of the reaction norm for TLBW and WT was high to moderate for the same breed, indicating that animals with higher breeding values in a standard environment responded better to improvement in environmental conditions or it indicates that individuals with a higher expression of WT and TLBW in a standard environment also tend to exhibit a more pronounced response (steeper slope) to environmental changes (Rauw and Gomez-Raya, 2015). Similarly, a higher correlation estimate between slope and intercept was noted by Carvalheiro et al. (2019) for body weight of beef cattle. Animals with high productivity levels may need to reallocate resources from other physiological processes to fulfill their requirements. This resource reallocation might render them less capable of coping with stress, ultimately increasing their vulnerability to environmental changes (Rauw, 2009). The correlation for LSW of Bonga sheep seems nonsignificant, as the standard error estimate indicates.

However, the correlations between the intercept and slope for LSB, ARR and LI of Bonga sheep and TLWW of Menz sheep are strongly negative, corroborating the result of Zhang et al. (2019) for fertility traits in Danish Holstein cattle. This negative genetic correlation estimate between slope and intercept implies that an animal with a high performance at the mean environmental descriptor may have a low genetic sensitivity to environmental variability (flatter slope) or a more stable trait expression across different environments. Besides, this negative relationship may suggest trade-offs, where animals excelling under standard conditions may not be as adaptable to environmental variations.

The very low and negative genetic correlation of intercept and slope for body weight in Menz sheep implies a very low association of performance at the mean environment and environmental responsiveness. This indicates that there is an opportunity to select for increased body weight and reduced sensitivity simultaneously. A moderate to high slope and intercept genetic correlation was observed for Menz sheep’s ADG, KR and TLBW, which suggests that a higher general productive performance is genetically associated with increased ADG, KR and TLBW in better production environments. Thus, individuals with a higher expression of these traits under a standard environment also tend to exhibit a more pronounced response (steeper slope) to environmental changes or a higher degree of plasticity of the trait to environmental changes. Regarding environmental adaptation, the moderate to high positive correlations may indicate that animals with high performance in standard environments are more adaptable to environmental changes. Nonetheless, the correlation for LSB, LSW, ARR and LI of Menz sheep was nonsignificant, as the standard error estimate indicates. The trajectories of various animals’ performances along the environmental descriptor are not parallel, so even with a very high correlation between the intercept and slope, re-ranking of the animals is still possible as long as the climatic variable range is sufficiently large and the slope variance is significantly different from zero (Zhang et al., 2019). The breeding goal and the stability of environmental conditions could determine the selection of steeper or flatter slopes. A steeper slope may be beneficial in changeable situations, whereas a flatter slope may be sufficient in stable ones. Besides, a steeper slope could be preferable for some traits, such as body weight or weight gain, to maximize production under various environmental conditions.

The genetic correlation between the intercept and slope of the reaction norm for productive and reproductive traits of Bonga and Menz sheep using rainfall as an environmental descriptor is shown in Table 2. A moderate to high (positive) genetic correlation of slope and intercept was observed for Bonga sheep’s body weight, whereas it was negatively high for TLWW. Similarly, a higher correlation estimate was noted by Carvalheiro et al. (2019) for beef cattle. The observed high genetic correlation of slope and intercept for body weight suggests that animals with better estimated breeding value for the average environment also tended to show higher sensitivity to environmental variation. The genetic correlation between slope and intercept for ADG, KR and TLBW of Bonga sheep was found to be low and positive, corroborating the report of de Paula Freitas et al. (2021) for cattle and suggesting that selection for these traits would result in low environmental sensitivity of animals. The body weight, ARR and LSB of Menz sheep had a high and positive genetic correlation between intercept and slope. Likewise, the genetic correlation for KR of Menz sheep was moderate and positive. These high genetic correlation estimates between the intercept and slope suggest that simultaneous selection for increased performance and reduced sensitivity may be limited. A similar observation has been made for the body weight of beef cattle (Carvalheiro et al., 2019). However, the interpretation of genetic correlation of other traits should be done with caution, as the standard error of the estimate failed to estimate. Moderate to high genetic correlation of slope and intercept for some traits in the current study suggests that individuals with a higher expression of these traits under a standard environment also tend to exhibit a more pronounced response (steeper slope) to environmental changes or a higher degree of plasticity of the traits to environmental changes. Likewise, Ribeiro et al. (2017) noted a high correlation for the body weight of Nellore cattle. Regarding adaptation, the moderate to high positive correlations may indicate that animals with high performance in standard environments are more adaptable to environmental changes. The negative genetic correlation between intercept and slope for TLWW in Bonga sheep suggests that selection for increased performance would result in increased susceptibility to rainfall variability.

The climate sensitivity of Bonga and Menz sheep is shown in Table 2. The ratio of variances of the slope and intercept of reaction norm (βg) for Bonga sheep’s WT, ADG, KR and TLWW using THI as an environmental descriptor was greater than one. On the other hand, other traits (LSB, TLBW, LSW, ARR and LI) for this sheep breed had a lower βg value. Likewise, the body weight and ARR of Bonga sheep had high slope and intercept variance ratios, using rainfall as an environmental descriptor. The variance ratio for other traits (ADG, LSB, LSW, TLBW, TLWW and LI) of Bonga sheep using rainfall as an environmental descriptor was found to be less than one, which is consistent with the report by de Paula Freitas et al. (2021) for the body weight of cattle. This low ratio indicates that Bonga sheep showed the lowest genetic variability associated with environmental sensitivity proportional to the intercept in terms of these traits. The high value of βg indicates substantial genetic variability in the responsiveness of the trait to environmental changes. However, a lower βg value indicates less genetic variability in how the trait responds to climate changes or individuals within the population exhibit more consistent responses to climate variation. The low to moderate βg value for reproductive traits of Bonga sheep indicates that Bonga sheep had less genetic variability in responsiveness to climate variability.

Menz sheep had a higher βg for ADG and KR and lower βg values for other traits, except LSB, using THI as an environmental descriptor. Using rainfall as an environmental descriptor, the βg for ADG and LSW was high, whereas it was low for other productive and reproductive traits of Menz sheep. The higher values of βg for ADG, KR and LSW in Menz sheep suggest a higher proportion of genetic variability in responsiveness compared to trait expression under standard environmental conditions. The value of βg for LSB indicates a balanced genetic influence on both the trait’s expression under the mean environment (THI) and its responsiveness to environmental change. However, the lower βg value for most of the traits of Menz sheep implies that genetic factors have a relatively greater impact on the trait’s expression under the standard environment than its responsiveness to the environment (THI and rainfall). This indicates that individuals within the population of Menz sheep exhibited more consistent responses to environmental changes (THI) for most of the investigated traits, except for ADG, KR and LSW.

The ability of an organism to develop different phenotypes in response to changes in its environment is known as phenotypic plasticity. Phenotypic plasticity allows individuals to adjust their physiological, morphological or behavioral traits to suit the prevailing environmental conditions better (Sommer, 2020), whereas phenotypes exhibiting small variation between environments are called ‘robust’ (de Jong and Bijma, 2002). There was a high incidence of robust animals (50.4%–81.1% of lambs and ewes) for all productive and functional traits of Bonga sheep using THI as an environmental descriptor (Figure 11). Likewise, the proportion of robust lambs and ewes for Menz sheep varied from 50.1% for TLWW to 95.4% for LSB. Consistently, the incidence of robust animals for both sheep breeds was high (52.9%–96.1%), using rainfall as an environmental descriptor (Figure 12). This implies that Bonga and Menz sheep can withstand climate change without exhibiting significant deviations from their normal productive and reproductive performance or a certain level of inherent stability and adaptability to adverse conditions. Likewise, Santana et al. (2013) noted that most animals had intermediate sensitivity. According to Strandberg et al (2000), an optimal breeding system would feature high-performing genotypes with a slope approaching zero, resulting in more robust animals that excel across diverse climatic conditions. The numbers of plastic and extremely plastic animals were relatively high for TLBW, LSB, ARR and LI of Bonga sheep compared to other traits, using THI as a climatic variable. The proportions of plastic and extremely plastic animals were also relatively greater (49.9%) for TLWW of Menz sheep compared to other traits, using THI as a climatic variable. Regarding adaptation, animals with higher plasticity may excel in adjusting to significant environmental changes, while robust animals might struggle to survive under specific conditions (Santana et al., 2013).

The genetic trend for intercept (level of production) and slope (environmental sensitivity) of productive and reproductive traits of Bonga sheep using THI as a climatic variable is presented in Figure 13. The regression coefficient of the intercept of the reaction norm for WT, ADG, KR, TLBW, TLWW, LSB, LSW, ARR and LI were 0.004 kg, −0.053 g/day, −0.0064 g/day/metabolic weight, −0.0019 kg, 0.00005 kg, −0.0005 lambs, −0.00001 lambs, −0.0006 lambs/ewe/year and 0.012 day/year, respectively. These values indicate that Bonga sheep tended to show higher performance for body weight and TLWW. The genetic trend of intercept was significant (p < 0.05) only for TLBW and LSB and the performance of these traits slightly decreased across the year.

The regression coefficients for slope were −0.0025 kg, −0.426 g/day, −0.0426 g/day/metabolic weight, −0.0002 kg, −0.0018 kg, 0.0002 lamb, 0.00008 lamb, 0.0002 lambs/ewe/year and −0.0047 day/year for WT, ADG, KR, TLBW, TLWW, LSB, LSW, ARR and LI, respectively (Figure 13). This result depicts that Bonga sheep moved toward sensitivity slightly in terms of LSB, LSW and ARR, indicating that these animals are becoming more responsive to environmental changes. However, Bonga sheep move toward less responsiveness to environmental change in terms of WT, ADG, KR, TLBW, TLWW and LI. The genetic trend of the slope of reaction norm for ADG, KR, TLBW and LSB of Bonga sheep was significant (p < 0.05). According to de Jong and Bijma (2002), the plasticity of production and health traits should be low, while high plasticity of resource (feed) intake is good. The higher and significant correlations were observed among feed efficiency and KR of small ruminants (0.96) and sheep (0.95), as noted by several authors (Tesema and Tiruneh, 2021 and Talebi, 2012), respectively. Thus, improving of this KR could enhance feed efficiency and, thereby, the coping ability of sheep for feed shortage resulting from climate change.

The genetic trends for intercept and slope of productive and reproductive traits of Menz sheep are presented in Figure 14. The regression coefficients of the intercept for WT, ADG, KR, TLBW, TLWW, LSB, LSW, ARR and LI were 0.0015 kg, −0.032 g/day, −0.0032 g/day/metabolic weight, −0.0001 kg, 0.0034 kg, −0.00005 lambs, −0.00002 lambs, 0.0009 lambs/ewe/year and −0.2315 day/year, respectively. However, the genetic trend of intercept was significant (p < 0.05) for LSB, ARR and LI. This indicates that the performance of LSB and LI tended to show a decreasing trend, while ARR tended to show high performance. Indeed, the improvement of litter size is determined by the area’s resources, feed availability and environmental conditions. The increment of WT and TLWW could be explained by more emphasis given to body weight in the existing selective breeding program.

The regression coefficients of slope for WT, ADG, KR, TLBW, TLWW, LSB, LSW, ARR and LI were −0.00002 kg, −0.0944 g/day, −0.0166 g/day/metabolic weight, −0.0001 kg, −0.0002 kg, 0.00001 lamb, 0.00004 lamb, 0.0004 lambs/ewe/year and −0.0797 day/year, respectively (Figure 14). This finding indicates that Menz sheep go to sensitivity in terms of TLWW, LSB, LSW and ARR. However, other traits move toward low sensitivity to THI variation. Nevertheless, the slope regression coefficient was insignificant (p > 0.05) for most of the traits investigated except for ARR and LI.

The genetic trend for productive and reproductive traits of Bonga and Menz sheep using rainfall as an environmental descriptor is shown in Table 3. Bonga sheep tended to show higher performance across years for WT, ADG, KR and LI, although it was not significant. However, LSB, LSW, TLBW and TLWW of Bonga sheep tended to show a significant performance reduction (p < 0.05) over the year. Bonga sheep exhibit sensitivity to performance variations due to rainfall, particularly in terms of WT, ADG, KR and LSW. However, other traits (LSB, TLBW, TLWW, ARR and LI) move toward low sensitivity to rainfall variability. Nevertheless, the trend for all traits except LSB and TLWW was nonsignificant (p > 0.05). The WT, KR, ARR, LSW and TLWW performance of Menz sheep increased genetically across years, while other traits (ADG, LSB, TLBW and LI) showed a decreasing trend across years. However, only WT was significant (p < 0.05). The ARR of Menz sheep goes to sensitivity and TLBW goes to less sensitivity to rainfall variation significantly (p < 0.05). The trend for other productive and reproductive traits of Menz sheep was found to be insignificant.

The genetic trend indicates that the performance of some traits under average environmental conditions has increased genetically over the years, but there is a slight increase in sensitivity to climate variation. Hence, selection for increased performance under optimal conditions is performed for some traits of both sheep breeds; there will probably be a loss in robustness of animals and vice versa. Similarly, Ribeiro et al. (2017) reported that as the animals are selected according to breeding values, the sensitivity to environmental changes also tends to increase due to high slope and intercept genetic correlation. The THI-based performance and sensitivity of LSB in Bonga sheep decreased and increased genetically, respectively. Likewise, the THI and rainfall-based performance and sensitivity of ARR of Menz sheep were increased genetically across years. For such situations, Usala et al. (2021) suggest that overall performance and tolerance for environmental stress should be combined in an economic selection index, with different weights depending on the likelihood of certain conditions occurring in a particular system.

The productive and reproductive performance of Menz and the reproductive performance of Bonga sheep were less sensitive to climate variation and there was a less uniform genetic response to varying climate conditions. The moderate heritability estimates for productive traits of Bonga sheep at high THI indicate a possible genetic improvement for heat tolerance by selecting the direct genetic component of these traits under heat-stressed conditions. Low and negative genetic correlations between low and high levels of THI for some traits of Bonga sheep and for most traits of Menz sheep suggest a less uniform genetic response to varying environmental conditions and some level of re-ranking of animals. However, the lower genetic variance ratio of slope and intercept for body weight and reproductive traits of Menz sheep and reproductive traits of Bonga sheep implies that genetic factors have a relatively greater impact on the trait’s performance expression under the standard environment than its responsiveness to the climate variation. Besides, high proportions of robust animals for both sheep breeds imply a certain level of inherent performance stability and adaptability to variable climatic conditions. Determining an animal’s sensitivity to environmental variation, as the findings of this study can support timely interventions to improve prductivity and prevent health and welfare issues related to climate change.

The authors are grateful to Accelerating the Impact of CGIAR Climate Research in Africa (AICCRA) project and CGIAR initiative on “Sustainable Animal Productivity for Livelihoods, Nutrition and Gender inclusion (SAPLING)” for funding this research. In addition, we would like to thank Debre Birhan Agricultural Research Center and Bonga Agricultural Research Center for data provision.

Zeleke Tesema: Writing – original draft, Software, Methodology, Investigation, Formal analysis, Data curation, Conceptualization. Tesfaye Getachew, Aynalem Haile, Berhanu Belay: Writing – review and editing, Methodology, Funding acquisition and Conceptualization. Yosef Amha, Tamrat Bekele, Teferi Demissie, Dawit Solomon, Mourad Rekik, Barbara Rischkowsky, Shanbel Besufkad, Zelalem Abate, Ebadu Areb: Writing – review and editing, Resources. All authors read, commented and approved the final manuscript.