3.1. AMMI analysis
3.1.1. Analysis of variance for the yield
Simple residual diagnostic plots revealed that grain yield values were normally distributed and had constant variance (Fig. 2). ANOVA was used to estimate the main effect and the interactions between and within the source of variations (Table 3 and Supplementary file 3). The results showed that the environmental effects accounted for 86.06% of the total variance, while the genotype and GEI explained 2.22% and 11.73% of the total variation, respectively (Table 3). The GEI effect was approximately five times greater than the G effect, showing a significant effect of GEI on the mean yield performance. This observation could be supported by the GEI effects on the mean yield, which ranged between 0.21 tons ha-1 for genotype G16 in E1 and 10.13 tons ha-1 for G15 in E5, whereas the genotype effects varied from 3.66 tons ha-1 for G24 to 5.10 tons ha-1 for genotype G18 (Table 4). A high level of variability was also detected for the mean grain yield across the environments, varying from 0.96 tons ha-1 in E1 to 9.37 tons ha-1 in E5. Notably, the environment (E1) was affected by severe drought stress due to the low level of precipitation and lack of supplemental irrigation applied during the growing season (Table 1), which explains its low yield performance. The remarkable variation defended by location suggested that the environments studied represent different agroclimatic conditions, which caused the most variation. The effects of the GEI magnitude on the mean grain yield were also observed in the box plot (Fig. 2).
Furthermore, the GEI consists of nine IPCA components whose weight is the sum of squares (SS) with decreasing significance. F-tests were used to estimate the significance of those components at the 0.01 probability level (Table 3). Splitting GEIs into multiplicative components revealed that the first four IPCAs were significant (p < 0.01) and increased the GE interaction variance by 75.32% (Table 3). The first three IPCAs, 1, 2, and 3, accounted for 33.81%, 17.91%, and 13.77% of the total variation, respectively (Table 3).
3.1.2. AMMI biplot
To identify the ideal genotypes for multi-environments, the AMMI biplot indicated that the genotypes could be discriminated based on their performance in different environments (Fig. 3A). Genotypes G10, G6, and G8 exhibited the highest values of IPCA1, whereas the nearest value of IPCA1 to zero was detected for genotypes G22 and G23. Likewise, the highest values of IPCA1 (Fig. 3A) and yield productivity (Table 4) were recorded for environments E2, E5, and E6, demonstrating their high impact on genotype discrimination. On the other hand, the closest value to zero of IPCA1 was detected for environments E1 and E4 (Fig. 3A), which had the lowest grain yield, with the minimum impact on GEI (Table 4). E3 and E7 were among the environments with the greatest contributions to the GEI, as indicated by their long vectors (Table 5). The environments that have acute angles between their vectors are most strongly associated with genotypes ranked in the same direction, while environments with an obtuse angle are least strongly associated with genotype rank, and those with right angles are independent and not associated. Positive associations were observed between environments E9 and E5 and between environments E5 and E6, and they were negatively correlated with E3 and E7. The genotypes G10, G6, and G8 interacted positively with E3 and E7 but negatively with the E2 environment (Fig 3A). According to our results, several genotypes, such as G13, G15, G18, and G21, exhibited general adaptation with high grain yield (Table 4), while most of them tended to exhibit specific adaptability.
3.1.3. Genotypes recommended by the AMMI model
According to the AMMI results, four promising genotypes for each environment were identified (Table 5). Genotype G10 was ranked as the major genotype in five environments: E3, E7, E8, E1, and E4. Genotype G15 was observed in three environments, with E6, E10, and E5 being the four most highly ranked and dominant genotypes (Table 4). Genotype G16 was observed in environment E9, and E2 was the genotype with the highest rank among the four genotypes in the three environments. G13 was found to be the dominant genotype in the E6, E10, and E5 environments and ranked among the top four genotypes in the three environments. Finally, G18 was dominant in the E1, E10, and E9 environments and was identified as one of the top four genotypes in the four environments. These results revealed that the selection of these genotypes for related environments by using the AMMI model is a suggestion of the best-adapted genotypes in those specific environments.
3.2. GGE biplot analysis
3.2.1. The winning genotypes for the studied environments
A highly significant result was observed for the AMMI’s residual mean square (Table 3), indicating the significance of creating a GGE biplot to visualize ‘which-won-where’ patterns of environments and genotypes and the environment's ability to discriminative and representativeness. Based on the ‘which-wins-where’ pattern, the environments and the genotypes were divided into five and six groups, respectively (Fig. 3B). A polygon based on the GGE biplot was created by associating the vertex genotypes with straight lines, while the rest of the genotypes were assigned inside the polygon. Accordingly, the vertex genotypes were G15, G10, G6, G5, G25, G16, and G19, which fall inside the relevant sector of the GGE biplot polygon because they are the best in the environment. Therefore, the performance of these genotypes is either the poorest or the best in one or more environments, as they have the greatest distance from the origin of the biplot. Genotypes inside the polygon and nearer to the origin of the axes are widely adapted and have a minimal response to environmental differences [44, 45]. The responsive genotype has either the poorest or the best performance in a single or all environments [46]. The G13, G15, and G18 genotypes exhibited high yields with overall adaptation, whereas G24 was considered the lowest-yielding genotype among the vertex genotypes. Furthermore, the environment fell inside the sectors of the vertex genotypes G13 and G15, which indicated that these genotypes were among the best genotypes in all the environments studied. The estimation of mega-environments with the winning genotypes is also an additional advantage of the GGE biplot. The results of the present study indicated that the meta-data could be divided into five mega-environments (ME1, ME2, ME3, ME4 and ME5) (Fig. 1 and Fig. 3B). Including the testing environments, E5, E6, and E10 fell into ME1; E3 and E7 fell inside ME2; E1, E4, and E8 fell inside ME3; and E9 fell inside ME4, while E2 fell inside ME5. The vertex genotypes in every sector are the best in the environments where their signs fall into the relevant sector. The environments inside the same sector contain similar winning genotypes and vice versa. Accordingly, G13, G15, and G18 were the winners and high-yielding performance genotypes in the ME1, ME2 and ME3 mega-environments, respectively. As a previous result, the GGE biplot polygon view is among the best methods for identifying the winning genotypes by displaying the patterns of interaction between genotypes and environments [28].
3.2.2. Selection for yield and stability performance
The stability vs. mean yield performance of the GGE biplot (Fig. 4A) displayed the ranking of genotypes beside the ‘average tester coordination’ (ATC) abscissa according to their mean performance. The results showed that G13 had the highest mean yield, followed by G18, G15, G14, G16, G21, G19, G10, G28, G32, G9, G2, G7, G12, G31, and G17, whose mean performance was greater than average (4.43 tons ha-1) (Table 3). In contrast, the G1, G4, G17, G13, G23, G26, and 27 genotypes had the lowest contribution to GI. Based on the results, the yield and stability of the G13, G15, G18, and G32 genotypes were the greatest.
3.2.3. Assessment of environments studied and genotype adaptation
The discrimination ability of the environment studied is relative to the length of the environment vector. Therefore, environments with prolonged vectors are considered to have a greater discriminating ability for genotypes because they provide information on genotype differences and vice versa [44]. In our study, environments E2, E5, and E6 exhibited the maximum discrimination power, which provided significant evidence of differences between genotypes, while environments E1, E8, and E4 had the lowest values (Fig. 4A). The representative environment depends on the size of the angle between the abscissa of the average environment coordinate (AEC) axis and the vector of an environment. The smaller the angle is, the stronger the representativeness of the test environment is [28]. The E9, E10, E6 and E5 environments exhibited reasonably strong representativeness (Fig. 4B). Accordingly, an ideal environment would have discrimination power (strength to define the studied genotypes) and representativeness (strength to represent the mega-environment). Among the environments studied in the present study, E5 was classified as an ideal environment that is powerful for discriminating genotypes in addition to being the most representative of all environments studied (Fig. 3A and Fig. 4B). However, discriminating ‘nonrepresentative’ environments such as E3 are suitable for selecting particular genotypes adapted when the target environments can be separated into mega-environments, or it could be valuable for discarding unstable genotypes in the case of a single mega-environment as a target environment. G13, G15, G18 and G21 exhibited the most adapted genotypes to these ideal environments. This can also be proven as shown in Fig. 4C, as it displays the genotype rankings relative to the ideal genotype. In that case, the ideal genotype is the one with the maximum yield and stable performance. Therefore, the G15, G13 and G32 genotypes could be the ideal genotypes among the genotypes G18 and G21 (Fig. 4C).
3.2.4. Assessment of genotypes according to stability indices
A total of 32 univariate and multivariate models were used to investigate the stability performance of the genotypes studied (Tables 6 and 7). According to the GY as a key parameter for genotype assessment, G13, G18 and G15 exhibited the greatest mean performance across environments; in contrast, the G24, G25 and G1 genotypes exhibited the lowest grain yield. These outcomes were confirmed graphically by the GGE biplot method (Figure 4B and Table 6), which demonstrated that the GGE biplot results were consistent with the yield meta-data.
According to parametric stability statistics (Table 6), the coefficient of variability (CV%) (Francis and Kannenberg, 1978) revealed that G17, G75 and G7 were considered the most stable and desirable genotypes with small CV% and high GY, while genotypes G25, G26 and G24 were unstable (Table 6). Based on the deviation means squares model (Eberhart & Russell, 1966), G27, G2 and G32 genotypes with bi ≥ 1 could be more adaptable to favorable conditions, while G15, G13 and G19 (bi < 1) could be adapted to unfavorable growing conditions. Genotypes G32, G1 and 27 exhibited the lowest variances in deviation from regression (S2di ≤ 0) and were the most stable across the environments, whereas G16, G25 and G6 (S2di ≥ 0) were unstable and exhibited lower stability. According to the coefficient of determination (R2) [16], G32, G13 and G18 had high R2 values and were the most stable genotypes, whereas genotypes G6 and G16 were unstable. Wricke’s ecovalence (W2i) and Shukla's stability variance (σi2) models showed that G32, G1 and G27 were stable and exhibited minimum stability variance across the environments, whereas G16, G25 and G6 were unstable genotypes. According to Perkins and Jinks's stability model, G27, G25, and G2 had the lowest regression coefficients (βi), while G32, G1, and G27 had minimum nonlinear sensitivity to environmental change (Dji). These genotypes are considered to be more stable across the environment. These outcomes were comparable to those obtained by deviation mean squares analysis. Using Lin and Binns's (1988) superiority index (Pi), genotypes G13, G18 and G21 expressed the lowest values of Pi and were identified as stable genotypes that exhibited high yield performance, whereas G5, G6, and G24 were unstable and low-yielding genotypes. According to Roemer’s environmental variance (Sxi2) model, genotypes G5, G17 and G4 exhibited minimum variance and were considered to be highly stable genotypes across the different environments, while G15, G13, and G24 exhibited unstable genotypes. According to the Tai stability statistics λ and α, genotypes G32, G1, and G27 were perfectly stable, as they displayed α values of −1 and λ close to unity (Fig. 5A). Based on this model, a stable genotype must have αi = -1 and λi = 1. However, a genotype with αi = 0 and λi =1 represents average stability, and those with α > 0 and λ = 1 indicate less average stability. Thus, genotypes G15, G13 and G19 exhibited average performance stability (Fig. 5A). According to Plaisted and Peterson’s low pairwise genotype-environment interaction (P59) and Hanson stability (Di2) models, G32, G27, and G1 and G32, G1, and G27 were the most stable genotypes and presented the lowest values of P59 and Di2, respectively, across the environments. According to Lewis’s stability factor (SF) (1988), the most stable genotypes were G5, G17 and G9, which exhibited high-yield performance and relatively little impact across environments. However, G16, G18, and G23 were the least stable genotypes, as the SF deviation was greater than unity.
Furthermore, for the nonparametric stability models (Table 7), the AMMI stability value (ASV) model indicated that G15, G16 and G26 were the best-ranked genotypes. These genotypes exhibited lower values of ASV and were considered to be more stable genotypes across the environments, while G8, G11 and G29 were unstable (Table 7). The yield stability index (YSI) identified G13, G2, G14, and G15 as desirable genotypes that combined high yield performance and stability (Figure 4C). According to Nassar and Huhn’s nonparametric stability model, G32, G1, and G27 were the most stable genotypes in the environment. These genotypes showed the lowest values for SI 1, with a yield greater than the overall mean (4.43 tons ha-1). According to the nonparametric stability model, G32, G1, G27, G23, G24 and G25 are stable genotypes. These genotypes displayed the lowest values for the four stable NPIs (1), NPIs (2), NPIs (3), and NPIs (4), which are based on adjusted means ranks of the genotypes in every environment. According to this model, the stable genotypes are those whose location is consistent in relationship to the others that remained unchanged in the set of environments evaluated based on adjusted means ranks of the genotypes in every environment. Ketata et al. (1989) proposed a ranking method (kr) by plotting the mean grain yield (gy) of the genotypes against the standard deviation (δr) across environments. A genotype is considered to be stable when its kr or gy value is high and relatively consistent with a low δr value across the environment. Accordingly, G32, G15, G7, and G9 in section 1 were the genotypes with the greatest mean and greatest stability performance. The grain yields of these genotypes were greater than the overall mean (4.43 tons ha-1) (Table 4), suggesting that they might be used for cultivation in favorable environments, while the genotype with the greatest mean performance and the most stable genotype in section 2 was G27. (Fig. 5B). However, the common adaptability of genotypes was detected in sections 1 and 2. The genotypes in sections 3 and 4 lacked mean yield and stable performance; however, their mean grain yield was moderate to high when single environments were considered. The last nonparametric model used in this study is the TOP-rank stability parameter of Fox et al. (1990). In this model, Fox et al. suggested the superiority of environments in which the genotype ranked in the top, middle, or bottom third of trial entries. Therefore, G7, G9, and G18 were the most stable genotypes, as they ranked in the top third of genotypes with a maximum percentage (60%) of environments. These genotypes had a high mean (4.43 tons ha-1) of 4.71 tons ha-1, 4.72 tons ha-1 and 4.86 tons ha-1, respectively (Table 4).
According to the two stability methods, G32, G27, and G1 were the most stable genotypes, with 17, 17, and 16 stability statistics, respectively, out of the 32 stability models used. Therefore, the genotypes G32, G27, and G1 are superior genotypes that might be utilized in barley breeding programs to improve barley production and stability across the environment. Furthermore, G13, followed by G14, G15, and G23, was the most stable genotype based on multivariate analysis only.
3.2.5. Selecting the best stability statistics
To identify the best stability parameters for the simultaneous selection of high yield and stable performance, a hierarchical cluster analysis [47] was performed. Based on the rank correlation matrix, the 32 stability measures were divided into four different clusters (Fig. 6A). The first cluster (CI) comprised the Pi and nonparametric statistics Kr and RSM. In the second cluster (C2), the following stability measures were measured: δr, NP I (1), NP I (2), SI 1, SI 2, SF, P59, σi2, W2i, Di2, λi, S2di and Dji. were classified together in the same subcluster. Cluster three (C3) consisted of nonparametric stability statistics of NP I (3), SI 3, SI 6, NP I (4), TOP and YSI along with the GY (Fig. 6A). Cluster four (C4) consisted of CV, bi, Bi, αi, Sxi2, δ gy, ASV, R2 and FW. Therefore, it could be safely suggested to use stability parameters from different groups to avoid the probability of assessing the same concept of stability. Furthermore, a pattern map was developed based on the squared Euclidean distance for the classification of the stability parameters as well as the 32 barley genotypes studied (Fig. 6B). The results indicated that nine distinct classes and subclasses of genotypes (C1 to C9) were present, revealing the genetic diversity among the genotypes studied. The first genotypic C1 comprised G32, G27, and G1, which were the most stable genotypes studied (Table 4). The second class (C2) included genotypes G23, G3, G20, G30, G29, and G22; the third class (C3) comprised genotype G16; the fourth class (C4) comprised genotypes G19, G14, and G13; G14 comprised 5.13 tons ha-1 and 5.06 tons ha-1; the fifth class (C5) comprised genotypes G2, G10, G12, G17 and G18; the sixth class (C6) comprised genotypes G31, G28, G21, G7 and G9; the C7 class comprised genotypes G15, G11, G8 and G24; the C8 class comprised genotypes G25, G6 and G26; and the C9 class comprised genotypes G4 and G5. In addition, six groups for stability parameters, I, II, III, IV, V and VI, were also identified and can be used for the characterization of genotypes based on stability groups (Fig. 6B).
The stability parameters R2, FW, and ASV are classified in the same group (I), showing that they are closely associated with genotype ranking (Fig. 6B). Group II comprised CV%, bi, Bi, αi, Sxi2 and δ gy, which were associated with a ranking of genotypes in related directions. The III group consisted of NPI (3), SI 3, SI 6, NPI (4), TOP, YSI and GY in a similar group, displaying the close correlation of these stability statistics with GY. Group IV consists of SF, P59, σi2, W2i, Di2, λi, S2di and Dji. The V group consisted of the nonparametric measures δr, NPI (1), NPI (2), SI 1 and SI. The VI group consisted of Kr, Pi, and RSM. The results indicate that most of the nonparametric statistics were identified in groups IV and V, indicating that using different statistics could be effective for the selection of stable genotypes for grain yield in barley.
3.2.6. Associations between the stability parameters
To understand the associations between parametric and nonparametric models with the concepts of stability, a graphical evaluation using Spearman's correlation matrix was performed (Fig. 7 and Supplementary file 4). The correlation matrix showed that the stability parameters were divided into five clusters (I - V), in which clusters II and III consisted of two subclusters (Fig. 7). The results indicated that grain yield was positively and significantly correlated (P < 0.01) with the bi, αi, Bi, and Sxi2 stability statistics (Supplementary file 4). These statistics are associated with parametric measures and are clustered into subgroup (d) of cluster III (Fig. 7). Furthermore, GY was positively and significantly correlated with the nonparametric measures YSI, NPI (4), TOP, NPI (3), SI 6, δ gy, and SI 3, which were clustered in subgroup (c) of cluster III (Figure 7). In contrast, GY was negatively correlated (P < 0.01) with the Pi, RSM and kr parameters. The coefficient of regression (bi) was significantly correlated with Sxi2, αi, and δ gy (P< 0.01), while the deviation from regression (S2di) exhibited a positive and significant correlation with σi2, Dji, W2i, λi, P59 and Di2 (P < 0.01). σi2 was positively and significantly correlated (P < 0.01) with Dji and W2i. W2i was positively and significantly correlated with λi, P59, and Di2 (P < 0.01). Pi was positively and significantly correlated with the RSM and kr (P < 0.05). Sxi2 was positively and significantly associated with αi and δ gy (P < 0.01). The λi was positively and significantly correlated with P59 and Di2 (P < 0.01). SI 1 and SI 2 were positively and significantly associated with NP I (1), NP I (2), and δr (P < 0.01).