High Sowing Densities in Rainfed Common Beans (Phaseolus vulgaris L.) in Mexican Semi-Arid Highlands under Future Climate Change

Mexico holds the largest single bean production area in the world that is vulnerable to drought. Using field data and two future climate scenarios (RCP4.5 and RCP8.5) for the period 2020–2039, this study evaluated three common bean (Phaseolus vulgaris L.) cultivars planted under rainfed conditions at different densities in two locations in the north-central Mexican semi-arid temperate highlands. The sowing densities were 90,000, 145,000, and 260,000 plants ha−1 established in single rows (SR), three rows (3R), and six rows (6R), respectively. The climate change scenarios were derived from an assembly model integrating 11 general circulation models (GCM) selected for Mexico with a 30” arc resolution. The baseline climate was for the period 1961–2010. The ALMANACMEX model (USDA-ARS-INIFAP, Temple, USA) was parameterized and evaluated and then re-run using the climate scenarios. Beans planted at 6R showed the highest increase in seed yield in both climate scenarios, although the response varied by cultivar and time periods. For the growth habit III cultivars, Flor de Mayo Bajio showed no difference in yield, while Pinto Saltillo, a drought-resistant cultivar, showed increases of 13% to 16% at 6R only until 2033. Growth habit I cultivar Azufrado 2 showed more than 60% increases at 6R in both climate scenarios for the full period 2020–2039. These results suggest that considering the projected climate conditions, high sowing densities may be a viable agronomic option for common bean production under rainfed conditions in semi-arid temperate regions, such as the highlands of Mexico, in the near future; however, the selection of the cultivar is a key element to consider in this regard.


Introduction
Common bean (Phaseolus vulgaris L.) is an essential food crop in Central and South America, Sub-Saharan Africa, and parts of Asia-regions that are predicted to be negatively impacted by climate change in the next few decades [1]. With the predicted decreases in water availability and increases in temperature in major agricultural production areas [2][3][4], sustainable production will be a central issue in agricultural and food systems in these areas. In this regard, common bean and other legumes could play an important role as they provide high-quality food and feed and help reduce the emission of greenhouse gases [5].
However, drought presents a risk to common bean production around the world, with Northeastern Brazil and the central and northern highlands of Mexico being the most vulnerable areas in Latin America [6][7][8]. Drought in relation to beans is defined as inadequate water availability (including precipitation and soil moisture storage capacity) in terms of quantity and distribution during the growing period that prevents the cultivar from achieving its full genetic potential [9]. The highlands of Mexico hold the largest single drought susceptible production area in the world, where more than 1,000,000 hectares of beans are cultivated [10]. Yields fall below 0.4 Mg ha −1 in dry years due to droughts and inadequate capital investment in irrigation in the dry areas, which are also poverty hotspots [10,11].
The use of high sowing density in improving bean production under future climatic scenarios needs to be further explored, especially in semi-arid production systems. Sowing density is considered a major management factor, affecting the growth and development of grain crops by modifying the canopy light environment and interplant competition for water and nutrients [12]. The effects of sowing density on the vegetative and reproductive growth traits of common beans have been explored previously [13][14][15] but mostly under non-limiting management conditions and without considering future climate change. The present study evaluated the use of high sowing densities on common bean cultivars under rainfed conditions in semi-arid temperate highlands in Mexico, considering predicted changes in precipitation and temperatures.

Study Area
The field study was conducted at two research stations of the Instituto Nacional de Investigaciones Forestales, Agricolas y Pecuarias (National Institute of Forestry, Agriculture and Livestock Research) (INIFAP) in Aguascalientes State in the semi-arid temperate highland common bean region of Mexico. The Pabellon de Arteaga research station, located at 22 • 11 latitude north and 102 • 20 longitude west and 1912 m a.s.l., has Calsicol soil with sandy-loam soil texture, 0.48 m depth 2.0% of organic matter and pH of 7.9 [16]; average annual rainfall is 460 mm [17]. The Sandovales research station, located at 22 • 09 latitude north and 102 • 17 west at 2000 m a.s.l., has Calcisol and Planosol soils with sandy-clay soil texture, 0.45 m depth, pH of 6.2, and 1% of organic matter [16]; the average annual rainfall is 400 mm [17].
Both sites have a semi-arid climate, with 308 to 349 mm of mean precipitation during the bean growing season (July to October), mean annual temperature of 16.2 • C, and mean maximum and minimum temperature of 20 • C and 7.1 • C, respectively [17,18].

Cultivars and Sowing Densities
The cultivars used in the study were Pinto Saltillo (PS), Flor de Mayo Bajio (FMB), and Azufrado 2 (A2). PS and FMB are growth habit III cultivars, i.e., indeterminate, with prostrate architecture, while A2 is a growth habit I determinate cultivar, with bush-type architecture [19]. These bean cultivars (Table 1) were from the Bean Genetic Improvement Program of INIFAP. The early maturing cultivars FMB and PS are generally preferred by Mexican farmers [20].
The three cultivars were established during three consecutive agricultural seasons (2011,2012,2013)  At the time of sowing, the seed was inoculated with the Glomus intraradices rhizobia strain at a dose of 350 g ha −1 of micro-rustic substrate. Before sowing, the land was prepared with a multi-plow to break the ground without inverting it, followed by harrowing. The sowing was done using an experimental prototype mechanical seeder designed by INIFAP [21] that enables the three methods of establishment, with distances of 0.76, 0.40, and 0.20 m between rows for SR, 3R, and 6R, respectively. Planting in Pabellon de Arteaga was done on 30 July in 2011 and 31 July in 2012 and 2013. Planting in Sandovales was on August 1 in 2012 and 2013.

Precipitation During the Years of Study
The precipitation during the growing season of rainfed beans (July-October) in the years of study was 175 mm, 208 mm, and 492 mm for the Pabellon de Arteaga research station during 2011, 2012, and 2013, respectively, and 168 mm and 364 mm in Sandovales during 2012 and 2013.
Considering that the July-October mean precipitation of the two semi-arid locations is 308 to 349 mm [17], it was possible to obtain information of the crop under two rainfall regimes: dry (<308 mm) in the two locations (175 and 208 mm in Pabellon de Arteaga during 2011 and 2012 and 168 mm in Sandovales during 2012), and wet (>349 mm) in both locations during 2013 (492 and 364 mm for Pabellon de Arteaga and Sandovales, respectively). Modeling the behavior of the crop under these two common water conditions in rainfed areas gives greater certainty to the parameters in the parameterization process and in the evaluation of the model.

Field Measurements
The leaf area index (LAI) was measured in Pabellon de Arteaga during the 2011 growing season, using a 0.8 m long LP-80 ceptometer (Decagon, Pullman, WA, USA). Three fixed points were randomly established per experimental unit consisting of four, six, and 12 rows of 30 m length per cultivar with three repeats for SR, 3R, and 6R, respectively. The readings were made between 10:00 and 13:00, during which time the incidence of solar radiation is relatively stable [22]. As the evaluated cultivars differ in development and precocity (Table 1), a calendar of measurements was established considering days after sowing (DDS); this made it possible to compare LAI in the same period of plant development. Samples were taken at 19,25,40,45,79, and 87 DDS.
Seed yield was evaluated at the end of the growing season of each cultivar. To obtain the yield data (Mg ha −1 ), plots of 1.5 m width and 2.0 m length were randomly selected from the two central rows of each experimental unit.

ALMANACMEX Model
The study used the ALMANACMEX model (USDA-ARS-INIFAP, Temple, USA) [23], which is a version of the ALMANAC model [24] that includes a Mexican interface with data on the climates and soils of Mexico. The model simulates processes of crop growth and soil water balance including light interception by leaves, dry matter production, and the partitioning of biomass into grain. Light interception is simulated by Beer's Law and considers the total leaf area and height of the canopy [24]. The water and nutrient balance subroutines are from the Erosion Productivity Impact Calculator (EPIC) model [25]. The simulated water demand is critical for yield simulation in water-limited conditions. Potential evaporation (E o ) is calculated first. Then, the potential soil water evaporation (E S ) and potential plant water transpiration (E P ) are derived from the potential evaporation and leaf area index (LAI). In addition, the model estimates stresses caused by water, nutrients (N and P), temperature, and aeration. The water stress factor is computed by considering supply and demand in the equation: where WS is the water stress factor, u is the water use in layer l, and Ep is the potential plant water use on day i. This is consistent with the concept that drought stress limits biomass production in proportion to transpiration reduction [26].

Model Parameterization
The first step in the modeling process was parameterization, i.e., finding parameters with the best correlation between the model and the reality being simulated [27]. The initial selection of parameters was refined during model calibration by running the model with the parameters and comparing the model results with observed data [28] so as to determine the best fit for the regional environment of interest [29]. This process has been used to derive parameters that are not available for important agricultural areas or crops [30][31][32].
In this study, the model was parameterized using two sources of crop parameters: (1) the bean crop parameters already available in the model (radiation use efficiency, plant height, etc.) and (2) the parameters reported for some common beans grown in the same study region (Aguascalientes) [33]: potential LAI, potential Harvest Index (HI), and Potential Heat Units (PHUs). Once the parameters had been incorporated into the model, it was run considering the climatic and edaphic characteristics of the study area of Pabellon de Arteaga. Daily weather data were obtained from the automated weather station closest to the plot (100 m). The soil data considered were soil depth, organic matter, soil texture, and pH as previously described in the Study Area section.
The model was run using these initial parameters together with observed daily climatic data (maximum and minimum temperature and precipitation) of the year 2011. The model results (simulated LAI and seed yield) were compared with the observed LAI and seed yield during the year 2011.

Model Evaluation
Model evaluation is defined as the process of demonstrating that, under certain conditions, the model is capable of performing simulations with sufficient precision [34].
Model evaluation was performed considering the seed yield and daily climatic conditions (maximum and minimum temperature and precipitation) observed during 2012 and 2013 in the study areas of Pabellon de Arteaga and Sandovales under rainfed conditions. The weather stations where daily weather data were gathered were located 100 to 150 m from the experimental plot on both sites.
The statistical indices mean bias error (MBE), root mean square error (RMSE), and the agreement index (d) [35,36] were used for the evaluation of the simulation performance of ALMANACMEX. The MBE of seed yield was calculated considering the difference between observed and simulated seed yield, using the following equation: where N denotes the observation numbers used for comparisons, and Si and Oi are the simulated and observed seed yield. The root mean square error (RMSE) was used to determine the difference between observed and simulated seed yield. The RMSE was computed by using Equation (3) to determine the degree of predictability [37]: where N denotes the observation numbers used for comparisons and Si and Oi are the simulated and observed seed yield.
The computed values of RMSE and MBE determined the degree of agreement between the predicted values and their respective observed values. To determine the degree of model prediction error, the index of agreement (d) [36] was used in this study, applying the following equation: where O and S are the observed and simulated yield, respectively. The numerator is the mean square error. The denominator is related to the variability in the observed and simulated values [35]. The d values range from 0 to 1; the closer the d index value is to unit (1), the better the fit and model simulation.
A linear regression analysis was also done using GraphPad Prism 5 ® (GraphPad Software, Inc. San Diego, USA) statistic software to describe, during model evaluation, the relationship between simulated and observed yields for the three bean cultivars grown in the two sites under study. In addition, the differences between observed and simulated seed yield were assessed using Fisher's paired t-test.

Future Climate Change Scenarios
Climate change scenarios for 2020-2039 were used in the analysis to determine the short-and medium-term impact of climate change. This period was selected because of the need to generate timely options that could be implemented early by Mexican producers to ensure production stability in the face of projected climatic changes in the next two decades.
For the estimation of future climate scenarios, we used a database with monthly anomalies of maximum temperature, minimum temperature, and monthly average precipitation derived from an assembly model that integrated 11 general circulation models (GCM) (BCC-CSM1-1, CCSM4, GISS-E2-R, HadGEM2-AO, HadGEM2-ES, IPSLCM5A-LR, MIROC-ESM-CHEM, MIROC-ESM, MIROC5, MRI-CGCM3, NorESM1-M), reduced in scale and calibrated for Mexico for the representative concentration pathways (RCP) of greenhouse gases (GHG) 4.5 and 8.5 for the period 2020-2039. The GCMs used are part of the climate modeling developed by the project on inter-comparison of coupled models Phase 5 (CMIP5) [38]. All the information was obtained from the WorldClim climate change data portal using raster images with a spatial resolution of 30"arc [38]. For the baseline climate (1961-2010), the monthly data of the minimum, maximum average temperatures, and the monthly precipitation were obtained for each experimental site from the reference climatology database [39].
The ALMANACMEX model contains a weather database with the monthly weather station parameters of more than 3000 stations located in Mexico [23]. This database is used by the weather generator component of the model to generate daily weather data. The weather parameters of average monthly maximum temperature, average monthly minimum temperature, and average monthly precipitation were manually updated by introducing the monthly anomalies projected for the period 2020-2039 for the weather stations located in Pabellon de Arteaga and Sandovales to reflect the predicted 20-year average from 2020 to 2039 for the RCP4.5 and RCP8.5 scenarios. All other weather and wind variables remained unchanged. The management and soil data input remained the same.

Model Application to Determine Sowing Density Impact
After the ALMANACMEX model had been parameterized and evaluated, the model was re-run using the climate scenarios while keeping constant the sowing methods and management conditions for SR, 3R, and 6R and the soil type of each site. A similar approach has been used in other climate impact studies [40][41][42].
The impact of sowing density in terms of increases or decreases in common bean seed yield during the period 2020-2039 was determined by calculating the difference between the baseline simulated yield, which was the yearly simulated yield under SR, and the yield under plant density 3R and 6R.

Crop Model Parameterization and Crop Parameters
In the parameterization process, the simulated LAI was not different (p ≥ 0.05) from that observed in SR, 3R, and 6R at 40, 45, 79, and 87 DDS for any of the three cultivars under study, as shown in Figure 1 and in the values presented in the supplementary (Table S1). Similarly, the simulated seed yield was not different (p ≥ 0.05) from the observed seed yield under the three sowing conditions ( Table 2).
The model parameterization process generated a set of key crop parameters for the three tested cultivars under study ( Table 3). The key parameters observed among the cultivars were the biomass-energy ratio; optimal, minimum, and maximum temperature; maximum height and root depth; dry matter; decline rate; and harvest index. Table 3 shows other parameters that are linked to LAI and that reflect the characteristics of the cultivar and its response to sowing density.
The PHUs responded to the cultivar cultivation growing cycle with 1100 to 1200 PHUs for the early maturing cultivars (FMB and PS) and 1600 for the intermediate cultivar (A2).
While ALMANAC simulates the grain yield using the Harvest Index (HI), it also simulates the decrease in HI with drought stress near anthesis. The model assumes that anthesis occurs the same date as when the LAI is at its maximum. This occurs in the simulation when the DLAI (fraction of growing season when leaf area declines) fraction of PHUs is reached. The model makes HI sensitive to drought stress when the stress occurs between 5% of the heat units before DLAI until 10% of the heat units after. Five days of severe drought stress during this period decreases the HI to the input value for minimum HI (WSYF). For this study, HI was assumed to be 0.36 for all varieties and the WSYF was 0.23 [33]. With DLAI set to 0.77, as it is for some of the varieties in this study, then 5 days of severe drought stress occurring between 0.72 and 0.87 of PHU decreases the simulated HI from 0.36 to 0.23. Thus, each day of severe drought stress during this interval decreases the simulated HI by 0.026.
Agronomy 2020, 10, x FOR PEER REVIEW 6 of 19 simulated yield, which was the yearly simulated yield under SR, and the yield under plant density 3R and 6R.

Crop Model Parameterization and Crop Parameters
In the parameterization process, the simulated LAI was not different (p ≥ 0.05) from that observed in SR, 3R, and 6R at 40, 45, 79, and 87 DDS for any of the three cultivars under study, as shown in Figure 1 and in the values presented in the supplementary (Table S1). Similarly, the simulated seed yield was not different (p ≥ 0.05) from the observed seed yield under the three sowing conditions ( Table 2).
The model parameterization process generated a set of key crop parameters for the three tested cultivars under study ( Table 3). The key parameters observed among the cultivars were the biomass-energy ratio; optimal, minimum, and maximum temperature; maximum height and root depth; dry matter; decline rate; and harvest index. Table 3 shows other parameters that are linked to LAI and that reflect the characteristics of the cultivar and its response to sowing density.
The PHUs responded to the cultivar cultivation growing cycle with 1100 to 1200 PHUs for the early maturing cultivars (FMB and PS) and 1600 for the intermediate cultivar (A2).
While ALMANAC simulates the grain yield using the Harvest Index (HI), it also simulates the decrease in HI with drought stress near anthesis. The model assumes that anthesis occurs the same date as when the LAI is at its maximum. This occurs in the simulation when the DLAI (fraction of growing season when leaf area declines) fraction of PHUs is reached. The model makes HI sensitive to drought stress when the stress occurs between 5% of the heat units before DLAI until 10% of the heat units after. Five days of severe drought stress during this period decreases the HI to the input value for minimum HI (WSYF). For this study, HI was assumed to be 0.36 for all varieties and the WSYF was 0.23 [33]. With DLAI set to 0.77, as it is for some of the varieties in this study, then 5 days of severe drought stress occurring between 0.72 and 0.87 of PHU decreases the simulated HI from 0.36 to 0.23. Thus, each day of severe drought stress during this interval decreases the simulated HI by 0.026.

Model Evaluation
When the optimized parameters were used to simulate the yield of each cultivar in the two sites under study and under the two conditions of precipitation, i.e., dry and wet years, the MBE of the model was −0.1 for both FMB and PS and −0.2 for AZ2 (negative values mean underestimation). The RMSE was 0.3 Mg ha −1 for FM and PS and 0.4 Mg ha −1 for AZ2. The d-index was 0.9 for the three cultivars. As shown in Table 4 and Figure 2, there was no significant difference (p ≥ 0.05) between the simulated and the observed yield of the three cultivars during the seasons evaluated. In addition, the slope and intercept of the regression line shown in Figure 2 were not significantly different from 1 and 0, respectively.

Model Evaluation
When the optimized parameters were used to simulate the yield of each cultivar in the two sites under study and under the two conditions of precipitation, i.e., dry and wet years, the MBE of the model was −0.1 for both FMB and PS and −0.2 for AZ2 (negative values mean underestimation). The RMSE was 0.3 Mg ha −1 for FM and PS and 0.4 Mg ha −1 for AZ2. The d-index was 0.9 for the three cultivars. As shown in Table 4 and Figure 2, there was no significant difference (p ≥ 0.05) between the simulated and the observed yield of the three cultivars during the seasons evaluated. In addition, the slope and intercept of the regression line shown in Figure 2 were not significantly different from 1 and 0, respectively.    (Table 4) for FMB, PS, and A2 cultivars, respectively, likewise show no significant differences (p ≥ 0.05) between the simulated and observed seed yield of the three cultivars.

Projected Changes in Temperature and Precipitation
It is predicted that climate change will result in rising temperatures and changes in precipitation patterns [43]. As the current and projected climatic factors are temporally and spatially variable in dryland cereal production systems in different countries [44], it is necessary to determine the impact of climate change in each region. Figure 3 shows the increases of maximum and minimum temperature and precipitation projected in the RCP4.5 and RCP8.5 scenarios for Pabellon de Arteaga and Sandovales. Under these scenarios, the increases in projected temperatures during the common bean growing season (July-October) are 1.5 • and 1.3 • C for the maximum and minimum temperature, respectively, for Pabellon de Arteaga and 2.5 • and 2.3 • C, respectively, for Sandovales. This means that the increases will be greater by 1 • C in Sandovales (Figure 3). The higher minimum temperature will be a cause of concern because common beans are sensitive to high night-time temperatures during reproductive development [2]. A 1 • C increase in temperature would be equivalent to a 150 km northward shift of isotherms (lines joining places with similar temperature) or about 150 m lower altitude [45,46]. The calculated paired t test of 1.42, 1.08, and 1.99 with 11 degrees of freedom and the Pearson correlation coefficients of 0.94, 0.91, and 0.86 (Table 4) for FMB, PS, and A2 cultivars, respectively, likewise show no significant differences (p ≥ 0.05) between the simulated and observed seed yield of the three cultivars.

Projected Changes in Temperature and Precipitation
It is predicted that climate change will result in rising temperatures and changes in precipitation patterns [43]. As the current and projected climatic factors are temporally and spatially variable in dryland cereal production systems in different countries [44], it is necessary to determine the impact of climate change in each region. Figure 3 shows the increases of maximum and minimum temperature and precipitation projected in the RCP4.5 and RCP8.5 scenarios for Pabellon de Arteaga and Sandovales. Under these scenarios, the increases in projected temperatures during the common bean growing season (July-October) are 1.5° and 1.3 °C for the maximum and minimum temperature, respectively, for Pabellon de Arteaga and 2.5° and 2.3 °C, respectively, for Sandovales. This means that the increases will be greater by 1 °C in Sandovales (Figure 3). The higher minimum temperature will be a cause of concern because common beans are sensitive to high night-time temperatures during reproductive development [2]. A 1 °C increase in temperature would be equivalent to a 150 km northward shift of isotherms (lines joining places with similar temperature) or about 150 m lower altitude [45,46]. As for the predicted precipitation, the analysis of monthly projections shows decreases and increases of precipitation during the growing season in the two locations. During the initial months of the growth period, i.e., July and August, precipitation will be reduced by 5 and 20 mm for Pabellon de Arteaga and 9.0 mm and 31 mm for Sandovales, under scenarios RCP4.5 and RCP8.5, respectively. However, toward the end of the growing season, i.e., September and October, increases of 4 mm and 40 mm are projected for Pabellon de Arteaga and 5 mm and 45 mm for Sandovales, under RCP4.5 and RCP8.5 climate scenarios, respectively. As for the predicted precipitation, the analysis of monthly projections shows decreases and increases of precipitation during the growing season in the two locations. During the initial months of the growth period, i.e., July and August, precipitation will be reduced by 5 and 20 mm for Pabellon de Arteaga and 9.0 mm and 31 mm for Sandovales, under scenarios RCP4.5 and RCP8.5, respectively. However, toward the end of the growing season, i.e., September and October, increases of 4 mm and 40 mm are projected for Pabellon de Arteaga and 5 mm and 45 mm for Sandovales, under RCP4.5 and RCP8.5 climate scenarios, respectively.

Impact of Sowing Density on Seed Yield
Simulated seed yield of the cultivars under study showed different responses to the sowing method in the two climatic scenarios. For the FMB cultivar (growth habit III), no difference was found in seed yield using the three sowing methods under the two climatic scenarios (Figure 4). The yields of SR versus 3R and SR versus 6R were similar (p ≥ 0.05) in both climate scenarios and in both locations (Figure 4). Table 5 shows that the average increments of seeds (%) in relation to SR projected for the period 2020-2039 for FMB were from 1% to 2% for 3R and from 1% to 4% in 6R in the two locations (Table 5).

Impact of Sowing Density on Seed Yield
Simulated seed yield of the cultivars under study showed different responses to the sowing method in the two climatic scenarios. For the FMB cultivar (growth habit III), no difference was found in seed yield using the three sowing methods under the two climatic scenarios (Figure 4). The yields of SR versus 3R and SR versus 6R were similar (p ≥ 0.05) in both climate scenarios and in both locations (Figure 4). Table 5 shows that the average increments of seeds (%) in relation to SR projected for the period 2020-2039 for FMB were from 1% to 2% for 3R and from 1% to 4% in 6R in the two locations (Table 5).    As for the PS cultivar (growth habit III), which is considered a drought-tolerant cultivar [47], seed yield in both locations increased when the sowing method was changed from SR to 3R and 6R ( Figure 5). As seen in Table 5, the potential increases in seed production within the period 2020-2039 were from 8% to 11% when changing from SR to 3R and from 13% to 16% when changing from SR to 6R in the two locations (Table 5). Rosales-Serna et al. [48] mention that with common bean cultivars such as PS that are resistant to drought stress, maturity acceleration coupled with a high seed filling rate helps lessen the impact of drought.
As for the PS cultivar (growth habit III), which is considered a drought-tolerant cultivar [47], seed yield in both locations increased when the sowing method was changed from SR to 3R and 6R ( Figure 5). As seen in Table 5, the potential increases in seed production within the period 2020-2039 were from 8% to 11% when changing from SR to 3R and from 13% to 16% when changing from SR to 6R in the two locations (Table 5). Rosales-Serna et al. [48] mention that with common bean cultivars such as PS that are resistant to drought stress, maturity acceleration coupled with a high seed filling rate helps lessen the impact of drought. PS cultivar production presented less variability in 6R than in 3R ( Figure 5). However, this was seen only in the first 11 and 14 years within the period 2020-2039 in Pabellon de Arteaga and Sandovales, respectively. Starting from 2030 and 2033, in Pabellon and Sandovales, respectively, both the 3R and 6R methods presented the same variability in production as the SR method ( Figure  5). This is in agreement with Folberth et al. [49], who mention that crop management may be the main factor in obtaining yield targets and sustaining agroecosystems but only until a certain time (2060s).
In order to analyze the response of the PS cultivar to the treatments, we made further analysis of the water available for cultivation in the study area of Pabellon de Arteaga for 2024 (the year with high production under 6R treatment compared to SR) and for 2030 (the year in which the treatments were similar in production) ( Figure 5) under the climatic scenario RCP4.5. For water balance (Figures 6 and 7), we used the information on rain (mm), evapotranspiration (mm), and surface runoff (mm) simulated daily by the model during the two years. Likewise, the simulated LAI was used daily for each treatment.
In Figure 6A,B, it is observed that during the 2024 growing season, there was an intermittent water deficiency (less than 0) with the most severe occurring during the period from 33 to 43 days after sowing (DAS) and from 49 to 51 DAS when the cultivar had not yet reached maximum LAI ( Figure 6B). However, if we compare the LAI during one of the critical days of water deficiency, on day 33 after planting, the SR treatment had 1.1 LAI while 3R and 6R had 0.5 and 1.40 LAI, respectively. Apparently, water stress had no effect on LAI development in the 6R treatment since it PS cultivar production presented less variability in 6R than in 3R ( Figure 5). However, this was seen only in the first 11 and 14 years within the period 2020-2039 in Pabellon de Arteaga and Sandovales, respectively. Starting from 2030 and 2033, in Pabellon and Sandovales, respectively, both the 3R and 6R methods presented the same variability in production as the SR method ( Figure 5). This is in agreement with Folberth et al. [49], who mention that crop management may be the main factor in obtaining yield targets and sustaining agroecosystems but only until a certain time (2060s).
In order to analyze the response of the PS cultivar to the treatments, we made further analysis of the water available for cultivation in the study area of Pabellon de Arteaga for 2024 (the year with high production under 6R treatment compared to SR) and for 2030 (the year in which the treatments were similar in production) ( Figure 5) under the climatic scenario RCP4.5. For water balance (Figures 6 and 7), we used the information on rain (mm), evapotranspiration (mm), and surface runoff (mm) simulated daily by the model during the two years. Likewise, the simulated LAI was used daily for each treatment.
during the year 2030. Starting from this year, the planting density using the cultivar PS did not present any advantage under projected changes of precipitation and temperature conditions for the highlands. Figure 7A shows that the cultivar had moisture deficiencies throughout the growing season affecting the LAI development ( Figure 7B) and seed yield. The total days with stress were 29, 26, and 53 for SR, 3R, and 6R, respectively, with seed yields (Mg ha −1 ) of 2.1 for SR, 2.5 for 3R, and 2.3 for 6R.  present any advantage under projected changes of precipitation and temperature conditions for the highlands. Figure 7A shows that the cultivar had moisture deficiencies throughout the growing season affecting the LAI development ( Figure 7B) and seed yield. The total days with stress were 29, 26, and 53 for SR, 3R, and 6R, respectively, with seed yields (Mg ha −1 ) of 2.1 for SR, 2.5 for 3R, and 2.3 for 6R.  In Figure 6A,B, it is observed that during the 2024 growing season, there was an intermittent water deficiency (less than 0) with the most severe occurring during the period from 33 to 43 days after sowing (DAS) and from 49 to 51 DAS when the cultivar had not yet reached maximum LAI ( Figure 6B). However, if we compare the LAI during one of the critical days of water deficiency, on day 33 after planting, the SR treatment had 1.1 LAI while 3R and 6R had 0.5 and 1.40 LAI, respectively. Apparently, water stress had no effect on LAI development in the 6R treatment since it showed a maximum LAI of 2.5 in contrast to 1.4 of SR. As mentioned earlier, one of the qualities of the model is that it calculates the days of water stress. During the whole cycle (2024), the total days of water stress were 39 for SR and 3R and 45 for 6R, and the seed yields were 1.9, 1.8, and 2.3 Mg ha −1 for SR, 3R and 6R, respectively. This demonstrates that high yields can be obtained without additional water and also that evapotranspiration does not change in higher plant populations if the water supply is adequate, as mentioned by some authors [50]. Figure 7A,B shows the water balance during the year 2030. Starting from this year, the planting density using the cultivar PS did not present any advantage under projected changes of precipitation and temperature conditions for the highlands. Figure 7A shows that the cultivar had moisture deficiencies throughout the growing season affecting the LAI development ( Figure 7B) and seed yield. The total days with stress were 29, 26, and 53 for SR, 3R, and 6R, respectively, with seed yields (Mg ha −1 ) of 2.1 for SR, 2.5 for 3R, and 2.3 for 6R.
In the case of the AZ2 (Figure 8), increases in seed yield were over 60% in 6R compared to SR in both locations for the full period of 2020-2039 under the two climate scenarios (Table 5). This indicates that AZ2 at 6R (260,000 plants ha −1 ) seems the most viable option for production under rainy conditions in the study areas, considering the future climatic scenarios. rainfed common bean cultivar Pinto Saltillo (growth habit III) under sowing densities of 90,000 (SR), 145,000 (3R) and 260,000 (6R) plants ha −1 in the semi-arid temperate highland of Pabellon de Arteaga in Aguascalientes, Mexico.
In the case of the AZ2 (Figure 8), increases in seed yield were over 60% in 6R compared to SR in both locations for the full period of 2020-2039 under the two climate scenarios (Table 5). This indicates that AZ2 at 6R (260,000 plants ha −1 ) seems the most viable option for production under rainy conditions in the study areas, considering the future climatic scenarios. Figure 9A,B shows the water balance in 2032, the year which presented the greatest increase in production for 6R ( Figure 8) compared to SR (2.7 versus 1.5 Mg ha −1 ). During almost the entire growing season, there were periods of water deficit ( Figure 9A). The total water stress days were 34, 41, and 46 for SR, 3R, and 6R, respectively. Figure 9B shows that the most critical water stress for the 6R treatment was recorded during days 37 to 45 after sowing, with available water values −9.6 to 0.0 mm, during which time the LAI reached its maximum (3.5) and subsequently declined. However, even with this abrupt stop to LAI development, the final production of 6R was 180% higher than that of the SR treatment.
The performance of this cultivar AZ2 under high sowing density could be associated with the high significant correlation between yield and seeds per plant, seeds per pod, and pods per plant, as has been previously reported [51,52].   Figure 9A,B shows the water balance in 2032, the year which presented the greatest increase in production for 6R ( Figure 8) compared to SR (2.7 versus 1.5 Mg ha −1 ). During almost the entire growing season, there were periods of water deficit ( Figure 9A). The total water stress days were 34, 41, and 46 for SR, 3R, and 6R, respectively. Figure 9B shows that the most critical water stress for the 6R treatment was recorded during days 37 to 45 after sowing, with available water values −9.6 to 0.0 mm, during which time the LAI reached its maximum (3.5) and subsequently declined. However, even with this abrupt stop to LAI development, the final production of 6R was 180% higher than that of the SR treatment. Our results show that the use of high sowing density in common beans as a means of adapting to climatic change in semi-arid temperate regions can be a viable agronomic option, but results may vary according to the growth habit of the cultivar. This is in agreement with the findings of Nienhuis and Singh [53]. Our results show that growth habit III cultivars such as PS, which was developed and released for drought-stressed environments in the northern highlands [49], will respond favorably to the change in sowing density only until a certain time. However, after 2033, with more extreme weather conditions projected in the localities under study, the use of a bush-type architecture cultivar (growth habit I), such as AZ2, at high sowing densities (260,000 plants ha −1 ) is viable for these semi-arid temperate regions.
The architecture of the cultivar and its response in terms of LAI and intercepted radiation will play an important role in future management strategies, such as sowing densities. Ricaurte et al. [12] reported that growth habit I cultivars had approximately 30% lower LAI in relation to growth habit III cultivars but only under lower densities. On the other hand, Gardiner et al. [54], reporting on the behavior of a bush cultivar in different sowing densities, found that leaf abscission and the reduction of LAI were greater for the 45.5-cm row spacing in relation to the 91-cm row spacing from 54 to 72 days after planting, but narrow spacing still intercepted a greater percentage of the photosynthetically active radiation (PAR) at 72 days.
The architecture of the root also plays an important role in the response of the cultivars under drought conditions. Early investigations at CIAT (International Center for Tropical Agriculture) have shown that root architecture differs among common bean cultivars exposed to drought stress under different soil conditions; deeper basal roots are considered advantageous under water stress conditions [2]. Despite being an essential component in combating stress, roots have not been given enough attention in bean breeding efforts. In the present study, the measurement of the root in the cultivars under study was not considered. Therefore, this component should be considered in future research work related to the responses of cultivars to high sowing densities.

Conclusions
The physiological parameters generated for the three bean cultivars under different sowing densities were adequate to simulate their behavior under rainfed conditions. For the growth habit III cultivars, Flor de Mayo Bajio obtained increases in seed yield of less than 2% under climatic scenarios RCP4.5 and RCP8.5, while Pinto Saltillo showed increases ranging from 8% to 16% when The performance of this cultivar AZ2 under high sowing density could be associated with the high significant correlation between yield and seeds per plant, seeds per pod, and pods per plant, as has been previously reported [51,52].
Our results show that the use of high sowing density in common beans as a means of adapting to climatic change in semi-arid temperate regions can be a viable agronomic option, but results may vary according to the growth habit of the cultivar. This is in agreement with the findings of Nienhuis and Singh [53]. Our results show that growth habit III cultivars such as PS, which was developed and released for drought-stressed environments in the northern highlands [49], will respond favorably to the change in sowing density only until a certain time. However, after 2033, with more extreme weather conditions projected in the localities under study, the use of a bush-type architecture cultivar (growth habit I), such as AZ2, at high sowing densities (260,000 plants ha −1 ) is viable for these semi-arid temperate regions.
The architecture of the cultivar and its response in terms of LAI and intercepted radiation will play an important role in future management strategies, such as sowing densities. Ricaurte et al. [12] reported that growth habit I cultivars had approximately 30% lower LAI in relation to growth habit III cultivars but only under lower densities. On the other hand, Gardiner et al. [54], reporting on the behavior of a bush cultivar in different sowing densities, found that leaf abscission and the reduction of LAI were greater for the 45.5-cm row spacing in relation to the 91-cm row spacing from 54 to 72 days after planting, but narrow spacing still intercepted a greater percentage of the photosynthetically active radiation (PAR) at 72 days.
The architecture of the root also plays an important role in the response of the cultivars under drought conditions. Early investigations at CIAT (International Center for Tropical Agriculture) have shown that root architecture differs among common bean cultivars exposed to drought stress under different soil conditions; deeper basal roots are considered advantageous under water stress conditions [2]. Despite being an essential component in combating stress, roots have not been given enough attention in bean breeding efforts. In the present study, the measurement of the root in the cultivars under study was not considered. Therefore, this component should be considered in future research work related to the responses of cultivars to high sowing densities.

Conclusions
The physiological parameters generated for the three bean cultivars under different sowing densities were adequate to simulate their behavior under rainfed conditions. For the growth habit III cultivars, Flor de Mayo Bajio obtained increases in seed yield of less than 2% under climatic scenarios RCP4.5 and RCP8.5, while Pinto Saltillo showed increases ranging from 8% to 16% when planted at high densities (140,000 or 260,000 plants ha −1 ), although the positive response to this agronomic management was only until 2033. On the other hand, habit I growth cultivar Azufrado 2 planted at the highest densities (260,000 plants ha −1 ) showed seed increases of 60% even in the most critical scenario, i.e., RCP8.5, during the full period 2020-2039. These results suggest that considering the projected climate conditions, high sowing densities may be a viable agronomic option for common bean production under rainfed conditions in semi-arid temperate regions, such as the highlands of Mexico, in the near future; however, the selection of the cultivar is a key element to consider in this regard.