Advances in Crop Science and Technology
Open Access

AMMI model; BLUP model; Genotypes; GGE biplot; GGE; stability analysis; Faba bean; R software

Introductions

Faba bean (Vicia faba L.) is an important cool-season grain legume in terms of its global area coverage and yearly production volume. It is being cultivated in more than 70 countries around the globe (FAOSTAT, 2021). Ethiopia is the second-largest faba bean (Vicia faba L.) growing country in the world next to China and the first producer in Africa, followed by Egypt, Sudan, and Moroko (FAOSTAT, 2021). According to the CSA (2021) report, the total area covered by Faba bean in the 2020 “Meher” production season was 504,569.99 ha, with a total annual production of 10,706,36.538 tons. The average national productivity of the crop was 2.122 tons per hectare. The crop is mainly cultivated in mid- and high-altitude areas of the country, with an elevation ranging from 1800 to 3000 meters above sea level. As a source of food and feed, the faba bean is important to the socioeconomic well-being of Ethiopia's small-holder farming communities. In addition, the crop provided the nation with foreign exchange and a reliable source of income for farmers. Combined with cereals like wheat and barley, it also serves as a good break crop for pests and benefits in restoring soil fertility [1].

The Faba bean crop is grown and adapted to different agroecological conditions and soil types in Ethiopia. Many scholars have stated that faba bean yield performance is unstable and affected by environmental variations (Cernay et al., 2015; Reckling et al., 2018). On the contrary, farmers need to use well-adapted and stable genotypes and good agronomic practices to boost productivity and get a high return from their lands. As the main goal of faba bean breeding programs is to increase the productivity of the crop by developing high-yielding, stable, and widely adapted varieties, breeders test large numbers of genotypes in various environments to evaluate the yield stability and wide adaptability of the genotypes. So, multi-environment trials are important in interpreting the genotype by environment interaction (GEI) effect and selecting superior genotypes at a later stage of variety developments. GEI arises due to the differences in the sensitivities of genotypes to different environmental conditions [2].

Therefore, to identify and select well-buffered and stable genotypes, the study of genotype by environment interaction (GEI) is very important for crop improvements. A better understanding of the level of GEI and performance stability in crops is a decision tool, mostly at the final stage of the variety development process. This helps the breeder generate essential information on the adaptation pattern in breeding lines and new varieties for release and determine the recommendation domains for released varieties (Yan W., 2011).

GEI can be computed using several procedures, all based on evaluating genotypes in multiple environments. The additive main effects and multiplicative interaction (AMMI) model is a multivariate parametric approach that is widely used to analyze and interpret GEI in METs (Gauch, 1992). The AMMI method combines variance analysis for additive or main effects and principal components (PCs) analysis for multiplicative effects to understand the patterns of GEI (Zobel et al., 1988). From a practical point of view, BLUP and AMMI can be seen as two distinct approaches to achieving the same goal. From a statistical point of view, these models are vastly different. The AMMI analysis retains most of the GEI pattern in the first interaction principal component axis (IPCA) resulting from the singular value decomposition (SVD) of the nonadditive effects matrix, while most of the random error is retained in the last IPCAs. The BLUP initially estimates the effects of the ANOVA model and then attributes weights to these effects; it could thus be considered a shrinkage estimator (Piepho 1994). This method estimates the mean yield of genotypes in mixed models with high efficiency. On the other hand, GGE biplot analysis is a beneficial graphical tool since it offers visual pictures and a clear summary of the main data and outcomes (Yan, W., 2015). The AMMI model is often used together with the GGE biplot graphical model to identify MEs as well as winning genotypes in each ME. The unique feature of the GGE biplot is that, based on the plots, it can be decided which genotype has the highest potential in which environment. These models are frequently used alone in the evaluation of METs. Some studies were successful in estimating genotypic values in MET using BLUP (Olivoto et al., 2017; Nardino et al., 2016), while others were successful in modeling GEI patterns using AMMI (Bocianowski et al., 2019; Veenstra et al., 2019). Combining the graphical tools of AMMI, GGE biplot, and the predictive accuracy of BLUP is very important in exploring GEI. Thus, this study aims to assess the genotype x environment interaction, apply the stability parameters, identify environments that are more suitable for faba bean growing, and identify varieties with a high and stable yield [3].

Materials and Methods

Experimental site description

The experiment was conducted during the main cropping season ("Meher") of 2018 and 2019 for two consecutive growing seasons at four locations, namely Kulumsa, Bekoji, Asasa, and Kofele experimental stations of the Kulumsa Agricultural Research Center in south-eastern Ethiopia. The testing sites are located in different districts of the Oromia Regional States of Ethiopia, which are characterized by mid- and high-altitude agroecology. The dominant crops grown within the experimental areas are wheat, barley, faba bean, and root crops like potatoes (Table 1) [4].

Table 1: Descriptions Geographical Descriptions of the Experimental Environments.

Experimental materials and design

Fourteen faba bean (Vacia faba L) genotypes obtained from Holeta Agricultural Research Center were evaluated under potential environments. The experiment was laid out in a randomized complete block design (RCBD) with four replications. Each genotype was planted on a plot size of 6.4 m2, with 4 rows at 40cm spacing between rows and 4 m length at all the testing sites. The central two rows were seed kg/ha and fertilizer followed uniformly for all plots throughout the experiment and location, as per the recommendation for faba bean (Table 2).

Table 2: Lists of genotypes used in the study.

Data collected

Data on grain yield, agronomic conditions, and disease reactions were collected from each plot. Days to flowering (DF) and days to physiological maturity (DM) were taken when each plot reached 50% of flower initiation and 90% of the pod attained physiological maturity, respectively. The days were calculated starting from the date of sowing. Plant height (cm) was taken at full maturity from five randomly selected plants in the central two rows, measured from the ground level to the top of the plant. The mean value is recorded as plant height per plot for analysis. Responses of genotypes to disease reactions like chocolate spot, ascochyta blight, and rust were recorded at late pod setting based on 1–9 scoring methods. Grain yield was measured on clean, dried seed, and plot yields were adjusted to 10% moisture level and converted to kilograms per hectare. Thousand seed weights (TSW) (gm) were counted and weighted. Data on the number of pods per plant and seeds per pod were also collected based on five plant bases and averaged for data analysis [5].

The analysis of variance of each location and combined data over location were performed using a mixed linear model to assess the differences among genotypes as per Gomez and Gomez (1984). R software version 4.4.0 with the packages “agricolae” and “metan” were used. Homogeneity of variance was tested and a combined analysis of variance was done using the Mixed Linear Model procedure to partition the total variation into components due to genotype (G), environment (E), and G × E interaction effects. The following individual and combined RCBD models were used for analysis.

where; Y_ij is the grain yield of the i^th genotype in the j^th environment, µ = the grand mean, 𝐺_𝑖 = the effect of the i^th genotype, 𝐸_𝑗 = the effect of the j^th location, 𝐺𝐸_𝑖𝑗 = the interaction of the i^th genotype with the j^th location, 𝛽(𝐸)_𝑗𝑘 = the effect of the k^th replication in the j^th location, and 𝜀_𝑖𝑗𝑘 = the error.

AMMI model analysis

The Additive main effects and multiplicative interaction (AMMI) model was performed for grain yield and thousand seed weight of 14 faba bean genotypes using performs ammi() function in metan packages of R software. Therefore, the estimate of the response variable for the i^th genotype in the j^th environment (𝑦_𝑖𝑗) using the AMMI model, is given as follows (Gauch, 1992).

where Y_ij = is the yield of the i^th genotype in the j^th environment; μ = is the grand mean; G_i and E_j are the genotype and environment deviations from the grand mean, respectively; λ_k = is the eigenvalue of the PCA analysis axis k; α_ik and 𝝲_jk = are the genotype and environment principal component scores for axis k; n is the number of principal components retained in the model, and e_ij is the error term.

AMMI Stability Value (ASV) which is the distance from the coordinate point to the origin in a two-dimensional of IPCA1 score against IPCA2 scores in the AMMI model was calculated using the formula developed by (Purchase et al., 2000). This weight is calculated for each genotype and environment according to the relative contribution of IPCA1 to IPCA2 to the interaction Sum of Squares as follows:

Where: IPCA1 = interaction principal component axis 1; IPCA2 = interaction principal component, axis 2.

GGE biplot model

The GGE Biplot model, as introduced by Yan et al. (2000), utilizes biplots, an effective tool for visualizing two-way data commonly conducted in MET data analysis. This model enables a simultaneous display of genotype main effects (G) and genotype × environment effects (GE) from a two-way data table (Yan et al., 2000). The first component of the GGE biplot, when closely associated with the genotype main effect (G), indicates the proportion of production solely attributed to the genotype, while the second component represents the proportion explained by genotype-environment interaction (GEI). Singular value decomposition (SVD) of the first two principal components was employed to fit the GGE biplot model [6].

where, Y_ij is the trait mean for genotype i in environment j, μ is the grand mean, β_j is the main effect of environment j, μ + β_j is the mean yield across all genotypes in environment j, λ_1, and λ₂ are the singular values (SV) for the first and second principal components (PC₁ and PC₂), respectively, ξ_i1 and ξ_i2 are eigenvectors of genotype i for PC1 and PC₂, respectively, η j₁ and η_j2 are eigenvectors of environment j for PC₁ and PC₂, respectively, ε_ij is the residual associated with genotype i in environment j. In GGE biplot analysis, scores of PC₁ were plotted against PC₂.

Best linear unbiased prediction (BLUP)

The best linear unbiased prediction (BLUP) model is used to assess the genetic merits of each genotype using restricted maximum likelihood (REML) for variance estimation components in R. It is used in linear mixed models to estimate random effects.   The Predictive accuracy of models can be assessed using cross-validation methods. These can be done by dividing the data into training and validation sets (Gauch et.al., 1988).

Results and Discussions

Analysis of variance

A combined analysis of variance (ANOVA) of 14 faba bean genotypes tested over four locations and two years (8 environments) showed that all traits studied are significantly influenced by both genotypes and environments at the 1% probability level. Specifically, traits such as days to 50% flowering, number of pods per plant, number of seeds per pod, thousand seed weight, grain yield, and severity of rust disease were notably affected by genotype by environment interaction (GEI) at a significance level of P≤0.01. This indicates that the faba bean genotypes exhibited inconsistent performance across different environments. The Considerable variations observed due to the environment in this study suggested that there were noticeable variations in the experimental environments. These results agree with previous findings of Tamne T. (2015), Teklay A. et al. (2015), and Mesfin T. et al. (2020) which reported that genotype, environments, and genotype by environmental interaction were significantly different for days to flowering, pods per plant and grain yield. Study conducted by Dereje A. et al. (2019) on 14 faba bean genotypes in the Kellem Wollega Zones of Western Oromia also found that genotype by environment interaction significantly affected the number of pods per plant and grain yield (Table 3) [7].

Table 3: Mean square of combined ANOVA for nine traits of 14 faba bean genotypes conducted at four locations for two consecutive years (8 environments).

Additive main effects and multiplicative interaction (AMMI) analysis

Highly significant variations due to environments, GxE interactions, and genotypes were observed by AMMI analysis. According to the AMMI analysis of variance, about 43.3% of the total sum of squares (SS) of variation for thousand seed weight accounted by environment (E) followed by 34.4% of genotype (G) and 7.5% of Genotype by environment interactions (GEI). Whereas, 19.1%, 21.5% and 15% of the total sum of squares (SS) of variation for grain yield was accounted by environment (E), genotype (G) and Genotype by environment interactions (GEI), respectively. The results indicate that environmental factors had the most significant impact on thousand seed weight, while genotype had a greater influence on grain yield. The variance due to GEI had the smallest impact on the phenotypic variance for both traits. Seven interaction principal component axes (IPCA) were identified through principal component analysis (PCA) that capture the variance in the interaction data between genotypes and environments. In this context, they represent how different interactions affect the traits of interest (TSW and GYLD). Of these seven, two IPCAs were significant at the 1% probability level for TSW. This means that for thousand seed weight, only two of the seven IPCA dimensions are statistically significant with a high level of confidence (1% probability level). This suggests that interactions affecting TSW are primarily captured by these two principal components, which explains a substantial portion of the variation in TSW due to interactions. The rest of the dimensions did not show significant contributions. Regarding GYLD five IPCAs were significant at the 5% probability level. This indicates that the interactions affecting GYLD are more complex and involve a larger number of principal components. These five dimensions together explain a significant portion of the variation in GYLD, showing that GYLD is influenced by a broader range of interactions compared to TSW. This complexity imply that multiple factors and their interactions influence grain yield more extensively than thousand seed weight (Table 4) [8].

Table 4: AMMI analysis table of thousand seed weight (TSW) and grain yield (GYLD) of 14 faba bean genotypes.

AMMI model analysis

The biplot analysis presented insights into how thousand seed weight (TSW) and grain yield (GYLD) vary. By examining the first two principal components (PC1 and PC2) from the AMMI model, the analysis revealed that PC1 and PC2 together accounted for 85.3% of the total variation in thousand seed weight (TSW). This indicates that the first two principal components together accounted for 85.3% of the variation in TSW. This substantial percentage implies that PC1 and PC2 capture most of the key interaction patterns and variability in TSW. Therefore, the biplot for TSW clearly illustrate the interactions between genotypes and environments, emphasizing the primary sources of variation. For GYLD PC1 and PC2 Explained 52.3% of the total variation. This is a lower percentage compared to TSW, suggesting that PC1 and PC2 capture only about half of the variation in GYLD. The biplot for GYLD will thus represent a less comprehensive view of the interaction patterns, with additional principal components potentially being necessary to fully explain the variation. According to AMMI 1 biplot analysis, Numan variety had the highest TSW followed by Dosha and Cool-0030. And also, Numan, variety yields the highest GYLD followed by Dosha, Cool-0018, Cool-0030, Cool-0024, Cool-0035, Cool-0031, and Cool-0034 genotypes (Figure 1,2,3 &4) [9].

GGE biplot analysis (which-won-where view of GGE Biplot)

Figure 1: AMMI1 biplot (A) and AMMI2 biplot (B) for thousand seed weight of 14 faba bean genotypes evaluated under eight environmental conditions.

Figure 2: A nominal yield of TSW describing the "which-won where" view for the 14 faba bean genotypes as a function of the environment scores of the first interaction principal component axis (IPCA1) (C) and Heat map showing the TSW variation of 14 faba bean genotypes across 8 environments (D).

Figure 3: AMMI1 biplot (E) and AMMI2 biplot (F) for grain yield of 14 faba bean genotypes evaluated under eight environments.

Figure 4: A nominal grain yield describing the "which-won where" view for the 14 faba bean genotypes as a function of the environment scores of the first interaction principal component axis (IPCA1) (G) and Line map showing the grain yield variation of 14 faba bean genotypes across 8 environments (H).

The polygon view of the “which-won-where” of the GGE biplot revealed the interaction patterns between genotypes and environments and highlighted the top-performing genotypes. In this GGE biplot, a polygon was formed by connecting the vertex genotypes which were positioned far from the origin with blue lines enclosing all other genotypes within this shape. For thousand seed weight (TSW), the vertex genotypes were Cool-0002, Cool-0011, Cool-0024, and Numan. For grain yield (GYLD), the vertex genotypes were Cool-0031, Numan, Dosha, Cool-0035, and Cool-0011. This indicates that these genotypes represent the extremes in performance (either highest or lowest) for each trait across specific environments. The GGE biplot analysis showed that PC1 and PC2 together explained 95.9% of the variation for TSW and 73.38% for GYLD. The Numan variety was identified as the top performer for TSW across all environments. For GYLD, the best genotypes were Cool-0035 at E1 and E4, Numan at E2 and E7, Dosha at E3 and E6, Cool-0018 at E5, and Cool-0031 at E8. Additionally, the graph also indicated that all environments have high potential for both traits as the quadrant I and IV environments better than the quadrant II and III environment (low potential environment) (Figure 5) [10].

Figure 5: Polygon view of biplot 3 (Which-Won-Where) for TSW (I) and GYLD (J) of 14 faba bean genotypes under 8 environments.

Discriminativeness versus representativeness of GGE biplot

Discriminativeness vs. representativeness’ view of GGE biplots for TSW and GYLD traits being investigated in this study were denoted as pattern K and L, respectively. The discriminative ability of a test location described by the length of its location vector, which serves as an estimate of the standard deviation for that location. Longer vectors indicate greater discriminative ability. From this study, for TSW E1 () was the most effective at distinguishing between genotypes, while E8 excelled in differentiating genotypes for GYLD. Conversely, E8 and E3 were the least effective for TSW, and GYLD, respectively. The GGE biplots revealed that E2 and E3 had the smallest angles from the average environment for TSW and E3 and E5 for GYLD, indicating that these environments are more similar to each other and to the average environment in their genotype differentiation capabilities. In contrast, E8 and E4 had the largest angle from the average environment for TSW and E1 and E8 for GYLD, suggesting greater variability in terms of genotypes performance (Figure 6) [11].

Figure 6: Discriminativeness versus representativeness of GGE biplot for TSW (K) and GYLD (L) of 14 faba bean genotypes under 8 environments.

Genotype ranking: best genotype assessment

The application of gge-biplot analysis enabled the determination of optimal and most desirable genotype from a set of genotypes. ideal genotype is consistently located within the middle region and in close proximity to the peak of the arrow within the circular band. Accordingly, from this study for TSW, it was noted that variety Numan was located within the inner circle and considered to be optimal. For GYLD, Dosha and Numan variety exhibited closeness to the inner circle. In contrast, Cool-0011 exhibited the greatest distance from the arrowhead in the plot for both TSW and GYLD (Figure 7) [12].

Figure 7: Ranking genotypes based on PC1 and PC2 of TSW (M Pattern) and GYLD (N pattern) showing G × E interactions of the 14 faba bean genotypes under 4 locations and two seasons (8 environments).

Best linear unbiased prediction (BLUP)

Overall performance of genotypes presented as BLUP values indicated that Dosha, Cool-0030 and Numan varieties were ranked the highest for TSW.  whereas, for GYLD, Numan, Dosha, Cool-0018, Cool-0030, Cool-0024, Cool-0035, Cool-0031, and Cool-0034 scored above average with better yield stability across environments. Conversely, Cool-0011 performed poorly in both TSW, and GYLD (Figure 8).

Figure 8: Best linear unbiased prediction (BLUP) for 14 faba bean genotypes evaluated under 8 environments for TSW on the left (a) and GYLD on the right side (b).

Comparison of AMMI family and blup models (Cross Validation)

Stability analysis was done using AMMI, and GGE biplot for the TSW and GYLD traits. Accordingly, the variability explained by the AMMI model for TSW was 85%, and GGE biplot was 95.9%. for GYLD, the variability explained by the AMMI model was 52.3% and GGE biplot was 73.38% (Figure 9) [13].

Figure 9: Predictive accuracy of the additive main effects and multiplicative interaction (AMMI) family and best linear unbiased prediction (BLUP) for TSW and GYLD.

Conclusion

Faba bean, a crucial cool-season grain legume, is grown in over 70 countries worldwide. Ethiopia is the second-largest growing country and the first producer in Africa. The crop is mainly cultivated in mid- and high-altitude areas, providing food and feed to small-holder farming communities and providing foreign exchange and income for farmers. However, faba bean yield performance is unstable and affected by environmental variations. To increase productivity, breeders test large numbers of genotypes in various environments to evaluate yield stability and wide adaptability. The study of genotype by environment interaction (GEI) is essential for crop improvements, helping breeders generate information on adaptation patterns and new varieties for release. GEI can be computed using the additive main effects and multiplicative interaction (AMMI) model, which combines variance analysis for additive effects and principal components analysis for multiplicative effects. BLUP and AMMI are two distinct approaches to achieving the same goal, but they differ statistically.

The experiment was conducted during the main cropping season in 2018 and 2019 at four locations in Ethiopia's Oromia Regional States. Fourteen faba bean genotypes were evaluated under potential environments using a randomized complete block design (RCBD) with four replications. Data on grain yield, agronomic conditions, and disease reactions were collected from each plot. Days to flowering and days to physiological maturity were taken when each plot reached 50% of flower initiation and 90% of the pod attained physiological maturity, respectively. Responses of genotypes to disease reactions like chocolate spot, ascochyta blight, and rust were recorded at late pod setting. Grain yield was measured on clean, dried seed, and plot yields were adjusted to 10% moisture level and converted to kilograms per hectare. Statistical analysis was performed using a mixed linear model to assess differences among genotypes. The Additive main effects and multiplicative interaction (AMMI) model was performed for grain yield and thousand seed weight of 14 faba bean genotypes.

The study analyzed 14 faba bean genotypes over four locations and two years, revealing that all traits were significantly influenced by both genotypes and environments. The genotype by environment interaction (GEI) significantly affected traits such as days to 50% flowering, number of pods per plant, number of seeds per pod, thousand seed weight, grain yield, and severity of rust disease. The AMMI analysis revealed that environmental factors had the most significant impact on thousand seed weight, while genotype had a greater influence on grain yield. The variance due to GEI had the smallest impact on phenotypic variance for both traits. Seven interaction principal component axes (IPCA) were identified through principal component analysis (PCA), representing how different interactions affect the traits of interest (TSW and GYLD). Two IPCAs were significant at the 1% probability level for TSW, indicating that interactions affecting TSW are primarily captured by these two principal components. The AMMI model analysis revealed that PC1 and PC2 together accounted for 85.3% of the total variation in thousand seed weight (TSW).

1. Bocianowski J, Niemann J, Nowosad K (2019) Genotype-byenvironment interaction for seed quality traits in interspecific cross-derived Brassica lines using additive main effects and multiplicative interaction model. Euphytica 215:1-13.

Google Scholar

2. Cernay C, Ben-Ari T, Pelzer E, Meynard JM, Makowski D, et al. (2015) Estimating variability in grain legume yields across Europe and the Americas. Scientific Reports 5:11171.

Indexed at, Google Scholar, Crossref

3. Dereje Abera (2019) Analysis of Genotype by Environmental Interaction and Stability for Grain Yield of Faba Bean (Vicia Faba L.) Genotypes in West and Kellem Wollega Zones of Western Oromia. International Journal of Research Studies in Agricultural Sciences (IJRSAS) 5: 37-43.

Indexed at, Google Scholar, Crossref

4. Nardino M, Baretta D, Carvalho IR, Olivoto T, Follmann DN, et al.(2016) Restricted maximum likelihood/best linear unbiased prediction (REML/BLUP) for analyzing the agronomic performance of corn. African J Agric Res 11:4864-4872.

Google Scholar

5. Olivoto T, Nardino M, Carvalho IR, Follmann DN, Ferrari M, et al. (2017) REML/BLUP and sequential path analysis in estimating genotypic values and interrelationships among simple maize grain yield-related traits. Genet Mol Res 16: gmr16019525.

Indexed at, Google Scholar, Crossref

6. Piepho HP (1994) Best linear unbiased prediction (BLUP) for regional yield trials: a comparison to additive main effects and multiplicative interaction (AMMI) analysis. Theoretical and Applied Genetics 89:647-654.

Indexed at, Google Scholar, Crossref

7. Reckling M, Döring TF, Bergkvist G, Stoddard FL, Watson CA, et al. (2018) Grain legume yields are as stable as other spring crops in long-term experiments across northern Europe. Agronomy for Sustainable Development 38:63.

Indexed at, Google Scholar, Crossref

8. Tamene T Tolessa (2015) Application of AMMI and Tai’s Stability Statistics for Yield Stability Analysis in Faba bean (Vicia faba L.) Cultivars Grown in Central Highlands of Ethiopia. Journal of Plant Sciences. 3: 197-206.

Google Scholar

9. Veenstra LD, Santantonio N, Jannink JL, Sorrells ME (2019) Influence of genotype and environment on wheat grain fructan content. Crop Sci 59:190-198.

Google Scholar, Crossref

10. Yan W (2015) Mega-environment Analysis and Test Location Evaluation Based on Unbalanced Multiyear Data. Crop Sci 55:113-122.

Google Scholar

11. Zobel RW, Wright MJ, Gauch G (1988) Statistical analysis of a yield trial. Agronomy Journal80: 388-393.

Indexed at, Google Scholar, Crossref

12. Purchase J Hesta & Deventer C (2000) Genotype × environment interaction of winter wheat (Triticum aestivum L.) in South Africa: II. Stability analysis of yield performance. South African Journal of Plant and Soil 17: 101-107.

Indexed at, Google Scholar, Crossref

13. Girmay, Zeyede, Mahari, Milite (2015) Knowledge, attitude and practices of peasants towards hyraxes in two selected church forests in Tigray Region, Northern Ethiopia. International Journal of Biodiversity and Conservation 7: 299-307.

Indexed at, Google Scholar, Crossref