Adaptation of Winter Wheat Cultivars to Di ﬀ erent Environments: A Case Study in Poland

: A proper understanding of cultivar adaptation to di ﬀ erent environments is of great relevance in agronomy and plant breeding. As wheat is the most important crop in Poland, with a total of about 22% of the total sown area, the study of its performance in environments with di ﬀ erent productivity levels for consequent cultivar recommendation is of major importance. In this paper, we assess the relative performance of winter wheat cultivars in environments with di ﬀ erent productivity and propose a method for cultivar recommendation, by considering the information of environmental conditions and drought stress. This is performed in the following steps: (1) calculation of expected wheat productivity, depending on environmental factors, (2) calculation of relative productivity of cultivars in the environments, and (3) recommendation of cultivars of a speciﬁc type and range of adaptation. Soil and weather conditions were conﬁrmed as the most important factors a ﬀ ecting winter wheat yield. The weather factors should be considered rather in shorter (e.g., 10 day) than longer (e.g., 60 day) time periods and in relation to growth stages. The ANCOVA model with genotype and management intensity as ﬁxed factors, and soil and weather parameters as covariates was proposed to assess the expected wheat productivity in particular environments and the expected performance of each genotype (cultivar). The recommendation of cultivars for locations of speciﬁed productivity was proposed based on the di ﬀ erence between the expected cultivar yield and the mean wheat productivity, and compared with the Polish o ﬃ cial cultivar recommendation list.


Introduction
Wheat is one of the most important crops in the world and the third largest crop in terms of production after sugarcane and maize and before rice [1]. It is the most important crop in Poland, where it occupies about 22% of the sown area [2]. In Poland, winter wheat is grown in environments of variable productivity, at both regional [2] and local scale. The variable productivity within environments is mainly conditioned by soil, which is relatively stable in time. However, some local soil variability may occur even within a single field [3]. Studies on a wide adaptation, even of the best varieties, show that it may not be sufficient in environments with soils of extremely variable productivity [4]. In such cases, it may be justified to select different cultivars of narrower adaptation to be grown on soils of contrasting productivity, e.g., one cultivar better adapted to fertile soil and another adapted to unfertile soil, and sow them on respective fields or parts of one field. Consequently, farmers need knowledge on which cultivars adapt to the soils of a given productivity: poor, medium or high.

Yield Dataset and Environmental Variables
The data used in this study comprises 19 experimental sites of the Research Centre of Cultivar Testing (COBORU) in Poland (Table 1), observed in the in the Post Registration Variety Testing System (PVTS) during the 2017-2018 cropping season [12,13]. The research stations are located in a warm temperate climate zone, between sea and continental transition. These areas fit the Köppen-Geiger climate classifications of Cfb (warm temperate climate, fully humid with hot summer) and Dfb (snowy climate, fully humid with warm summer) [23].
The winter wheat grain yield in each location calculated as the mean yield of cultivars considered and two crop management levels is shown in Table 1. The environmental, weather and soil variables are shown in Tables 1 and 2. Arable land suitability groups (agricultural soil suitability complexes) were transformed to points in the way described by Witek et al. [24] and used in the paper to represent soil conditions (Tables 1 and 2). Soil pH values were within the range of 4.9-7.0 ( Table 1).
The winter wheat grain yield data covered a subset of 29 winter wheat cultivars, which were evaluated in the growing season 2017-2018 at each trial location at two levels of crop management intensity (for the cultivar names see Table 3). The lower, a moderate-input management system (MIM) included standard fertilization adapted to local conditions, seed treatment, and herbicide and insecticide applications, according to necessity. The high-input management system (HIM) included nitrogen dose increased by 40 kg ha −1 , an additional applications of foliar fertilizers, fungicides and growth regulators [25,26]. Each field experiment was conducted according to a two-factor strip-plot design with two replications and the grain yield for the cultivars was calculated. The mean air temperature (T) and precipitation (P) were provided by the COBORU experimental station and used to compute the Selyaninov's Hydrothermal Coefficient (HTC) in the 10-day periods (each single record covered 10 days in each location). In this paper, the HTC was considered to be the sum of air temperatures higher than 0 • C, while in the original version [27] the HTC was considered to be the sum of air temperatures over a given minimum air temperature (e.g., 5 • C or 10 • C; Radomski [28]). The climatic water balance (CWB) was collected from the Agricultural Drought Monitoring System for Poland (ADMS) (http://www.susza.iung.pulawy.pl), provided by the Institute of Soil Science and Plant Cultivation-State Research Institute (IUNG-PIB). If the value of CWB in the district of the trial location was lower than the critical value of CWB determined for a particular soil texture grouping (agronomic category, Jadczyszyn et al. [29]), the occurrence of drought was recognized. Drought length was determined on the base of the successive ADMS reports with drought occurrence (Tables 1 and 2).
The dates of winter wheat principal growth stages were provided by the COBORU experimental station and rated according to Meier [30], considering only the principal growth stages (e.g., 5), according to the BBCH (modified Zadoks) scale [31]. Due to the type of data and very similar results, the Selyaninov hydrothermal coefficient (HTC) was identical to the Aridity Index (AI) [34,35], which was not included in the paper. In this paper, the sum of air temperatures >0 • C was considered, differently than in Selyaninov's original HTC formula [27].

Statistical Analyses
In this paper, the assessment of cultivar adaptation to environments of different productivity was performed in the following steps: (1) calculation of expected wheat productivity, depending on environmental factors; (2) calculation of cultivar productivity in location that were treated as separate environments; and (3) recommendation of cultivars of a specific type and range of adaptation.

Calculation of Expected Productivity
First, the analysis of covariance (ANCOVA) was performed by considering the genotype and management intensity levels as fixed factors and the soil and weather parameters as covariates as follows: where: m is the overall mean, a i is the ith cultivar main effect, b j is the jth management intensity level main effect, ab ij is the cultivar-management interaction, HTC is the HTC index of nth period denoted in kth location and c n is the coefficient for this quantitative variable, DL, LS and pH are the drought length, arable land suitability and soil reaction, respectively, and e ijk is the error term. The location ( Table 1) was not treated directly as a factor but its influence was included by the location's quantitative parameters (Table 2 and [13]). In this way, the estimated yields were calculated based on all available data such as soil conditions and weather. Thus, the obtained yield values should contain a smaller portion of random noise than the observed data, the influence of the weather conditions being reduced.
Then, the stepwise procedure based on the Akaike information criterion was used for covariate selection [36]. The use of too many variables leads to overestimating the reduction of variables needed to increase the model prediction ability. Finally, selected covariates were used to fit a linear model to estimate the expected yield. In this way, the location estimated mean yield across cultivars (environmental mean productivity in kth location and jth management intensity; expected grain yield, Y Ejk ) can be obtained from the model in Equation (2), which is a transformation of model in Equation (1), by removing the terms associated with the cultivars (a i , ab ij ), the error term (e ijk ) and the unjustified (according to the Akaike information criterion) quantitative variables. Y Ejk = m + b j + c 1 × HTC 1,k + . . . + c 10 × HTC 10,k + c 11 × D Lk + c 12 × LS k + c 13 × pH k (2)

Calculation of Relative Productivity of Cultivars (Genotype Specific Reaction)
The simple linear regression model was used to calculate single cultivar yield depending on the expected mean wheat yield (independent variable), according to general equation: where: Yc ijk is cultivar yield for the kth location under jth management level for ith cultivar, g i is the intercept for ith cultivar, h i is the coefficient describing the relationship between the ith cultivar and the environmental expected mean productivity, similar to Finlay and Wilkinson [17], Y Ejk is the environmental mean productivity in kth location under jth management level and e ijk is the error term. This regression model in Equation (3) (one for each cultivar) describes the cultivar response to generally considered environmental conditions, especially to the water stress described by the ten HTC indexes. The cultivar (genotype) specific relative reaction (GSR) was calculated as a difference between estimated cultivar yield Equation (3) and estimated environmental productivity Equation (2). Due to the fact that average coefficients across cultivars average across cultivars (i index) for h i in Equation (3) equals one, the GSR can be simplified to: where g i and h i are as defined in Equation (3) and Y Ejk is as defined in Equation (2). This parameter, defined in Equation (4) was used for cultivar recommendation.

Recommendation of Cultivars
The GSR and the mean yield of wheat expected in a given environment (location) were the unique criteria used for selection of cultivar to be recommended. The relative specific reaction (GSR) retains the order of estimated cultivar yield in a particular environment (location). For this reason, the GSR values and wheat productivity levels (between 4 and 10 t/ha) was used for cultivar recommendation. All cultivars were ranked according to the GSR value at a given yield level, i.e., low (4 t/ha, corresponding approximately to national average, GUS [2]), medium (7 t/ha), and high value (10 t/ha). The five cultivars with highest GSR at each level were chosen for recommendation.
The selected top five cultivars were also compared with the official list of recommended cultivars by the COBORU for cultivation within the particular province (voivodeship) in 2019 [37].
The statistical analyses were performed using the R software [38]. The analysis of covariance and regression analyses were done with the use of the lm function and the orthogonal contrasts were chosen for factors. The selection of independent variables contained in the model was done according stepwise selection method based on Akaike information criterion (AIC) and the calculation was performed with the 'stepAIC' function from 'MASS' package [39]. Table 3 presents the estimates for the main parameters in Equation (1). The overall mean of winter wheat yield across experimental locations (m) was estimated at 7.70 t/ha. Cultivar main effects (a i ) show yield relative to the overall mean [17]. Positive effects of a i inform about a higher yield in relation to the overall mean and vice versa. The best performing cultivars, outyielding the overall mean by 0.2 t/ha or more, were RGT Bilanz, Artist, Bonanza and Frisky. The cultivars with the lowest yield were Ostroga, KWS Ozon, Owacja, Arkadia and Pokusa. Table 3. The effects of cultivar (a i ), management intensity level (b j ) and their interaction (ab ij ), the AMMI (additive main-effects and multiplicative interaction) stability measure and the Finlay-Wilkinson regression coefficient on winter wheat yield. The main effects of the management intensity level (b j ) indicates a mean difference between moderate (MIM) and high (HIM) management intensity levels of about 1 t/ha ( Table 3). The effects of the cultivar-management interaction (ab ij ) indicate the response of cultivars to the level of management intensity level (MIM or HIM). The effect of this interaction (Table 3) was sometimes: (i) around zero (Pokusa and RGT Kicker cultivars had estimated differences between management levels around 0.006 t/ha); (ii) positive for HIM (and negative for MIM at the same time), which indicated positive effect of improved management (0.186 t/ha KWS Ozon); or (iii) negative for HIM, which showed lower positive effect (than occurred for most others cultivars) of management level on yield for particular cultivar (Owacja, HIM main effect equal to 0.504 t yield/ha was reduced by 0.154 for this cultivar). Thus, the results showed which cultivars were more sensitive to the management level, and for those cultivars, the higher management intensity is especially important.

Cultivar
Stability indexes such as the Finlay-Wilkinson regression (FW) coefficient and the AMMI stability value (ASV) were considered for comparison purposes (Table 3). For the FW coefficient, the most stable cultivars are those with regression coefficient equal to one, being a possible measure for instability as the absolute distance from the unit coefficient [17]. As for the ASV, since it is computed as the Euclidian distance between the scaled genotype scores and the origin of the biplot in a two-dimensional space formed by the first two interaction principal components, a lower absolute value of this index indicates greater cultivar stability [40,41].
Based on a preliminary analysis of Table 3, by considering the FW stability index, the top five most stable cultivars are Bonanza, Hondia, LG Jutta, Patras and Rotax. When considering the AMMI stability value, the most stable cultivars are Formancja, Hondia, KWS Ozon, Linus and Rivero.

The Impact of Environmental Variables (Covariates)
The data set that was used in the analysis was characterized by a large interdependence between variables. Although HTC in 10-day periods was calculated based on precipitation and air temperature, it was most frequently correlated with rainfall [13].
Both the sum of squares in ANCOVA connected to the independent quantitative variables and the signs of the coefficients were different for each sub-model ( Table 4). The sum of squares and AIC result in the same ranking of factors influencing winter wheat yield, meaning that both criteria chose the same factors as main source of yield gaps. The signs of coefficients determine if the influence of the quantitative variable was positive or negative on the yield. The strongest impact on the winter wheat yield was obtained by HTC in the first 10-day period in June and in the third 10-day period in May, both with positive signal ( Table 4). The drought stress generated yield gaps. The next factors were drought length (DL) with negative impact, i.e., longer water deficit decreased yield, and soil quality (LS) with positive impact. The influence of HTC was positive or negative, independent of the period, so the water deficit is beneficial in some growing stages (Table 4 and Figure 1).
The 10 points of LS value increased yield by 0.93 t/ha. Greater soil pH results in lower yield, and an increase of 1 pH represents a loss of about 0.9 t/ha. Also, increasing the DL parameter by one (water stress time is longer by 10 days) results in a loss yield equal to 1 t/ha.
The response of winter wheat to the HTC index depended on the time period. A negative relationship between yield and HTC was observed in the second 10-day period in April, the first 10-day period in May and between the second 10-day period in June and the second 10-day period in July. This relationship was particularly strong in the first 10-day period in May and first 10-day period in July. In contrast, a positive relationship between yield and HTC was found in the third 10-day period in April and between the second 10-day period in May and the first 10-day period in June, being especially strong in the third 10-day period in May ( Figure 1 and Table 4). It is important to mention that CWB also had a significant effect on winter wheat yield, however, weaker in comparison with the HTC impact. For this reason, CWB was not included in the final model. The 10 points of LS value increased yield by 0.93 t/ha. Greater soil pH results in lower yield, and an increase of 1 pH represents a loss of about 0.9 t/ha. Also, increasing the DL parameter by one (water stress time is longer by 10 days) results in a loss yield equal to 1 t/ha.
The response of winter wheat to the HTC index depended on the time period. A negative relationship between yield and HTC was observed in the second 10-day period in April, the first 10day period in May and between the second 10-day period in June and the second 10-day period in July. This relationship was particularly strong in the first 10-day period in May and first 10-day period in July. In contrast, a positive relationship between yield and HTC was found in the third 10day period in April and between the second 10-day period in May and the first 10-day period in June, being especially strong in the third 10-day period in May ( Figure 1 and Table 4). It is important to mention that CWB also had a significant effect on winter wheat yield, however, weaker in comparison with the HTC impact. For this reason, CWB was not included in the final model. Based on the results described above, the considered independent quantitative variables (covariates) in the linear model ANCOVA were: DL, LS, pH and the set of HTC values in decades without precipitation and air temperature (both were included in HTC index). A large number of covariates was avoided to improve parsimony and avoid overfitting during estimation [42][43][44]. Akaike information criteria confirmed the importance of all covariates in the analysis of covariance ( Table 4).
The model (2) describes the mean yield under each environmental condition (YEjk) based on the studied weather (mainly HTC index) and soil parameters, and allows the calculation of the location means according to the analysis of covariance (Figure 2). The location Tomaszów Bolesławiecki was characterized by the lowest yield, whereas Czeslawice was the most productive location. The estimated productivity decreased from Czeslawice to Tomaszów Bolesławiecki and seems to be a continuous variable. The arable lands in Poland are not clearly separated in terms of environmental conditions so the continuous increasing in the adjusted mean yield (Figure 2) instead of a discrete behavior was expected. Based on the results described above, the considered independent quantitative variables (covariates) in the linear model ANCOVA were: DL, LS, pH and the set of HTC values in decades without precipitation and air temperature (both were included in HTC index). A large number of covariates was avoided to improve parsimony and avoid overfitting during estimation [42][43][44]. Akaike information criteria confirmed the importance of all covariates in the analysis of covariance ( Table 4).
The model (2) describes the mean yield under each environmental condition (Y Ejk ) based on the studied weather (mainly HTC index) and soil parameters, and allows the calculation of the location means according to the analysis of covariance (Figure 2). The location Tomaszów Bolesławiecki was characterized by the lowest yield, whereas Czeslawice was the most productive location. The estimated productivity decreased from Czeslawice to Tomaszów Bolesławiecki and seems to be a continuous variable. The arable lands in Poland are not clearly separated in terms of environmental conditions so the continuous increasing in the adjusted mean yield (Figure 2) instead of a discrete behavior was expected.  Figure 3 shows the genotype specific reaction (GSR) of each cultivar in environments of determined productivity ranging from 4 to 10 t/ha, based on the simple linear regression according to equation [4], for each single cultivar. The GSR gives the relative difference between expected yield of each cultivar and mean yield of all cultivars. The row position in the figure indicates relative productivity of a particular cultivar in comparison to the other cultivars. The intensity of shading informs the type and degree of winter wheat cultivar adaptation to environment productivity. The green color indicates high and positive reaction of cultivar while red means low and negative reaction in the environment of the specified productivity. The response of a single cultivar to the productivity of the environment is seen in each row.   Figure 3 shows the genotype specific reaction (GSR) of each cultivar in environments of determined productivity ranging from 4 to 10 t/ha, based on the simple linear regression according to equation [4], for each single cultivar. The GSR gives the relative difference between expected yield of each cultivar and mean yield of all cultivars. The row position in the figure indicates relative productivity of a particular cultivar in comparison to the other cultivars. The intensity of shading informs the type and degree of winter wheat cultivar adaptation to environment productivity. The green color indicates high and positive reaction of cultivar while red means low and negative reaction in the environment of the specified productivity. The response of a single cultivar to the productivity of the environment is seen in each row.

Recommendation of Cultivars
Cultivars of wide adaptation are marked as green lines in the full range of productivity. For example, Artist and Bonanza produce higher yield than the environmental average in both low and high productivity environments (GSR values from 0.20 to 0.53 t/ha). Cultivars of narrow adaptation are characterized by green color in a narrower range of productivity of environments. In environments with high productivity (from 7 to 10 t/ha) the cultivars RGT Bilanz and Frisky are recommended as producing higher yield (GSR values from 0.15 to 0.64 t/ha) than average wheat yield in these environments. Cultivars LG Jutta, Rotax and Belissa are characterized by good adaptation (GSR values from 0.17 to 0.53 t/ha) to environments with low and medium productivity (from 4 to 7 t/ha), similarly to Artist and Bonanza with wider adaptation to productivity of environments. Thus, our method recommends the following, top-five highest yielding cultivars: Belissa, Bonanza, LG Jutta, Artist, Rotax and RGT Bilanz in environments of low productivity (4 t/ha), three of them being the same as those suggested by the FW stability index. The cultivars Artist, Bonanza, LG Jutta and Rotax (two also recommended by the FW stability index) are recommended for environments of medium productivity (7 t/ha), and the cultivars RGT Bilanz, Frisky, Dolores, Hybery and KWS Kiran are recommended for locations of high productivity (10 t/ha).
Agronomy 2020, 10, x FOR PEER REVIEW 12 of 20 4, 5, 6 denotes environmental average yield in t/ha as it was calculated by Equation (2), the figure is sorted in descending order for yield at the highest management intensity, or most productive environments, this is at 10 t per ha.
Cultivars of wide adaptation are marked as green lines in the full range of productivity. For example, Artist and Bonanza produce higher yield than the environmental average in both low and high productivity environments (GSR values from 0.20 to 0.53 t/ha). Cultivars of narrow adaptation are characterized by green color in a narrower range of productivity of environments. In environments with high productivity (from 7 to 10 t/ha) the cultivars RGT Bilanz and Frisky are recommended as producing higher yield (GSR values from 0.15 to 0.64 t/ha) than average wheat yield in these environments. Cultivars LG Jutta, Rotax and Belissa are characterized by good adaptation (GSR values from 0.17 to 0.53 t/ha) to environments with low and medium productivity (from 4 to 7 t/ha), similarly to Artist and Bonanza with wider adaptation to productivity of environments. Thus, our method recommends the following, top-five highest yielding cultivars: Belissa, Bonanza, LG Jutta, Artist, Rotax and RGT Bilanz in environments of low productivity (4 t/ha), three of them being the same as those suggested by the FW stability index. The cultivars Artist, Bonanza, LG Jutta and Rotax (two also recommended by the FW stability index) are recommended for environments of medium productivity (7 t/ha), and the cultivars RGT Bilanz, Frisky, Dolores, Hybery and KWS Kiran are recommended for locations of high productivity (10 t/ha).

Main Effects of Cultivars on Wheat Productivity
The model used in the study takes into account fixed factors depending on farmers (cultivar, management intensity level and their interaction, Table 3) and covariates (justified weather and soil factors, Table 4). Among the fixed factors, management had a greater impact than cultivars and interaction between both factors, as it was observed before by Mądry et al. [45] and Rozbicki et al. [46].
Mean difference between moderate (MIM) and high (HIM) management intensity levels was of about 1.0 t/ha (Table 3) and it was similar to the respective difference reported by Rozbicki et al. [46] (0.9 t/ha).
The cultivars incurred different effects on winter wheat productivity (a i varying between-0.472 for Ostroga and 0.360 for RGT Bilanz, Table 3). A significant effect of winter wheat cultivar on yield was reported by Mądry et al. [45] in contrast to Rozbicki et al. [46]. However, both studies differ partially in terms of the growing season, selected cultivars and location.
The cultivar-management interaction gives the response of the cultivars to the level of management intensity (MIM or HIM). The highest and positive ab ij values at HIM (and, consequently the lowest negative ab ij values at MIM) distinguish cultivars with the highest requirements regarding agrotechnology level. Examples of such cultivars were KWS Ozon, KWS Firebird, Arkadia and Bonanza. The cultivars Owacja, Hybery, Delawar and Formacja, characterized by lower requirements regarding agrotechnology, actually responded negatively on improvement of management intensity (Table 3). This fact is quite surprising and requires explanation by further and more detailed studies. However, according to Praczyk [47], the growth regulators used for stem shortening, which were applied at HIM, and not at MIM level, may have negative effect on crops during the years with longer periods of dry and hot weather, as it actually occurred in 2018.
The effects of cultivars and cultivar by management interaction were not quantified in previous studies but it may be useful for recommendation of cultivars depending on technology level.

Soil
In the model, land suitability was the fourth most important quantitative (environmental) factor determining winter wheat yield, after HTCs in the first 10-day period in June and in the third 10-day period in May and drought length in the dry year 2018. The importance of soil as the main or one of the most important factors affecting crop yield is widely known and has been continuously confirmed in previous and recent studies [10][11][12][13]48]. However, according to Iwańska and Stępień [13], based on the same yield, soil and weather dataset, the arable land suitability had stronger effect on winter wheat yield than any weather factor. At the same time, while in this study the HTC values for the third 10-day period in May, and first 10-day period in June as well as drought length had stronger and more significant effect on wheat yield than arable land suitability (LS). This may result from the use of a simple Spearman's rank correlation coefficient by Iwańska and Stępień [13], while in this study more advanced, multivariate statistical analyses were performed.
In the model developed in this study, the 10 points assigned to the LS group corresponded to 0.93 t/ha of yield increase (Table 3). Actually, the number of points assigned to each land quality group depended on yield of the main cereals obtained in 1970s [24]. In those times the maximum average cereal yield obtained on the best land suitability groups (1-94 points and 2-80 points) were of about 4.8-5 t/ha. Currently, the progress in technology used in agriculture led to an increase yield in winter wheat to about 8-10 t/ha and even more (Table 1) in these suitability groups. Thus, the valuation of land suitability in points developed in 1970 s [24] is still valid and well related to yield.
The increase of soil pH by one reduced average winter wheat yield by 0.9 t/ha. In this study the soil pH values varied between 4.9 and 7.0, and, most frequently, it was above six. Fotyma and Zięba [49], Farhoodi and Coventry [50] and Miller [51] reported a possibility of wheat yield decrease at such pH values. Also, a Spearman rank correlation analysis performed by Iwańska and Stępień [12] indicated a significant and negative effect of soil pH on wheat yield in 2015 and lack of any effect in the years 2016-2018.

Weather
The impact of some weather factors on winter wheat yield was greater than the impact of land suitability (Table 4). These factors were the HTC in the last ten days of May and the first 10 days of June and the drought length. The 10 days of DL decreased yield by one t/ha and confirmed the occurrence of drought, as reported by the agricultural drought monitoring system for Poland (ADMS) [52].
The analysis of the same data set using the Spearman correlation rank coefficient showed a negative relationship of winter wheat yield with drought length [13]. The impact of HTC on winter wheat yield were positive or negative depending on period similarly to Babushkina et al. [9] and Iwańska and Stępień [12]. As the temperature and precipitations vary across years, so does the HTC values. Consequently, the effect of HTC, as well as T and P, on wheat yield, is year-specific [12].
A negative coefficient for HTC indicates a decrease in yield along with an increase in HTC and as the amount of water increases, which might indicate the occurrence of excess of water in this period. A negative effect of HTC on wheat yield was found in the second 10-day period in April, which corresponded to tillering (BBCH 2), with particularly strong HTC in the first 10-day period in May, which corresponded to shooting (BBCH 3) and between the second 10-day period in June and the second 10-day period in July, which corresponded to milk maturity (BBCH 8) and ripening (BBCH 9) with particularly strong HTC in the first 10-day period in July ( Figure 1 and Table 4). During this period there were two highest values of HTC significantly reducing wheat yield by an amount ranging from 0.0 t/ha (Czesławice) to 5.0 t/ha (Pawłowice) and from 0.0 t/ha (Nowa Wieś Ujska) to 2.7 t/ha (Pawłowice). Similar results as for the final winter wheat growth stage under the strong influence of water availability in June and July were reported by Wójcik-Gront [11].
In contrast, a positive coefficient in the ANCOVA analysis for HTC indicates a decrease of yield, along with a low amount of HTC, indicating water deficit and hence drought ( Table 4). The lower the value of the HTC index, the greater the decrease of yield in the respective period. A positive impact of HTC was noted in the third 10-day period in April, which corresponded to tillering (BBCH 2) and between the second 10-day period in May and the first 10-day period in June, which was related mostly to heading (BBCH 5) and flowering (BBCH 6), with particularly strong HTC in the third 10-day period in May ( Figure 1 and Table 4). This period is critical regarding water supply for winter wheat, as it reported previously by Hanson and Nelsen [6], Rane et al. [53], Podolska [8] and Senapati et al. [54]. During the third decade of May the lowest values of HTC significantly reduced wheat yield by an amount ranging from 0.0 (Tomaszów Bolesławiecki, the most pronounced water deficit of all locations, HTC 0.0, data not shown) to 4.4 (Zybiszów, the best water supply of all locations, HTC 1.6, data not shown). In this study, the effect of the climatic water balance on winter wheat yield was significant (data not shown). The importance of CWB for winter wheat yield was reported by Wójcik-Gront [11] using the data from seven growing seasons (from 2009-2010 to 2015-2016). However, CWB did not show drought in Tomaszów Bolesławiecki, while the analysis of HTC in decades showed a strong negative impact of drought on winter wheat yield in the third 10-day period in May ( Figure 1, Table 4).
Moreover, HTC proved to be more useful than CWB, which is determined in 60-day periods and difficult to relate to grow stages. Probably, CWB reported for shorter periods, would lead to similar conclusions as HTC [12]. DL was calculated using CWB and was included in the model. Consequently, HTC and DL are also more convenient to be used in yield modeling and plant breeding programs.

Recommendation of Cultivars and Validation
In the current study, cultivar recommendation was based on genotype specific reaction methodology taking into consideration the criteria of winter wheat yield for the selected cultivars, which were those occurring in each of 19 locations during the 2017-2018 growing season. The top-five highest yielding cultivars were selected for good adaptation to specified productivity levels of the environment such as low (4 t/ha), medium (7 t/ha) and high (10 t/ha) (see Section 3.2) and also because of their presence in the list of recommended cultivars by COBORU for cultivation within the province in 2019.
These top five cultivars were then compared with the list of recommended cultivars by the COBORU (Research Centre for Cultivar Testing) for cultivation within the provinces (voivodeships, or administrative areas) of Poland in 2019 ( Table 5). Five of the top-yielding cultivars for the low value of wheat productivity (4 t/ha) were found on the list of recommended cultivars by COBORU [37]: Artist (in 16 provinces), Rotax (in eight provinces), Belissa (in seven provinces), Bonanza and LG Jutta (in five provinces). Table 5. The selected TOP five highest yielding cultivars for the specified value of Genotype Specific Reaction (GSR) and their presence in the list of recommended cultivars by the Research Center of Cultivar Testing (COBORU) for cultivation within the province in 2019.

Low (4 t/ha) Medium (7 t/ha) High (10 t/ha)
For the medium level of wheat yield (7 t/ha), Artist (in 16 provinces), Rotax (in eight provinces), Bonanza and LG Jutta (in five provinces) and RGT Bilanz (in three provinces) were also recommended by COBORU. The top-five provinces for the high yield level (10 t/ha) that were selected in this study were Hybery (in four provinces), RGT Bilanz (in three provinces), KWS Kiran (in three provinces), Frisky (in two provinces) and Dolores (in one province), and also recommended by COBORU.
Cultivars, which we recommend for lower productivity, are recommended by COBORU in more provinces, and those that we recommend for higher productivity, are recommended in fewer provinces (Table 5). It is also worth to note, that some cultivars recommended by COBORU in more than 10 provinces (RGT Kilimandżaro, Linus and Hondia) were not found within our top five cultivars selected for recommendation in low-, medium-and high-producing environments.
The differences in recommendations made by COBORU and those resulting from the current study may be explained as follows: (1) Our criterion is based on the expected yield, whereas COBORU does not regard yield as the only criterion; (2) Our recommendations are based on a one-year study, since the main goal is to develop a methodology, and COBORU considers the results from several years; (3) At COBORU experimental sites, it is possible to provide generally higher agrotechnology even at MIM level, as the experiments are carried-out on small plots, than in production conditions [20,22]; (4) At COBORU experimental sites, wheat is usually grown (on average) on better kinds of soils than (on average) under production conditions (see Table 1); (5) The detailed algorithm used by COBORU for cultivar recommendation is not available to the public; (6) Not all cultivars recommended by COBORU were included in this study, as they were not investigated in all locations in 2018.
The ranking for the top five most stable genotypes was computed for each stability measure, including our proposal, the FW stability index and the AMMI stability value, and Linus only for the AMMI stability value. Two genotypes, Formacja and KWS Ozon, were in the top five for our method and for the AMMI stability value. Genotype Hondia was on the top five for both FW stability index and AMMI stability measure.
The Spearman correlation of rankings for these three stability measures (Table 6) showed that the AMMI stability value was significantly correlated with the other two, and compliance with them was comparable. However, no significant correlation was found between our method and the FW stability index. Since only the proposed method takes into consideration weather conditions and soil yield-forming potential, it is expected to result in different ranking than other methods that do not take into consideration this information. The main limitation of our study is that it was based on single year. However, the main objective of this paper was to propose a new methodology for cultivar recommendation by considering the information of connected with environmental conditions and drought stress. The reasons for this approach to be based on one single year are: (i) often, the data available to researchers includes cultivars that are not repeated in consecutive years (this is the case for the COBORU research testing for variety testing); (ii) by including more than one year, many missing values would be included in the data and further strategies would need to be considered to generalize our proposal; and (iii) we were interested in proposing a relevant methodology and conclusions/suggestions based on recent data that can be useful for practical use in Polish plant breeding, instead of considering older data and a limited number of cultivars in a data set without missing values. However, we do believe that this article introduces the basis for a long-term analysis taking into account weather conditions, and, therefore represents a valuable step in this direction. Nevertheless, as a future development, we intend to develop a new methodology to analyze data collected along the time that can also account for missing values.

Conclusions
In this paper, we studied the relative performance of winter wheat cultivars in environments with different productivity and proposed a method for cultivar recommendation that takes into consideration the information of environmental conditions and drought stress. This was performed in three steps: (1) calculation of expected wheat productivity, depending on environmental factors, (2) calculation of relative productivity of cultivars in the environments, and (3) recommendation of cultivars of a specific type and range of adaptation.
This study confirmed the importance of soil and weather conditions as important factors to determine wheat yield. However, the weather factors should be considered in shorter time periods such as 10-day periods and in strict relationship with crop growth stages. Even during dry years, the crop may suffer not only from water deficit, especially in the period during or near to heading, but also water excess, especially in early spring, i.e., during tillering, and during ripening.
From the weather variables affecting crop yield, the Selyaninov's hydrothermic coefficient (HTC), which is calculated from precipitation and temperature, was more useful in wheat yield modeling than precipitation and temperature treated as separate variables. HTC was also more strictly correlated with wheat yield than the climatic water balance in 60-day periods, which is freely available in Poland, without any charges.
In the current study, we propose a genotype-specific reaction methodology for the recommendation of a cultivar for environments of a determined wheat productivity level. The method should be further developed for the analysis of long-term and to account for missing values, in order to increase its precision. Cultivars of wide adaptation should be recommended especially for fields with high soil variability, i.e., with large variations in productivity within the field. Cultivars of narrow adaptation should be recommended for fields with low soil variability.
This approach can be especially useful for farmers having long term experience, whose knowledge allows to assess the wheat yield level expected on their farms or even in particular places of one field, and then select the best-adapted cultivar for the respective wheat productivity level and also for their specific agrotechnology level.