Genotype × environment interaction and yield stability of Arabica coffee ( Coffea arabica L . ) genotypes

Lack of suitable varieties that exhibit stable yield performances across wide ranges of environments is the major factor among several production constraints contributing to low productivity of Arabica coffee in Ethiopia. Eleven advanced Limmu coffee genotypes were evaluated in eight environments (four locations over two years) to determine the existence of GEI and yield stability performances. The experiment was laid out in a randomized complete block design of two replications under all locations. Combined analysis of variance showed a highly significant effect of genotype by environment interaction indicating the differential yield response of genotypes across different environments. The major proportion of the variation explained by environments was 42.74% of the total variation. Nevertheless, the contribution of the genotypes to the total variance was much smaller than the environments, and the genotype by environment interaction. Different stability models such as additive main effect and multiplicative interaction (AMMI), AMMI stability value, cultivar superiority index and yield stability index were used for stability analysis. The first two Interaction Principal Component Axis (IPCAs) of AMMI exhibited a highly significant effect and cumulatively contributed about 63.21% of the total interaction sum of squares. Two high yielding genotypes, namely (L52/2001) and (L55/2001), on average, showed stable performance across environments. On the other hand, the study also illustrated the presence of location specific high yielding coffee genotype such as L56/2001. Regarding the test environments, Gera 2015/16 (E5) is considered as a more stable site over the rest environments, while Agaro 2015/16 (E7) was considered to be the most interactive environment. Based on the result of the study, coffee breeders or farmers would be recommended for wise selecting either for location specific or wider adaptable coffee genotypes leading to substantial yield increase under Limmu coffee growing areas.


INTRODUCTION
Arabica coffee is the most widely consumed and highly preferred international beverage mainly for its best quality and is also one of the most important agricultural commodities in the world contributing to more than 60% *Corresponding author. E-mail: lbeksisa@gmail.com. Tel: +251910822464.
Author(s) agree that this article remain permanently open access under the terms of the Creative Commons Attribution License 4.0 International License of the world coffee production (Van der Vossen and Bertrand, 2015). Particularly in Ethiopia, coffee cultivation plays a fundamental role both in the cultural and socioeconomic life of Ethiopians. It represents the major agricultural export crop, providing 20 to 25% of the foreign exchange earnings (ECFF, 2015). The coffee sector contributes about 4 to 5% to the country's Gross Domestic Product (GDP) and creates hundreds of thousands of local job opportunities (EBI, 2014).
Ethiopia is the largest producer of coffee in sub-Saharan Africa and is the fifth largest coffee producer in the world next to Brazil, Vietnam, Colombia and Indonesia, contributing to about 7 to 10% of total world coffee production (Gray et al., 2013). The total area coverage of coffee in Ethiopia is estimated to be about 800,000 ha of land with an annual production capacity of 500,000 tons of which about 95% is produced by 4 million small scale farmers (Berhanu et al., 2015). Despite the high genetic diversity of Arabica coffee and naturally suitable climate condition in Ethiopia, the coffee production and productivity is not yet fully improved. Lack of high yielding improved varieties for each agro ecological zones and lack of suitable varieties that exhibit stable performance across wide ranges of environments are the major constraints in coffee production and productivity in Ethiopia (Bayetta, 2001;Yonas and Bayetta, 2008).
Crop performance is the product of the interaction between the genotype and the environment in which the crop is grown (Acquaah, 2007). The differential response of genotypes to environmental changes refers to genotype x environment interaction (GEI) (Crossa et al., 1990). Significant G x E interaction represents a major challenge to plant breeders to fully understand and obtain the genetic control of variability (Luthra and Singh, 1974). Measuring and understanding G x E interaction and stability performance of genotype should be an essential component in plant breeding programs for the decision making process such as identification of the most relevant testing environments, allocation of resources within a breeding program, and choice of germplasm and breeding strategies (Lean et al., 2016).
In fact, Jimma Agricultural Research Center (JARC) has been conducting coffee research work for about five decades to improve production and productivity in the country. As a result, about 34 improved pure lines and six hybrid coffee varieties were released for the various major coffee growing agro-ecologies of the country. Although, many high yielding and disease resistance varieties were developed and released, the released varieties did not relate to the whole coffee growing region due to the existence of vast and divers agro ecologies in the country. It was shown that, coffees grown under these environments are different in quality, disease resistance and yield potential. Therefore, development of varieties that have the potential for wider adaptation would be of paramount importance to overcome the Beksisa et al. 211 shortage of improved varieties in the potential coffee growing regions of the country. In Ethiopia, to increase production and productivity of coffee using stable varieties, the first adaptation tests across different environments was carried out by Mesfin and Bayetta (1987). They reported that some of the genotypes showed poor adaptation to major coffee production areas outside the high land forest which showed location specific nature of Arabica coffee genotypes. Similarly, Yonas and Bayetta (2008), Meaza et al. (2011) andYonas et al. (2014) also confirmed that varieties that exhibit better adaptation at one location did not perform well at other locations. On the other hand, some reports at the same time stated the presence of high yielding genotypes with regular responses to changes in environment. However, in spite of the G x E interaction impacts on production and productivity of coffee yield, due emphasis have not been given on this particular area of investigation. Based on this evidence, testing coffee genotypes across environments has vital prominence before deciding either for specific or extensive uses of genotypes. Therefore, this study was undrtaken with the objective to determine the existence of G x E interaction and stability performance of Arabica coffee genotypes for bean yield.

Description of experimental site
The trials were conducted in different major coffee producing agro ecological zones of southwestern Ethiopia, Oromia Regional state at four specific places: Jimma, Agaro, Manna and Gera for two consecutive cropping seasons (2014/15 and 2015/16). The first three locations represent mid altitude, while Gera represents high land area. Description of the testing locations with some of their climatic and soil characteristics are presented in Table 1.

Experimental materials
The experimental materials used in this study comprised 11 Arabica coffee (Coffea arabica L.) genotypes. The genotypes were the only common genotypes at all locations; obviously the trials consisting of thirteen genotypes and therefore, the genotypes which did not exist at all locations were not incorporated in this study. The genetic materials were selected from Limmu Kossa and Limmu Seka collection of 2001 based on their yield, cup quality and disease resistance during the initial investigation at Gera and Agaro research centers. The geographical origin and description of the experimental materials are presented in Table 2.

Experimental design and management
Randomized complete block design (RCBD) was used and treatments were replicated two times in each location. The experimental plots consisted of ten trees 2 m x 2 m between rows and between plants have a plant density of 2500 per hectare. All field management practices were conducted according to the recommendations for the crop in the region (Endale et al., 2008).

Data collected
Total fresh cherry yield was harvested and recorded in grams from ten trees in a plot and used to compute mean yield per each tree. Clean coffee yield in kg/ha was obtained by multiplying the yield of the fresh cherry by percent out-turn.

Statistical analysis
Analysis of variance (ANOVA) was done for each location separately based on the the standard procedure developed for a randomized complete block design. Bartlett's (1974) test was used to determine the homogeneity of error variances between environments. Comparison of treatment means was done using least significant difference (LSD). A combined analysis of variance was done to determine the significant effects of the genotypes, environments and their interactions. The SAS version 9.2 (SAS, 2008) statistical software was used for statistical computations and estimation of differences among genotypes. The effects of the genotypes and environments as well as their interactions were determined from ANOVA. Analysis of genotype stability across eight environments (locations and years) was computed using Additive Main effect and Multiplicative Interaction (AMMI), AMMI stability value (ASV), cultivar superiority index (Pi) and yield stability index. The detail of each stability model is separately presented as follows: AMMI model which combines standard analysis of variance with principal component analysis (PCA) analysis was used to investigate genotype x environment interaction. To show a clear insight into specific GE interaction combinations and the general pattern of adaptation, a biplot of genotype and environment interaction (Kempton, 1984) was developed. In the AMMI1 biplots, the first IPCA was used as ordinate (Y-axis) and the main effects or means of genotypes and environments represented abscissa (Xaxis). AMMI2 biplot is generated using genotypic and environmental scores of the first two AMMI components. The AMMI analysis was done using GenStat version 16 th software according to the model suggested by Crossa et al. (1990).
Where: i=1, 2…….11; j=1, 2….8; Yij is the performance of the i th genotype in the j th environment; µ is the grand mean; Gi is additive effect of the i th genotype (the genotype deviation from the grand mean); Ej is additive effect of the j th environment (the environment deviation from the grand mean); Kn is Eigen value of the IPCA axis n; Uni and Snj are score of genotype i and environment j for the IPCAs; Qij is residual for the multiplicative components; eji is random error. AMMI stability value (ASV), which is stability value based on the AMMI model's IPCA1 and IPCA2 values for each genotype and each environment was calculated as suggested by Purchase (1997). The larger the ASV value either negative or positive, the more the genotypes specifically adapted to certain environments. Conversely, lower ASV values indicate greater stability of genotypes to different environments (Purchase, 1997).
Where, is the weight resulting from dividing the sum of IPCA1 squares by the sum of IPCA2 squares.
Cultivar superiority index (Pi) was done using GenStat version 16 th software as described by Lin and Binns (1988). Mathematically, the value of Pi was obtained as follows: Where, Y is the yield means of the i th genotype in the j th location, Yij max is the yield mean of the genotype with maximum yield in the j th environment and nl is the number of environments. A small Pi value indicates a better fit of a genotype to this stability concept.
The new approach known as yield stability index (YSI) that was developed by Mahmodi et al. (2011) is recommended as a measure of genotype stability. YSI incorporate both mean yield and stability Beksisa et al. 213 in a single criterion. Low value of this parameter shows desirable genotypes with high mean yield and stability.YSI is calculated as: Where, RASV is the rank of AMMI stability value and RY is the rank of mean yield of genotypes across environments.

RESULTS AND DISCUSSION
The combined analysis of variance for yield of the tested coffee genotypes is shown in Table 3. The result showed that testing environments testing environments were significantly different at P<0.05. The significant difference observed between environments indicates that the mean bean yield value of genotypes differed from one environment to another due to dissimilarity of environments. In this case, unpredictable variations such as fluctuating features of the environment for instance; rainfall, relative humidity, temperature and soil characteristics might cause differential performances of genotypes from one environment to another environment. Fehr (1993) also reported that every factor which is a part of the environment has the potential to cause differential performances of the genotypes that is associated with genotype x environment interaction. Similar findings on the existence of genotype x environments were also reported by many authors (Mesfin and Bayetta, 1987;Meaza et al., 2011;Yonas et al., 2014). The genotypes also revealed highly significant difference (P<0.01). The significant difference among the genotypes demonstrated the presence of variability in the inherent genetic constitute of the Coffea arabica L. genotypes tested. The result of the study also indicated the existence of genotype by environment interaction and the interaction was highly significant (P<0.01) reflecting the differential response of genotypes in various locations and seasons. This variation could be attributed to differences in climatic and edaphic conditions at different testing environments.
In the presence of the G x E interaction, the phenotypic expression of one genotype might be superior to another genotype in one environment but inferior in another environment. Hence, such presence of a significant G x E interaction complicates breeding strategy because superiority of genotypes across environments cannot be identified by considering their mean performance and the need to develop genotypes that are adapted to specific environmental conditions or the need to identify genotypes that are exceptional in their stability performances across environments. Similarly, the significant effect of G x E interaction in Ethiopia in different types of the quantitative traits of Arabica coffee was reported by previous researchers (Mesfin and Bayetta, 1987;Meaza et al., 2011;Yonas et al., 2014). Agwanda and Owuor (1989) and Agwanda et al. (1997) with Arabica coffee and Montagon et al. (2000) with Coffea canephora also reported the presence of  significant genotype × environment interaction in bean yield and yield components.

Additive main effects and multiplicative interaction (AMMI) analysis for bean yield
In this study, the estimated magnitude of different variance for yield variation of the genotypes tested across environments showed that the largest portion of variation was accounted to environment contributing to the bean yield of Arabica coffee with 42.75% (Table 4). The large sum of square and highly significant mean square of environment indicated that the environments have significant influence on bean yield performances of the genotypes tested. Genotypes and G × E interaction accounted for 9.31 and 32.32% of the total variation explained, respectively (Table 4). The current finding indicated that bean yield of genotypes was found to be significantly affected by changes in the environment, followed by G × E interaction and genotypic effect. Thus, the large differences among environmental means causing most of the variation in bean yield of Arabica coffee was mainly due to environments. The result is in agreement with the discoveries of Yonas and Bayetta (2008), Meaza et al. (2011) and Yonas et al. (2014) who reported the significant influence and/or largest portion of environments on Arabica coffee bean yield performances.
Genotype by environment interaction effects were further partitioned into six possible interaction principal component axes (IPCA) along their contribution of sum of squares with decreasing importance (Table 4). Among these, the first four IPCAs exhibited highly significant difference (P<0.01). The first and second interaction principal component axis (IPCA) explained 38.15 and 25.06% of the total variation accounted by the G × E interaction sum of squares, respectively (Table 4). The third and fourth interaction principal component axis (IPCA3 and IPCA4) explained 17.73 and 11.43% of sum of squares of G × E interaction, while the first four IPCAs cumulatively explained 92.37% of sum of squares of G × E interaction (Table 4). Meaza et al. (2011) also reported significance of first four IPCAs for Arabica coffee bean yields which is in agreement with the current finding. The cumulative sum of squares of the first two IPCAs accounted for a total of 63.21% of the interaction. The first six interaction principal component axis (IPCA 1-6) accounted for 98.95% of total G × E interaction, leaving 1.05% of the variation in the residual. The first two principal components showed sum of squares greater than half of all and evaluation using F-test revealed highly significant P<0.01, indicating the capability of the first two principal components axis for cross-validation variation explained by G × E interaction (Zobel et al., 1988;Gauch and Zobel, 1996).

AMMI 1 biplot for yield
In Figure 1, the mean yields of the genotypes grown across different environments, the environment means and the first IPCA scores can be clearly understood. Based on the biplot analysis, genotypes or environments with large IPCA1 scores, either positive or negative had large interactions, whereas genotypes with IPCA1 score of zero or nearly zero had smaller interactions and was considered as stable over wide range of environments (Crossa et al., 1990;Gauch and Zobel, 1996). Accordingly, genotypes G4 (L45/2001), G9 (L55/2001), G7 (L63/2001) and G11 (L67/2001) had low IPCA1 value closest to zero score indicating that these genotypes were more stable than other C. arabica L. tested ( Figure  1). However, for genotypes to be stable or generally adaptable to all environments, the genotypes should attain above average mean performance and the IPCA score would be nearly zero. Therefore, genotypes G9 (L55/2001), G7 (L63/2001) and G4 (L45/2001) registered above average yield together with the IPCA1 score close to zero, whereas, G1 (L68/2001) with low average mean performances was the most unstable as its IPCA1 score is largest when compared with the others, while G9 (L55/2001) showed the most stable performance than the rest C. arabica L. genotypes tested.
Environments with IPCA score located farther away from the origin in the biplot interacted more with the genotypes and made the selection difficult. In this study, E7 was high yielding environment but the most interactive as its IPCA1 score is largest when compared with the others. On the other hand, the environments E8, E5 and E4 had IPCA1 score close to zero (Figure 1), but low Beksisa et al. 215 yielding environments, except E5 which have above average yield. E5 (Gera) is characterized by high altitude, high rainfall, cool temperature and long maturity period; thus, genotypes constantly exploit their genetic potential giving average mean yield under this particular location. Moreover, this biplot also indicated E7 as the highest yielding environment and E8 as the lowest yielding environment. In general, high yielding environments were sparsely distributed in quadrant II (E1, E5 and E6) and III (E2 and E7), while the lower yielding environments were sparsely distributed in quadrant I (E3, E4 and E8); but none of the environment was plotted in quadrant IV (Figure 1). In addition, genotypes plotted in quadrant II G4 (

AMMI 2 biplot for yield
The AMMI 2 biplot (Figure 2) was generated using the genotype and environment scores of the first two AMMI components (Vargas and Crossa, 2000). The first interaction principal component captured 38.2%, while the second interaction principal component captured 25.1% of the total G x E interaction sum of square. The first two IPCAs cumulatively captured 63.3% of sum of square of the G × E interaction of tested coffee genotypes. From earlier yield trial of G × E interaction in Arabica coffee, Yonas et al. (2014) reported that the first IPCA alone accounted for 36% of the total interaction sum of square and Meaza et al. (2011) reported that AMMI with the first two IPCAs explained 74% of the total interaction sum of squares.
According to Purchase (1997), the genotypes and environments that are located far away from the center are more responsive or unstable, while genotypes that are closer to the center of biplot have higher stability performance. Hence   G1=L68/2001, G2=L01/2001, G3=L52/2001, G4=L45/2001, G5=L54/2001, G6=L03/2001, G7=L63/2001, G8=L32/2001, G9=L55/2001, G10=L56/2001, G11=L67/2001, E1=Gera 2014/15, E2=Jimma 2014/15, E3=Agaro 2014/15, E4=Manna 2014/15, E5=Gera 2015/16, E6=Jimma 2015/16, E7=Agaro 2015/16 and E8=Manna 2015 plotted relatively close to the center designating their minimum involvement in the total G x E interaction sum squares and considered as stable genotypes. However, for genotypes to be considered as stable, it should attain high mean performance having greater than grand mean. Therefore, only G7 (L63/2001) and G3 (L52/2001) could be considered as most stable genotypes with their high bean yield performance and being closer to the origin as compared to the others. Whereas, genotypes G1 (L68/2001) and G4 (L45/2001) were farthest from the center of biplot having substantial involvement in G x E interaction sum squares. Therefore, these genotypes were considered as unstable genotypes. Similarly, E4 and E8 can be considered as stable environments due to closeness of its vector end points to the center of biplot. In contrast, the farther away from the center of biplot for the environments, the more interaction the environment has with genotypes. As was already identified, E7 was the most interactive environment on AMMI1 biplot, AMMI2 biplot also identified E7 and E5 as the most interactive environments as it was farthest from the center of biplot. In addition, the response of the locations for the genotypes performance was different from one season to the other. Such imbalance genotypes performance over the two seasons is largely attributed to the very conducive or unconducive environment prevailing at all locations during the experimental period.

Cultivar superiority index (Pi)
The superiority index (Pi) values ranged from 204137 to 1202144 (Table 5) which indicated large differences among tested genotypes for this stability model. Abtew et al. (2015) also reported the large differences among tested wheat genotypes under their investigation. Regarding superiority index (Pi), the genotypes with minimum Pi-value could be considered as stable (Lin and Binns, 1988). Accordingly, the high yielding genotypes, namely G3 (L52/2001), G9 (L55/2001) and G7 (L63/2001) displayed the greatest yield performance and the lowest Pi-values. On the other hand, genotypes with maximum Pi-value could be considered as most unstable. Therefore, genotype G8 (L32/2001) but with maximum Pi-value; could be considered as most unstable genotype. The strong association between mean bean yield and Pi was expected because the values of this stability parameter were high for highyielding genotypes, that is, the top-ranking in bean yield could also be the top-ranking in this stability parameter. In the dynamic concept of stability such as in cultivar superiority model, it is not required that the genotype response to environmental conditions should be equal for all genotypes (Becker and Léon, 1988). This type of stability parameter is preferable for commercial farming with high mean yields and the potential to respond to agronomic inputs or better environmental conditions rather than for resource poor farmers who prefer lower but stable optimal environmental conditions and inputs.

AMMI stability value (ASV)
The ASV parameter is used to quantify and classify the genotypes according to their stability performance. In this model, genotypes with least ASV or have smallest distance from the origin are considered as the most stable, whereas those which have highest ASV are considered as unstable (Purchase, 1997). Accordingly, genotypes G7 ( (Table 5). AMMI stability value which defines stable genotypes by the distance of the genotypes from the zero point of the IPCA1 vs. IPCA2, is consistent with the AMMI2 model but have a little relationship with AMMI1 model which has only genotype G7 (L63/2001) in common. Genotypes, G3 (L52/2001) and G9 (L55/2001) which ranked first and third in their bean yield performance according to their order, are considered as moderately stable by this stability model.

Yield stability index (YSI)
Yield stability index (YSI) proposed by Mahmodi et al. (2011) incorporates both stability and yield performance in one criterion. AMMI stability value (ASV) takes into account both IPCA1 and IPCA2 that justify most of the variation in the G × E interaction. The rank of ASV takes the rank one, at the same time the highest yield mean takes the rank one and then the ranks are summed in a single simultaneous selection index of yield and yield stability called yield stability index (YSI). The genotypes with low YSI would be considered as high yielding and stable genotypes. Hence, YSI identified G3 (L52/2001) G7 (L67/2001) and G9 (L55/2001) as most stable genotypes, whereas G1 (L68/2001) was identified as least stable (Table 5). This stability parameter was also used by Tadesse and Abay (2011) in sesame genotypes to describe the stability performance of genotypes studied.

Conclusion
In Ethiopia where coffee production plays a major role in the national economy, yield fluctuation and yielding pattern of coffee being varied with small geographic variation and thus, attributed to low productivity. To this effect, assessments of the stability as well as the performance of coffee genotypes across diverse environmental conditions are important for selection of wider adaptable or superior genotypes for the target environments before variety release. In this study, eleven Arabica coffee genotypes which were common at all locations were evaluated at different agro-ecologies of southwestern Ethiopia; at eight environments (four locations for two cropping seasons) to determine the existence of G × E interaction and yield stability performances.
Combined analysis of variance exhibited highly significant difference among the genotypes. The finding showed significant effects of both environments and G x E interaction. The major proportion of the total variation in bean yield was explained by environments (42.75%) followed by G x E interaction (32.32%) and genotypes (9.31). The finding indicated that the genotypes G3 (L52/2001) and G9 (L55/2001) with high mean yield of 1558 and 1473 kg/ha -1 , respectively proved to be the best in stability among the studied genotypes. On the other hand, environment, Gera 2015/16 (E5) showed average response to all genotypes, while Agaro 2015/16 (E7) exhibited non-additive behavior. Therefore, this study clearly indicated the possibility of exploiting the yield potential of Limmu coffee genotypes under its growing conditions either by using wider adaptable coffee types or location specific high yielder genotype under favorable environmental condition.