Climate change impact on wheat and maize growth in Ethiopia: A multi-model uncertainty analysis

Ethiopia’s economy is dominated by agriculture which is mainly rain-fed and subsistence. Climate change is expected to have an adverse impact particularly on crop production. Previous studies have shown large discrepancies in the magnitude and sometimes in the direction of the impact on crop production. We assessed the impact of climate change on growth and yield of maize and wheat in Ethiopia using a multi-crop model ensemble. The multi-model ensemble (n = 48) was set up using the agroecosystem modelling framework Expert-N. The framework is modular which facilitates combining different submodels for plant growth and soil processes. The multi-model ensemble was driven by climate change projections representing the mid of the century (2021–2050) from ten contrasting climate models downscaled to finer resolution. The contributions of different sources of uncertainty in crop yield prediction were quantified. The sensitivity of crop yield to elevated CO2, increased temperature, changes in precipitations and N fertilizer were also assessed. Our results indicate that grain yields were very sensitive to changes in [CO2], temperature and N fertilizer amounts where the responses were higher for wheat than maize. The response to change in precipitation was weak, which we attribute to the high water holding capacity of the soils due to high organic carbon contents at the study sites. This may provide the sufficient buffering capacity for extended time periods with low amounts of precipitation. Under the changing climate, wheat productivity will be a major challenge with a 36 to 40% reduction in grain yield by 2050 while the impact on maize was modest. A major part of the uncertainty in the projected impact could be attributed to differences in the crop growth models. A considerable fraction of the uncertainty could also be traced back to different soil water dynamics modeling approaches in the model ensemble, which is often ignored. Uncertainties varied among the studied crop species and cultivars as well. The study highlights significant impacts of climate change on wheat yield in Ethiopia whereby differences in crop growth models causes the large part of the uncertainties.


Introduction
In large parts of the world, local economies and food security are dependent on agriculture which is very sensitive to climate variations [1,2]. Both shifts in the mean climate as well as changes in the variability can considerably affect agricultural production [3]. There is now insurmountable evidence linking anthropogenic greenhouse gas emissions to climate change [4]. If global warming continues at the current pace, the global average temperature is likely to increase by 1.5˚C by 2050 compared to pre-industrial levels [3]. This will undermine agricultural production in many parts of the world, in general, and Ethiopia, in particular. It must be seen as a very severe threat to meeting the growing world food demand [5,6]. A higher and faster temperature increase of 2˚C is projected for Africa [7]. This would lead in a major food security risk by 2050 [8]. Therefore, adapting agriculture to the changing climate is an additional and daunting problem for the continent [9]. Consequently, climate change impact assessment has developed as an important area of research in Africa [6,10], aiming at identifying adaptation options.
The impact of climate change on crop production has been studied extensively using statistical models [5,11,12], process-based crop models [13][14][15], and a suitability evaluation approach [16]. Simulations with process-based crop models driven by climate projections from dynamical downscaling have dominated the recent literature. Crop models can describe crop responses to different climate scenarios and variable agronomic management options. Globally, crop model structures exhibit quite a large variability enabling more robust projections and quantifying associated model uncertainties. This will improve the relevance of crop model outputs in different environments [17]. However, not accounting for the related uncertainties may jeopardize decision-making in the agricultural sector. There are many sources of uncertainties including those related to differences in the governing process model equations and parameterizations, climate model projections, and methods of climate downscaling. In order to help decision makers to decide between plausible decisions options, Challinor et al. [18] recommended that these sources of uncertainties should be distinguished and communicated.
Considering the complex and non-linear influence of climate on crop development and growth, quantifying, and dissecting the various sources of uncertainties in the prediction of the impact of climate change is a challenge. Recently, the multi-model ensemble approach has been introduced to study uncertainties in projected climate change impact on crop production [14,19,20]. By this, yield projections could be improved by calculating the ensemble means and quantifying the related uncertainties in comparison to the common single crop model-single climate projection simulations. One major finding from the previous studies was that differences in the crop models were the main source of uncertainty and usually dominated over the uncertainties resulting from the climate model scenarios [14,17,20].
Crop models differ not only in the equations describing crop growth, but also in the level of detail of boundary conditions needed for running them, and the equations representing related processes in the soil [21]. The role of different approaches in modelling these processes determining crop growth has so far not been analyzed in previous studies. For instance, modeling soil processes is an important component of crop modeling. However, considerable uncertainties remain to be addressed regarding the quality of predictions among different soil models of despite their significant progresses in soil models [22]. In previous multi-model studies such as that by Asseng et al. [23], the participating modeler groups contributed only with complete models, i.e. with a given combination of model equations that describe the soilplant-atmosphere system. Ceteris paribus analyses of how single model components affect the outcome of simulations were not performed, which limits assigning uncertainties to individual processes.
Studies focusing on the influence of climatic variations on crop productivity are of particular importance for regions such as Ethiopia, where only few studies have been conducted so far. Ethiopia's national economy is dominated by the agricultural sector, which accounts to 42% of the GDP, employs 85% of the population, and contributes about 90% to total export earnings [24]. Crop production is the main agricultural activity, and it is mainly practiced under rain-fed conditions. This can exacerbate the negative effects of climatic variations [25]. Since cereals are the main source of food, failures in cereal production will almost inevitably jeopardize food security. In 2018, 87% of total crop production were cereals grown on 80% of the total cultivated land [24]. Of the many cereals grown in Ethiopia, maize (Zea mays L.) and wheat (Triticum aestivum L.) are two of the three most important crops with a combined share of 42.6% of the total crop production [24].
Significant increasing trends in frequencies and magnitude of climate extremes have been documented in Ethiopia [26][27][28]. In the face of a changing climate, the likely increase of extreme events will have an adverse impact on cereal crop production in Ethiopia [10,29,30]. However, the results have shown large discrepancies in the magnitude and sometimes in the direction of the impact. For example, some studies projected maize yield to decline by 20-43% [31][32][33], other studies show that maize yield could increase by 2-51% [10,29,31]. The possible reasons of these inconsistencies might stem from limitations related to either the low number of crop models used, or the few numbers of climate model projections considered in each of the studies, or both.
Thus, the main objectives of this study were i) to quantify the impact of projected climate change on maize and wheat yields in Ethiopia, and ii) to dissect the contributions of different sources of uncertainties from the overall uncertainty in crop yield projections. To achieve the objectives, we applied a multi-model ensemble approach by running 48 crop growth and soil submodel combinations in conjunction with climate projections of ten different climate models from the Coordinated Regional Climate Downscaling Experiment (CORDEX). To the authors' best knowledge, no crop modeling study has been conducted so far using dynamically downscaled climate projections in East Africa in general and Ethiopia in particular. Models were calibrated for the phenology and growth across 3 locations representing the major wheat and maize growing areas in Ethiopia. To identify the most influential factors that affect grain yield, a sensitivity analysis was also conducted by discretely varying selected input variables and analyzing their sensitivity on modeled yield for the 30 baseline years . The potential impact of climate change on yield was computed by comparing the relative changes between the baseline and future climate simulations. Finally, comparison was made for the contribution of the different sources of uncertainty in projected impact on grain yield.

The study sites
Our study was based on experimental field data for two wheat and two maize cultivars, which were collected over several years at three research sites with differing soil and climate: Adet and Sinana for the two wheat cultivars and Kulumsa for the two maize cultivars. Adet lies at an altitude of 2170 m with a long-term (average) annual rainfall between 1000 and 1800 mm ( Fig  1). The dominant soil types in Adet are Fluvisols and Vertisols. Wheat is widely grown as the major crop following teff. Sinana is one of the major wheat producing regions in the country. It is located at the foot of the Bale Mountain west of the Great Rift Valley at an altitude of 2460 m. Its average annual rainfall varies from 600-1600 mm with two distinct rainy seasons; the main one from July-October and the minor one from March-June. The research site is characterized by a deep, fine textured and aggregated soil structure. Predominant soil types are Phaeozems and Cambisols. Kulumsa is in the central part of the country at an altitude of 2260 m. It represents a highland area with a cool climate receiving an annual total average rainfall of 700-1100 mm rainfall and a rainy season from May to October. Most of the area is a flat plain dominant by Haplic Xerosol.

Observed climate and crop data
Observed daily weather data (rainfall, solar radiation, maximum temperature, and minimum temperature) were obtained from the National Meteorological Agency of Ethiopia for the period of 1981 to 2017. Datasets from the Agriculture Modern-Era Retrospective Analysis for Research and Applications [34] were used both for gap-filling and as source to estimate unavailable data such as wind speed and relative humidity. Field data for wheat was obtained from Adet and Sinana agricultural research centers while the data for maize was obtained from field experiments conducted by the International Maize and Wheat Improvement Center (CIMMYT) in collaboration with Kulumsa agricultural research center. The experimental data we used were gathered in parts during genotype by environment experiments and extracted from the agricultural research centers' field books. The records contain summary of measured data, published in annual reports for each season and stored in the research centers' libraries. The reports include records of crop development and growth observation data as well as soil and crop management information for each field experiment. Crop data include maturity date, anthesis date and grain yield, whereas crop management data include planting date, planting density and fertilizer application dates and rates. The reports include only the average values across the replications of the experiment. The four crop cultivars considered in the study were 'Shina' (in Adet) and 'Medawolabu' (in Sinana) for wheat as well as 'Jibat' and 'Wenchi' (in Kulumsa) for maize (Table 1). Available soil data included soil texture as well as physical and chemical properties ( Table 2). The soil data were extracted from the report that match the beginning of our experimental crop data wherever available. For maize (i.e., Kulumsa site), soil data from the 2012 report were used whereas for wheat the data were extracted from the 2002 (i.e., Adet site) and 2006 (i.e., Sinana site) report.

Climate projections
Future climate change data of ten contrasting Global Climate Models (GCMs, S1 Table) were accessed from CORDEX [35] which had been downscaled to 44 km horizontal resolution in Africa (available at https://esgf-node.llnl.gov/search/esgf-llnl/). Sponsored by the World Climate Research Program, CORDEX generated high-resolution future climate projections by downscaling GCM projections using different regional climate models. A 30-year window period representing climate projections by the mid of the century (2021-2050) was considered in this study. A single Representative Concentration Pathway (RCP4.5) emission scenario was selected as the radiative forcing does not diverge in the different RCPs by the selected time [36]. The data were bias-corrected using the quantile-quantile mapping procedure [37,38], in order to match distributional characteristics of historical observations of precipitation and temperature.
Compared to the baseline climate (1981-2010), a consistent temperature increase is expected by 2050 for all study sites whereas the change in precipitation is expected to differ, both seasonally as well as between sites. Averaged over the climate models, both annual maximum and minimum temperature are projected to increase by about 1.2-1.4˚C by 2050 where the projected change in the minimum temperature would be slightly higher than that of the maximum temperature (Fig 2). Averaged across climate models, a slight increase in total annual precipitation is projected by 2050 in Sinana (+6%) and Kulumsa (+7%), whereas a decrease (>25%) is projected in Adet. However, the projected changes are highly variable among the climate models as well as between seasons ( Fig 2C). For the growing season only, an increase by 17% and 29% in Kulumsa and Sinana, respectively, and a -2% in Adet are projected. Generally, a negative relationship was evident between the projected changes in temperature and precipitation with a majority of GCMs predicted hot and dry condition ( Fig 2D).

Model description
The Expert-N modelling package. The agro-ecosystem modelling package Expert-N provides process-based models and submodules for the description of coupled processes in the soil-plant system (Fig 3). This enables the simulation of saturated and unsaturated water flux, soil heat transfer, C and N turnover in the soil, evapotranspiration, transport of mineral nitrogen, and plant processes, given abiotic boundary conditions as well as soil and crop management [39]. Its modular design facilitates multiple model configurations by combining various submodels. The crop growth submodel enables the dynamic simulation of phenological development, crop growth, canopy photosynthesis, biomass accumulation, leaf area and root distribution, root water and nitrogen uptake, and yield.
Model ensembles. There exists a multitude of crop models and modelling approaches for simulating processes in agroecosystems. In many cases, there is no clear reason for preferring one process representation over another. It also depends on the data situation whether it is possible to decide a priori whether certain process formulations are to be preferred. In fact, many submodels were developed specifically for limited data availability. In our study, we followed a multi-model approach whenever it was not possible to decide a priori which approach would be superior given our data. In the present study, a multi-model ensemble (n = 48; Fig 3) was formed by coupling four crop growth submodels (CERES, SUCROS, SPASS, and GECROS), two soil organic matter turnover submodels (SOILN and DAISY), two soil water transport submodels (based on Richards equation and the tipping-bucket (TB) approach) and three soil heat transfer submodels (LEACHN, DAISY and SHAW).
The four selected crop growth submodels took part in the world-spanning agricultural model intercomparison and improvement project (AgMIP) for wheat [15,40] and maize [13,41] under contrasting conditions and environments across the globe. The models covered a good space of approaches and methods used in many common models. Most state-of-the-art crop models originate from two crop modelling schools, the DSSAT models and the Wageningen models. In this study, we used a further developed version of the CERES model based on an earlier version of DSSAT, two models from Wageningen (SUCROS and GECROS), and the hybrid model SPASS. Table 3 presents short summaries of the main features of the four crop growth models and their differences. The models vary in the level of detail and complexity of growth process representation [42]. The CERES was designed to simulate management and environmental effects on crop growth and development [43] and tested in different environments around the world [44]. SUCROS originally developed to simulate potential crop production with unlimited supply of plant nutrients [45], was combined in Expert-N with submodules for nitrogen uptake from WAVE [46] to take into account possible N-limitation effects on photosynthesis and crop growth. SPASS was developed by adopting approaches mostly from CERES and SUCROS with some new own functions [47]. GECROS is relatively a recent model. It is based on the SUCROS model and uses a set of new algorithms that enhanced the description of genotype-environment interaction [48]. Photosynthesis is calculated based on biochemical approach. This sets GECROS apart from the rest of the models that use less advanced approaches. Whereas SPASS was originally developed in Expert-N [49,50], the other crop growth models were implemented from freely available code or implemented based on their documentations. Therefore, they may slightly differ from their original version with regard to the coupling to soil modules of Expert-N (water uptake, nutrient uptake, litter transfer to soil organic matter pools). Short summaries of the main features of the four plant models and their differences are provided by previous studies [42,51,52].
The two soil organic matter turnover submodels included here were based on the SOILN model [53] and DAISY [54,55] model. Both use a multi-pool structure of soil organic matter. SOILN distinguishes three pools with distinct turnover rates (humus, litter, and manure) whereas DAISY considers six of them. In DAISY, the soil organic matter is divided into three main pools namely: added organic matter, soil microbial biomass and native soil organic matter where each pool is further divided into two parts, fast and slow. Turnover is determined by the availability of fractions for the soil organisms [56]. Descriptions of the implementations of SOILN and DAISY in Expert-N can be found in Berkenkamp et al. (2002).
The three soil heat transfer modules included were taken from LEACHN [57], DAISY [54] and SHAW [58]. These models differ in the terms considered in the general heat transfer equation, the numerical solution of the applied equation and in the way the upper and lower boundary conditions are defined. In LEACHN, the upper boundary condition is given by interpolation of average weakly air temperatures and average weakly temperature amplitudes. It is assumed that there is no heat flux across the lower boundary in 2 m soil depth. Convective heat flow with water and vapor is neglected in LEACHN; however, it is considered in the other models. In addition, SHAW considers the effects of soil frost and evaporation on the heat storage in soil. In both DAISY and SHAW, it is assumed that soil temperature at zero soil depth equals air temperature. At the bottom boundary, an analytical solution of the heat transfer equation is used assuming that frost/thaw and convective heat transfer can be neglected.
The two approaches used to simulate soil water dynamics were based on a numerical solution to the Richards equation, as implemented in the HYDRUS 1D model [59], and alternatively, on the tipping-bucket (TB) approach of CERES [60]. The former approach is based on the numerical solution of the governing water flow equation using an implicit finite element scheme, which can be used in unsaturated, partially saturated, or fully saturated porous media. It considers the gradients of soil water potential to simulate soil water movement. The tipping- bucket (TB) approach divides the soil into a cascade of ''buckets" and calculates the water balance processes. Drainage occurs only when the field capacity of the overlying bucket is reached.
Since the matric potential is neglected, upward soil water movement (capillary rise) is not simulated. Potential evapotranspiration was simulated based on the standard Penman-Monteith FAO56 method with a single crop coefficient [61] for all model ensemble members. The majority of the crop models participated in previous AgMIP studies were based on Penman-Monteith method [41,62]. Webber et al. [63] showed that the difference in methods to calculate the reference evapotranspiration (ET0) can result in significant uncertainty in yield estimates. However, the Penman-Monteith method is physically the most comprehensive method for calculating evapotranspiration. It consistently provides better estimates than other methods, in particular including the Penman equation [64] provided all input variables are available.

Model calibration
The crop growth models were calibrated by minimizing the sum of squared errors between the model simulated and the observed experimental data. The data used for the calibration were flowering and maturity dates, and crop yields for six years for each cultivar ( Table 1). The calibration procedure involved two steps. First, the genetic parameters (S2 Table) controlling simulated phenology were estimated based on observed flowering and maturity dates. In a second step, the parameters from step 1 were fixed, and the genetic parameters which control crop growth were estimated. For parameter estimation, we used the freely available inverse modelling code UCODE [65].

Statistical analysis
Performance metrics. The performance of the calibration was evaluated using most commonly used statistical methods such as Pearson correlation coefficient (r), Root Mean Square Error (RMSE), the percent bias (PBIAS) and Nash-Sutcliffe efficiency (NSE). As a measure of agreement between simulated (sim) and observed (obs) data, we computed r (Eq 1). A value of r close to 1(-1) shows a perfect positive (negative) correlation between the simulated and observed data.

PBIAS ¼
The NSE (Eq 4) measures the magnitude of the residual variance relative to the measured data variance [67]: where obs and sim refer to observed and simulated data whereas, obsand simrefer the average values respectively. NSE can range from −1 to 1 where NSE = 1 corresponds to a perfect simulated data match, whereas NSE < 0 occurs when the mean of observed data describes the data better than the model simulation.
In addition, the calibration results were visualized in Taylor diagrams [68]. A Taylor diagram summarizes how well model simulation match the observed system behavior to facilitate the relative assessment of different models. For this, it combines three statistics in a 2-D graph: the normalized standard deviations (SD) of both simulated and observed data, Pearson's correlation coefficient between simulated and observed data, and the centered RMSE. The radial distance (blue contours) from the point on the x-axis identified as "observation" is a measure of the centered RMSE between the simulated and observed. The normalized SD of the simulated pattern are proportional to the radial distance from the origin to the points and the cosine of the polar angle (pink) is equal to the correlation between simulated and observed data. The standardization of SD is done using the SD of the respective observational data so that the reference has a standard deviation of 1. The centered RMSE refers that the simulated and observed data are normalized by subtracting the respective means before calculating the RMSE.
Sensitivity analysis. To identify the most influential factors that affect grain yield, a sensitivity analysis was conducted by discretely varying selected input variables and analyzing their impact on modelled yield for the 30 baseline years (1981-2010). We considered five levels of fertilizer applications: 0, 40, 80, 120 and 160 kg ha -1 ; five [CO 2 ] levels: 360, 450, 540, 630 and 720 ppm; four levels of temperature increase: 0, 2, 4, and 6˚C; and five levels of precipitation change: -50, -25, 0, 25, and 50%. The perturbations are applied as absolute offsets for temperature treatments from daily maximum and minimum temperature data whereas precipitation perturbations are applied as fractional changes to daily data. The sensitivity to changes in [CO 2 ], temperature and precipitation were evaluated under three (0, 80, 160 kg ha -1 ) N-fertilizer levels to evaluate their combined effects on grain yield. Here, the main purpose of the sensitivity analysis was to sample over a wide range of treatment space to understand the influential factors that drive grain yield changes. The levels were selected one to encompass the current situation (represented by 360 [CO 2 ], 0˚C temperature increase, and 0% precipitation increase from the baseline climate) and second to approximate as much as possible to the levels used in previous studies [13,41,69] for comparison purpose. Grain yield was simulated for 30 years � 141 perturbations � 48 model ensembles resulted in 67,680 simulations per cultivar. We analyzed the effect of the treatments using Student's t-test and linear modelling with analysis of variance (ANOVA) and post-hoc Tukey test.
Climate change impact calculation. The potential impact of climate change on crop growth (yield and grain [N]) were computed by comparing the relative changes in the simulations driven by baseline (1981-2010) climate and future (2021-2050) climate projections. The relative percentage change (ΔZ) due to climate change was calculated as: where Z future and Z baseline represent the simulated median values for the future and baseline climates, respectively. Dissecting sources of uncertainty. In predicting climate change impact on crop growth, uncertainties may arise for a variety of reasons. In this study, we compared the contributions to the uncertainty in the projected impact of climate change on yield from i) climate models, and ii) Expert-N models and submodules selected to represent soil and crop processes. As in Asseng et al. [17], the uncertainties were quantified by calculating the respective coefficient of variation (CV) for each crop cultivar. For instance, to calculate the CV due to climate models for a certain crop cultivar, first the standard deviation of yield changes over the ten climate models was calculated for each combination of crop and soil models: The same procedures were repeated to calculate the CVs resulting from the use of different plant models, soil water models, SOM models and soil heat models.
While the trials on the two maize cultivars were conducted at one site, the two wheat cultivars were grown at two different sites. Since the maize and wheat sites are not the same, we present our study per cultivar. Specifically, this means that the sensitivity analysis, the climate impact study, and the uncertainty decomposition analysis were conducted per cultivar.
Crop management assumptions. Throughout our analysis, the same crop cultivars and management were assumed to evolve in both baseline and future periods (S3 Table). The planting window was calculated by taking the mean +/-2SD of the planting dates used for model calibration. Planting was triggered on the first occasion from the beginning of planting window with 20 mm or more within a 3-day period and no dry spell exceeding 10 days in the following 30 days [70]. The last day in the particular window will be set as a planting date if the condition was not met. In the present study, the models were run with same initial conditions each year. In addition, the pressure from pest and disease were not considered here. The effects of CO 2 fertilization were not included in our climate change impact assessments.

Model performance
The detail statistics (r, NSE, RMSE and PBIAS) of individual models' performances are presented in S4 and S5 Tables. In the case of phenology (flowering and maturity dates), only the crop growth submodels affect the calibrations (S4 Table and S1 Fig). The RMSEs for the flowering dates were around two and five days for wheat and maize cultivars, respectively, and the corresponding correlation coefficient was higher than 0.90. The NSE values were also reasonably good ranging from 0.62 to 0.90 and from 0.6 to 0.94 for wheat and maize cultivars, respectively, where the PBIAS was below 3.2% in both crops. Similarly, close agreements were observed between the simulated and observed maturity dates. The models simulated maturity dates with RMSE lower than 2 and 5 days for wheat and maize cultivars, respectively, with r values greater than 0.97. The NSE values ranged from 0.53 to 0.96 and 0.8 to 0.98 for wheat and maize cultivars, respectively with a PBIAS value of below 2.0% in both crops. Overall, all simulated phenological stages (flowering and maturity dates) were in good agreement with the observations of both crops and in all cultivars.
The standard deviation, correlation coefficient, and RMSE (relative to the reference observation) for grain yield is shown Fig 4 and S5 Table. The variations due to soil heat transfer submodels were minimal and are not shown in the diagram. Each symbol and color combination represents a plant-soil water-soil organic matter model combination. For example, red circles represent simulation results obtained by running the CERES crop growth model together with the tipping-bucket approach for soil water movement and humus mineralization by DAISY. Other than with the phenological stages, the different ensemble model members showed quite large spread in simulating crop yields. The resulting ensemble members simulated grain yields of wheat and maize with RMSE values ranging from 0.13 t ha -1 to 0.54 t ha -1 and 0.36 t ha -1 to 1.01 t ha -1 , respectively. The NSE values ranged from 0.21 to 0.75 and 0.15 to 0.95 for wheat and maize cultivars, respectively, where PBIAS values range approximately from -2% to +2% and -1.2% to +3.9%, respectively. The respective coefficient of correlation ranges from 0.51 to 0.90 and 0.4 to 0.97 for wheat and maize cultivars, respectively. Among the model ensemble, members that involved the GECROS plant model coupled with the tipping bucket approach perform less in simulating wheat grain yield. In contrast, model ensemble members based on the Richards equation reproduced the observed grain yield significantly better than those based on the tipping bucket approach. This holds for both crops.

Sensitivity analysis
The 'probability-of-exceedance' plots in Fig 5 show the spread of the ensemble runs for the simulated effect on grain yield, dependent on the different nitrogen fertilization application rates. Maize yield continuously increased with increasing fertilizer level, whereas the response of wheat grain yield was pronounced only until 80 kg ha -1 N, after which the effect of additional fertilizer was small. Based on the ensemble median, increasing fertilizer application rates from 0 to 160 kg ha -1 could increase maize grain yield by 4 to 4.4 t ha -1 while the corresponding range of increase in wheat grain yield was 1.0 to 1.7 t ha -1 . The variability in the ensemble was consistently higher in the lower fertilization levels across cultivars. For instance, at zero fertilization level, the standard deviation (SD) of the yields was more than 1.5 t ha -1 for maize and 0.8t ha -1 for wheat. However, the SD were less than 0.7 t ha -1 at the highest fertilization rate (160 kg ha -1 ) for both crops.
The results of the sensitivity analysis of grain yield to different levels of [CO 2 ], temperature and precipitation are presented in box plots under three fertilizer levels (Fig 6). The associated multi-comparison ANOVA test result is also appended (S6 Table). The pattern of the effects of elevated CO 2 were much clearer for wheat than for maize. The response of wheat grain yield to increased [CO 2 ] were positive at all treatment levels and the effects were enhanced at higher fertilizer amounts. For instance, increasing [CO 2 ] levels from 360 to 720 ppm increased median wheat grain yield by 1.2 t ha -1 under high fertilizer rates (160 kg ha -1 ) while the effect was only 0.5 t ha -1 for unfertilized treatment. However, the maize grain yield responded to a change from 360 to 450 ppm [CO 2 ] level, beyond which no further response was observed. The corresponding maize grain yield increase was only 0.2 t ha -1 and the effects also show not much difference between different levels of fertilizer treatments.
Temperature had a clear effect on grain yield and the effect was stronger for wheat. Wheat grain yield decreased by 1.6 to 1.8 t ha -1 for a 6˚C increase in temperature, which conforms with a relative decrease in yield of 47-57%. However, the effects on maize yield were modest. The decrease was less than 0.2 t ha -1 for a 6˚C increase in temperature, which is less than 10%. In particular, there was no considerable difference in maize grain yield for 2˚C temperature increase. The temperature effect on wheat yield was higher when N fertilization was high. A 6˚C increase in temperature decreased the wheat yield by 50-60% under high N fertilization (160 kg ha -1 ) while without fertilizer, the corresponding decrease was 35-50%.
The simulated grain yields reacted rather insensitive to changes in precipitation. In some cases, a negative relationship between increases in precipitation and yield was simulated. A 50% decrease of precipitation resulted in a decrease of median maize yield by around 0.5 t ha -1 (p<0.05) while a 50% increase in precipitation led to no significant change. Wheat yield, on the other hand, show no significant response for precipitation change. However, despite its weak relative strength, wheat yield decreased by about 0.1 t ha -1 for a 50% decrease in precipitation under high N fertilization. There was also an instance where wheat yield decreased (increased) by 0.3 (0.2) t ha -1 for a 50% increase (decrease) in precipitation under unfertilized treatment. These results were opposite in tendency to those observed for the temperature and [CO 2 ] treatments in the sense that the response decreased as N fertilization increased. Further investigation revealed that much of the effects of precipitation perturbation were attributed to N leaching which was relatively higher for wheat cultivars (S2 Fig). An increase in precipitation resulted in higher N leaching where the effects were different at different level of N input. The effects of increased precipitation on N leaching were larger at higher N fertilization rate (160 kg ha -1 ) for wheat cultivars while maize cultivars were more affected when N was limiting (0 kg ha -1 ). The results clearly show an interaction effect of nitrogen and precipitation. The pattern of precipitation effects was consistent between the maize cultivars but not between wheat cultivars.  Fig 7(A) presents the percentage change in grain yield between the baseline and future climate scenraios by the ensemble. Wheat grain yield is expected to decrease by 36 to 40% by 2050 whereas the impact on maize yield is minimal (~-2%). However, the variability of the projected change in wheat grain yield was much higher SD~18-21%. The results were consistent among the cultivars. The grain N concentration is expected to decrease by about 27% for wheat (Fig 7B). In contrast, the simulations suggest that maize grain N concentration will increase by 12 to 14%. The variability associated with the projected change of grain N was higher for wheat cultivars (SD: 17-28%) than for maize cultivars (SD: 7.5-15%). We analyzed the effect on the growing season length, as a definition for the operation window for the farmers. Due to climate change, a substantial shortening of crop life cycle (LGS) is projected for both crops (Fig 7C). A median shortening in LGS of 12 days was found for maize and between 8 to 9 days for wheat. The uncertainty in projected grain yield change caused by the variability in crop growth models was larger compared to the rest of the model components considered in this study. The corresponding median CVs were 40-80% (crop growth submodels), 11-25% (climate models), 3-22% (soil water flow submodels), 1-15% (soil organic matter submodels) and less than 2.5% (soil heat submodels). The uncertainty due to crop growth submodels was larger for wheat than maize. The uncertainty in the wheat grain yield changes was largely dominated by the variations in the crop growth models (with 71-80% CV) followed by climate models (with 11-25% CV). Meanwhile, the sources of uncertainty in the maize grain yield changes were fairly distributed among all model groups except for the soil heat submodels. Accordingly, 39-47% of the variations were attributed to the differences in crop growth submodels followed by climate models, soil water flow and soil organic matter submodels with CVs 19-23%, 19-22% and 10-16%, respectively. The contributions in the uncertainty were largely consistent within the crop cultivars but not across crop species. Among the model components considered, the uncertainty due to variations in the soil heat transfer submodels were minimal (CV<2.5%).

Yield response to elevated CO 2 , warming atmosphere, change in precipitation and fertilizer application
It is established knowledge that elevated CO 2 affects crop production. However, the strength of the response varies widely among crops and environment [71][72][73][74][75][76]. Compared to current yields, our study suggested that elevated CO 2 to 720 ppm would increase wheat grain yield by 35-37% while maize grain yield would increase only by 3 to 5%, without considering the increase by the breeding progress or new technologies. Based on 59 peer-reviewed papers, Wang et al. (2013) conducted a meta-analysis on the impacts of elevated CO 2 and found that wheat yield increased by 20-28% as [CO 2 ] increased from 450 to 800 ppm. Depending on the environment, some studies even reported much higher yield increases of up to 44% [77] and 70% [73]. The effects of elevated CO 2 were lower under low N fertilization treatments. The response of grain yield of C 3 plants for elevated CO 2 was higher under non-limiting N and H 2 O [74]. Regarding maize yield response, our results are consistent with previous modeling studies [13,78]. A doubling of [CO 2 ] from 360 to 720 ppm resulted in an average 7.5% increase in maize grain yield [13]. The lower response of maize yield to increasing [CO 2 ] in our study can be related to the site climate conditions, where rainfall is high. In C 4 plants such as maize, the response of grain yield to increased levels of [CO 2 ] is more pronounced under water stress conditions [72,75], due to water saving by a reduced stomatal conductance. Temperature, along with photosynthetic radiation and humidity are important agroclimatic factors that affect crop growth and development [79]. It affects various photosynthesis-related reactions and hence crop yield [80]. The present study highlights the negative effects of increasing temperature on grain yield. This is in-line with previous findings [81]. Higher temperature generally decreases yield by shortening the grain filling phase. Like the effects of elevated CO 2 , however, the level of the impact on wheat yield was much higher than that on maize yield. A temperature increase by 6˚C from current climate reduces wheat grain yield by 47-57%, but the same temperature increment results in a reduction of only 10% in maize yield. The relative reduction in wheat yield was higher than that presented in Pirttioja et al. [82]. These authors reported 28% to 37% reduction in spring wheat yield over Europe in response to a 6˚C temperature rise. However, in their global analysis, Asseng et al. [17] also found that wheat yield can be decreased by up to 56% in the tropics due to a 4˚C increase. For maize, our results are comparable to the study of Bassu et al. [13]. Evaluating maize crop response to different climate factors using multi-crop models, they reported a modest reduction of grain yield.
The impact of precipitation on grain yield was the weakest in our sensitivity analysis. The effects of a 50% increase/decrease of current precipitation resulted in less than a 0.5 t ha -1 and 0.3 t ha -1 yield difference for maize and wheat. In some cases, precipitation increase had also a small negative effect on grain yield, particularly in case of wheat. These results have to be seen in the light that perturbation of precipitation did not affect the frequency of rainy or dry days but the quantity of rainfall. Guan et al. [83] also found similar decline in crop yield with increasing rainfall under high mean annual precipitation scenario in west Africa. Major cereal production including maize, wheat and sorghum showed no significant correlation with precipitation and should be studied together with other agronomic and climatic factors [84]. On the other hand, the effects of precipitation were smaller when considering increases in N fertilization which is in-line with the fact that application of N increases water use efficiency [85]. The simulated effects of precipitation changes were rather more pronounced in N leaching compared to grain yield. A possible explanation for the negligible effect of precipitation perturbation on grain yield could be attributed to the associated high soil organic carbon content at the study sites [86]. The high water holding capacity of the soils due to high organic carbon contents may provide the sufficient buffering capacity for extended time periods with low amounts of precipitation. More importantly, the low sensitivity could be attributed to the fact that the higher amount of precipitation at the study sites would be satisfactory even with -50% perturbations. Similar results were reported for sites in Africa (like Benin and Ethiopia) with high seasonal precipitation [41].
Results of the sensitivity analysis showed that nitrogen fertilization played a crucial role in modelled crop productivity in the study area. Compared to the unfertilized treatment, application of 160 kg ha -1 nitrogen increased maize and wheat median grain yield by 4.0-4.4 t ha -1 and 1.0-1.7 t ha -1 , respectively, corresponding to a relative increase of 90 to 100% and 40 to 65%, respectively. In their field experiment conducted in a maize growing area in Ethiopia, Abera et al. [87] found up to 254% increase in maize grain yield with 110 kg ha -1 nitrogen application. For wheat, studies highlight that in the highlands of Ethiopia nitrogen application has even higher impact on grain yield [88].

Possible effects of climate change on crop growth
Climate change will negatively affect crop production in Africa [6,89,90]. Our study suggests that future wheat yield will be a major challenge in Ethiopia (Fig 7A). However, the associated variability in the prediction of the grain yield was high. Averaged over all climate simulations, annual temperature is projected to increase by 1.2 to 1.4˚C by 2050 in the study area (Fig 2A  and 2B). This entails a tremendous impact on a cool season crop like wheat through the acceleration of phenological development [62]. Compared to the baseline period, the median length of the wheat life cycle is expected to shorten by eight to nine days by 2050 (Fig 7C). The exposure to climate extremes particularly during the grain filling phase [79] could probably explain the large yield decrease despite modest change in growing season duration which clearly worth further research. In their review of climate change impacts on major crops in eastern Africa, Adhikari et al. [91] highlighted wheat as the crop most vulnerable to climate change in the region. However, the projected yield losses could partly be offset if we included CO 2 fertilization in our simulations, which have been found to have positive effects in C3 crops [76].
Studies suggest that atmospheric CO 2 exposure also affects the availability of N and N use efficiency [80]. Our study indicate that the impact of changing climate also manifested itself in terms of grain quality where a decrease of up to 28% in wheat grain N concentration is projected by 2050 (Fig 7B) but with a high associated variability. Asseng et al. [15] evaluated the effect of climate change on wheat protein around the world using 32 different wheat crop models with combination of 5 GCMs. They found that grain protein concentration would decline as a result of climate change in most parts of Africa. The study suggested that reductions in grain quality were attributed to nitrogen availability limiting growth stimulus from elevated CO 2 particularly in low-rainfall regions [15]. A meta-analysis study supported the result that elevated CO 2 reduces grain quality albeit it stimulates wheat grain yield [92].
In contrast, our estimates showed maize production would not suffer as much as wheat crop with a projected grain yield reduction of less than 2% by 2050. This was even though maize growing season was projected to decrease by up to two weeks under future climate ( Fig  7C). It suggests that changes in the critical development stages like the timing of flowering and length of grain filling period could be more important than the total growing season in determining the final grain yield. The projected increase in seasonal total precipitation for the site (Fig 2C) can also partly explain the minimal impact. The estimated relative change in maize grain yield is in line with the study conducted by Tesfaye et al. [10] who reported a relatively small yield reduction of 3-9% projected for 2050s for the highlands of Eastern Africa. Likewise, in their recent working paper, Timothy et al. [29] pointed out that the impact of climate change on crop yield including maize will be minimal. Some studies even reported that climate change would increase maize yield by 2050s. Araya et al. [30] studied the effects of climate change on maize yield in southwestern Ethiopia and found that median grain yield would increase by up to 4.2%. Generally, the projected impact of climate change on maize production would be minimum in east Africa [10,91,93]. In the lowland part of the country, however, much higher maize yield reduction (up to 20%) was projected by 2050 [32]. In contrast to the grain yield, climate change affected maize grain quality positively, where grain N concentration increased with median projection ranging from about 12% to 14% (Fig 7B).

Sources of uncertainty
Understanding and quantifying different sources of uncertainties has become the core interest among scholars in multi-model climate change impact studies [13,20,62,78,94]. In this regard, sources including the choice of crop models and climate model scenarios contributed to the uncertainty of the simulated results. Comparisons of uncertainties among the different model components showed that differences in crop growth models resulted in largest uncertainties in the crop yield change estimates (Fig 8). Our result is in line with similar previous studies [17,[95][96][97]. In the largest global crop model intercomparison exercise, Asseng et al. [17] compared 30 different crop models in simulating wheat growth around the world under rising temperature. They found that the uncertainty in simulated grain yield due to differences in crop models were more than the uncertainty caused by the variations in locations as well as the inter-annual variability. Similarly, in their climate change impact assessment, Araya et al. [95] reported the choice of crop models to cause much of the variations in projected maize yield.
The present study also shows that the uncertainties vary among crop species and cultivars as well (Fig 8). The uncertainties in the projected yield change due to crop growth models were higher for wheat than maize. This is related to the larger sensitivity of the wheat plant models and the magnitude of impact of climate change on wheat. The uncertainty in simulated wheat yield increased with increasing temperature [17]. The relative inconsistent pattern among the wheat cultivars also suggested that the contributions of different sources of uncertainties might also be site dependent [98]. The variation due to soil water dynamics submodules could be related with the different performance of the two approaches suggesting the importance of soil processes in such studies [99]. Our study showed that soil water dynamic submodule based on Richards equation approach performed better than the tipping bucket based approach in simulating yield. The result suggest that Richards based soil dynamics modeling approach would help in capturing the high variability in the observation data like in the case of ours. Richards equation based approach has a potential for better simulation of soil water dynamics provided that the necessary data are available [100]. Our study highlighted that a considerable amount of uncertainty originated from the different soil water flow and soil organic matter models. However, the variance in projected yield changes due to variations in crop growth models remain the main source of uncertainty across examined crop species and cultivars. This is mainly due to the fact that crop growth models vary largely in their approach to represent key dynamic processes [42] and overall, their predictive skill might still be limited. Model improvement has been hindered by unavailability of high-quality data as well as uncomplete understanding of some plant processes [101].

Conclusions
In the face of climate change, agricultural decision making is and will be a challenging task due to the uncertainties associated with impact assessments. Different sources of uncertainties, including choice of crop models and climate model projections, challenge results of climate change impact assessment. Our results indicate that the use of multiple dynamic crop models provides a better understanding of the uncertainties as a consequence of model structural differences, leading to more robust analyses. The current study also included the uncertainties that could arise from future climate prediction through involving a variety of distinct climate models that were downscaled to finer resolution. In such an effort to quantify, dissect and compare uncertainties, the Expert-N agro-ecosystem modeling package was an ideal framework to facilitate multiple model configurations. Importantly, our results identified crop growth model associated uncertainty to be larger than the rest of model components considered in this study. Uncertainties varied between the crop species and cultivars as well. Incorporating different soil water flow and soil organic matter modeling approaches in multi-model ensemble uncertainty studies could also reveal a considerable portion of uncertainty in grain yield changes. Under the current crop management scenario, our results suggest that maintaining current wheat yield will be a challenging task for the region due to the projected greater impacts by 2050. To abate the associated problems, proper management strategies will be required. In this regard, proper application of N fertilization would alleviate part of the projected impact. Regarding precipitation, our sensitivity study was based on simple perturbations of the quantities to which the crops did not show a significant sensitivity. Possibly, the distribution of rainy days within the year might matter more than the difference in rain amounts in the region. In summary, our study quantifies the impact of climate change and demonstrates the importance of multi-model ensemble approach in understanding the associated uncertainties.