Stability, the Last Frontier: Forage Yield Dynamics of Peas under Two Cultivation Systems

The stability of performance may be proved to be the last frontier for adopting certain genotypes in various cultivation systems and environments. The main objective of the present study was to analyze the forage yield stability performance of seven pea (Pisum sativum L.) genotypes based on various stability indices. The genotype behavior was studied based on the yield of peas under both conventional and low-input cultivation systems. Five cultivars of peas (broadly distributed) and two lines were used in a strip-plot design. Significant positive correlations were detected between forage yield and some other traits. This way, forage yield stability may be indirectly improved by improving certain traits showing qualitative inheritance. Comparisons revealed that genotypes exhibited stable performance, even in low-input farming systems. AMMI analysis, GGE biplot, and analysis of variance (ANOVA) combination showed statistically significant differences between genotypes and environments and the farming system. Our analysis depicted specific cultivars of peas for different areas and farming systems to attain highly stable performance. Vermio was confirmed to be a stable genotype for forage yield performance in low-input farming in Trikala and Kalambaka areas, while Pisso was indicated as the best in Florina and Giannitsa areas in low-input farming. The two pea lines exhibited stable performance in Giannitsa and Florina areas, especially in low-input conditions. The stable behavior of some genotypes in these conditions may be useful for farmers that raise livestock in mountainous areas. The genetic parameters show that the selection for fresh forage yield and dry matter yield in breeding programs is expected to be effective.


Introduction
The pea (Pisum sativum L.) is a very useful crop for supporting livestock [1,2]. Pea cultivations can be found in a wide range of agro-climatic zones; the potential of this crop and its high nutritional value are referred to in many related publications [3,4]. Field pea is often used as a main protein source, as it has a rich and unique protein profile, different from other natural protein sources [5]. On the basis of the aforementioned points, our study focused on estimating the genetic potential needed to support the growth and yield performance of pea cultivations in various environmental conditions and, subsequently, to propose the best cultivars.
The stability of performance could be crucial for adopting certain genotypes in various cultivation systems and environments [6]. Stability may be dependent on cultivar tolerance Peas are an interesting cultivar for animal feeding; thus, the highly productive varieties combined with a high quality of forage are preferable. The exploitation of the genetic parameters of the traits for breeding purposes is desirable. This knowledge is fundamental for effective breeding programs. Therefore, in order to initiate any breeding program, the exploitation of suitable parameters, such as the genotypic coefficient of variation GCV, heritability in a broad sense (H 2 ) is necessary.
In the present study, the main scope was to determine the forage yield stability of pea genotypes along with various correlated traits based on the innovative approach of estimating the stability index, with the specific aim of studying pea genotypes' behavior under both conventional and low-input cultivation systems. Greveniotis et al. [23,39] used a stability index, based on Fasoulas [7] and Fasoula [8] remarks, as an estimation of the heritability of various traits, leading to clear discrimination between qualitative and quantitative traits. Our approach includes stability performance analysis and reveals the stability performance and the kind of heritability of traits.

ANOVA and Descriptive Statistics on the Stability Index
Regarding the ANOVA table (Table 1), the main effects for all traits expressed significant differences. The G × E interaction showed significant differences for all traits. Multiple interactions involving genotypes, the environment, and the cultivation system were found to be very significant, especially for forage yield, and these data must be analyzed in combination with the genotype performance within each environment and cultivation system in order to define the best genotype for specific conditions. Days to flowering showed no multiple interactions. Main stem length showed no cultivation × environment interaction. Days to flowering showed no multiple interactions. Table 1. Mean squares (m.s.) from analysis of variance over environments and cultivation methods for tested traits: days to 50% flowering, main stem length (cm), main stem thickness (mm), fresh forage yield (kg ha −1 ), dry matter yield (kg ha −1 ), forage dry matter crude protein content (%), and ash content % of dry matter. To better analyze the performance of the genotypes in different environments and estimate the stability of each genotype for all traits, we used AMMI and GGE analysis as the most appropriate tools. For all traits, the GE was much more significant than G effects, and additional AMMI analysis was needed.

Source of Variation
Stability estimations based on the calculation of the stability index of each trait are presented in Tables 2-4. The tabulated stability index data across environments for the seven characteristics under study are listed in Table 2. Days to 50% flowering and main stem length showed generally high indices (over 1000) in many cases. Fresh forage yield showed low indices. The low-input farming system seems to improve stability indices in many cases for forage yield.
Days to 50% flowering showed some extreme values for certain areas; for example, it was found to be 38,704 for cv. Pisso in Kalambaka (cv. Livioleta was 27,429). Generally, it was the most stable trait, with values over 10,000 in many cases. The main stem length showed a very stable behavior, with values over 1000. For some cultivars or lines, there were a few extreme values over 10,000 (line Zt2, cv. Olympos and cv. Vermio)-but also line Zt2, cv. Vermio in Kalambaka, cv. Olympos in Trikala, and cv. Livioletta in Giannitsa and in Florina, respectively, depending on the environment or cultivation method. Evidently, interactions led specific genotypes to exhibit varying behaviors for stability according to the environment or cultivation method.

The AMMI Tool for Multi-Environment Evaluations
The AMMI model is a widely used statistical tool in the analysis of multi-environment experiments. The purpose of the tool is to understand the complex GEI. In the AMMI model, the data are represented by a two-way table of GEI means. In complete tables, least-squares estimation is equivalent to fitting an additive two-way ANOVA model for the main effects and applying a single value decomposition to the interaction residuals [40].
Using this statistical tool, AMMI software can generate the adaptation map and AMMI1 biplot, where one axis is the axis of the factor, and the other is the PC1 value. When the PC1 value and its distance from the X-axis are close, the factor analyzed is stable. Regarding the AMMI1 biplot, the desirable genotypes were those with a high value on the axis of trait performance (x-axis, right position) and close to the center of the PC1 axis (near zero).
GGE stands for the genotype's main effect (G) plus the genotype by environment interaction (GE), which is the only source of variation that is relevant to genotype evaluation. Mathematically, GGE is the genotype by environment data matrix after the environmental means are subtracted. A GGE biplot is a biplot that displays the GGE of a genotype through two-way environmental data. The GGE biplot methodology originates from the graphical analysis of multi-environment genotype trials (MET) data but is equally applicable to all other types of two-way data.
Regarding the GGE biplot for environments, the most stable environment was that placed close to the dot of the ideal and average environment and in the concentric area of the ideal environment dot. In terms of the GGE biplot for genotypes, the desirable genotypes (stable and productive) were those placed near the ideal genotype and in the concentric area of the ideal genotype dot.
The AMMI1 and G × E biplot analysis created biplots depicting the performance of the genotypes in different environments. The biplots created serve as a simple tool that can easily characterize each genotype for performance and stability.
The stability analysis using both AMMI and GGE biplot for days to 50% flowering is depicted in Figure 1a-d.  The stability analysis using both AMMI and GGE biplot for main stem length is depicted in Figure 2a-d.  Data from the main stem thickness used in AMMI and GGE biplot analysis (Figure 3a-d) show that this trait was very environmentally dependent.  The stability analysis using both AMMI and GGE biplot for fresh forage yield is shown in Figure 4a-d.  The stability analysis of dry matter yield using AMMI and GGE biplot is depicted in Figure 5a-d.  The stability analysis using AMMI and GGE biplot for forage dry matter crude protein content is presented in Figure 6a-   The stability analysis for the ash content % of dry matter trait, using the AMMI and GGE biplot, is presented in Figure 7a-  For AMMI analysis, as visualized by the adaptation map figure, the most desirable genotypes were those placed high on the axis of trait performance, showing a nearly parallel line to the PC1 axis, which was an indication of stability in different environments. For AMMI analysis, as visualized by the adaptation map figure, the most desirable genotypes were those placed high on the axis of trait performance, showing a nearly parallel line to the PC1 axis, which was an indication of stability in different environments.
For the AMMI1 biplot, the desirable genotypes were those placed high on the axis of trait performance (x-axis, right position) and close to the center of the PC1 axis (near the zero point).
Regarding the GGE biplot for environments, the most stable environment was that placed close to the dot of the ideal and average environment and in the concentric area of the ideal environment dot.
Concerning the GGE biplot for genotypes, the desirable genotypes (stable and productive) were those placed close to the ideal genotype and in the concentric area of the ideal genotype dot.

Exploratory Data Analysis of Peas
In order to estimate the phenotypic distances among genotypes, the clustering method of Ward was performed and formed clusters based on the traits tested. The clusters were formed based on the fresh forage yield and dry matter yield and the relations among them.

Genotypic and Phenotypic Coefficients of Variation and Heritability
In Table 5, estimations of genetic parameters for the traits are tabulated. The genetic parameters, along with the heritability in a broad sense, were estimated for all traits except the trait of main stem thickness. The parameters show that there is enough phenotypic variability for all traits. Furthermore, a large portion of phenotypic variability was genotypic, and this is desirable for geneticists in order to select superior genotypes for all traits. The heritability for all traits ranged from 99.4% to 83.8%. These estimates of heritability combined with the high percentage of genetic variability to the phenotype and the high diversity for all traits indicates that the selection of new varieties would be effective. Table 5. Estimations of genetic parameters for tested traits: days to 50% flowering, main stem length (cm), main stem thickness (mm), fresh forage yield (kg ha −1 ), dry matter yield (kg ha −1 ), forage dry matter crude protein content (%), and ash content % of dry matter. sd-standard deviation, σ 2 g -genotypic variance, σ 2 p -phenotypic variance, GCV-genotypic coefficient of variation, PCV-phenotypic coefficient of variation, and H 2 -broad sense heritability (%).

Correlations between All Characteristics
In Table 6, correlations between all traits are tabulated. Many correlations were statistically significant, especially between forage yield and traits such as the main stem length (r = 0.203), dry matter yield (r = 0.974), and forage dry matter crude protein content (r = 0.100). Table 6. Correlations between all traits measured: days to 50% flowering, main stem length (cm), main stem thickness (mm), fresh forage yield (kg ha −1 ), dry matter yield (kg ha −1 ), forage dry matter crude protein content (%), and ash content % of dry matter.

Discussion
Farmers and breeders need both high and stable performance regarding forage yield. In our work, the two cultivation systems (conventional and low-input) displayed differences in genotype-yielding performance, but overall estimations on various pea characteristics seemed to be unaffected. In combination with GGE biplot analysis, the two farming systems revealed the most stable genotypes across all environments, as well as those more stable in specific environments and farming systems. Additionally, some genotypes exhibited stability in low-input conditions. Generally, very significant GGE interactions were recorded. Sayar and Han [41], based on ANOVA findings, state that G × E interaction is the most important concept to deal with. In our work, G × E interaction was revealed due to multiple interactions recorded for many traits. Sayar's work [27] was based on AMMI analysis in order to define cultivar interactions with the environment. We described the interactions of each trait of pea cultivars and lines across different environments as follows:

Days to 50% Flowering
Regarding days to 50% flowering, AMMI analysis produced the figures adaptation map (Figure 1a) and AMMI1 biplot (Figure 1b). Both figures explained a portion of the total variability (71.5%), which is high enough for the genotype × environment (Gx) variation. Both the adaptation map and AMMI1 figures show that the most stable genotypes for environments E1 (Gianitsa) and E2 (Florina) were G6 (Zt1) and G5 (Dodoni), the late genotypes, whereas, for E3 (Trikala) and E4 (Kalambaka), the most stable genotypes were the early ones, G1 (Olympos), G2 (Pisso) and G7 (Zt1). The GGE analysis explained a total variability of 98.8% (PC1:96%, PC2: 2.8%), which was very high. The GGE biplot of the environment (Figure 1c) shows that all environments were quite similar and in the concentric circles of the ideal environment. The GGE biplot for the genotype view (Figure 1d) shows that all genotypes were very stable in all environments; the early genotypes were G1 (Olympos) and G2 (Pisso), and the late genotypes were G6 (Zt1) and G5 (Dodoni). The ideal for cultivation genotypes depends on what is desirable among early and late ones.

Main Stem Length (cm)
Regarding the main stem length, AMMI analysis produced the figures adaptation map (Figure 2a) and AMMI1 biplot (Figure 2b). Both figures explained a portion of the total variability (62.5%), which is high enough for conclusions. Both the adaptation map and AMMI1 figures show that the most stable genotypes were G2 (Pisso) and G1 (Olympos), where G2 (Pisso) had the highest performance for the main stem length trait. The GGE analysis explained a total variability of 90% (PC1:66.2%, PC2: 24.6%), which is very high. The GGE biplot for the environment view (Figure 2c) shows that all environments were very diverse, where the E1 (Giannitsa) environment was very close to the average environment.
The GGE biplot for the genotype view (Figure 2d) shows that the most stable genotype and identical to the ideal genotype was G2 (Pisso), followed by the G1 (Olympos) genotype, which was very stable but with lower performance for this trait.

Main Stem Thickness (mm)
The AMMI analysis via the adaptation map ( Figure 3a) and AMMI1 biplot expressed the PC1: 48.4% of the total variability. In both figures, there no clear pattern was found for stability, but the genotypes G3 (Livioletta), G4 (Vermio), and G2 (Pisso) were relatively stable. The GGE biplot analysis explained 70.0% (PC1:38.5%, PC2:31.5%) of the total variability. The GGE biplot for the environment view (Figure 3c) shows that all environments were very diverse, and no environment was placed near the average environment. The GGE biplot for the environment view (Figure 3d) shows that relative stable genotypes were G3 (Livioletta), G5 (Dodoni), and G7 (Zt2), but all were placed out of the concentric circles of the ideal genotype.

Fresh Forage Yield (kg ha −1 )
The AMMI analysis explained a portion (57.1%) of PC1's total variability. Both the adaptation map ( Figure 4a) and AMMI1 biplot (Figure 4b) show that the relatively stable genotypes were G2 (Pisso), G6 (Zt1), and G1 (Olympos), while G2 (Pisso) had the highest fresh forage yield. The GGE biplot analysis expressed 93.3% (PC1:79.3%, PC2:14.0%) of the total variability. The GGE biplot for the environment view shows that E1 (Giannitsa) and E4 (Kalambaka) were close to the average environment, and all environments were very diverse. The GGE biplot for the genotype view shows that the most desirable genotype was G2 (Pisso), followed by G6 (Zt1), which was less stable, and G1 (Olympos), which was less productive than the other two but very stable. AMMI analysis assisted Sayar [27] in recommending the best cultivars for fresh forage yield in certain cultivation areas.

Dry Matter Yield (kg ha −1 )
The AMMI analysis as presented from the adaptation map and AMMI1 biplot figures explained the PC1: 62.7% of the total variability. Both the adaptation map ( Figure 5a) and AMMI1 biplot (Figure 5b) show that the most stable genotypes were G2 (Pisso), G6 (Zt1), and G1 (Olympos). The most productive genotype was G2 (Pisso), followed by G6 (Zt1) and G1 (Olympos). The GGE biplot analysis explained 95.5% (PC1:77.9%, PC2:15.6%) of the total variability. The GGE biplot for environment view shows that E1 (Giannitsa) and E4 (Kalambaka) were close to the average environment, and all environments were very diverse. The GGE biplot for the genotype view shows that the most desirable genotype was G2 (Pisso), followed by G6 (Zt1) which was less stable, and G1 (Olympos), which was less productive than the other two but very stable. Acikgoz et al. [14] investigated the dry matter, yield relations, and G × E interactions and concluded after a comparison of cluster and stability analyses that the stability analysis failed to recommend cultivars for different regions.

Forage Dry Matter Crude Protein Content %
The AMMI analysis explained the PC1: 78.2% of the total variability, which is quite high. Both the adaptation map ( Figure 6a) and AMMI1 biplot (Figure 6b) figures show that the most productive genotypes were G2 (Pisso), G6 (Zt1), and G5 (Dodoni), which showed relatively low stability. The GGE biplot analysis explained 97.2% (PC1:79.9%, PC2:17.3%) of the total variability. The GGE biplot for the environment view shows that E2 (Florina), E1 (Giannitsa), and E4 (Kalambaka) were close to the average environment. The GGE biplot for the genotype view shows that the most desirable genotype was G2 (Pisso), G6 (Zt1), and G5 (Dodoni), which showed relatively low stability. Only the G6 (Zt1) genotype was placed in the concentric region of the ideal genotype, which indicates relatively acceptable stability and performance for this trait.

Ash Content % of Dry Matter
The AMMI analysis explained the PC1: 68.8% of the total variability, which is quite high. Both the adaptation map ( Figure 7a) and AMMI1 biplot (Figure 7b) figures show that the most productive genotypes were G3 (Livioletta), followed by G6 (Zt1) and G2 (Pisso). The most stable one was the G3 (Livioletta) genotype. The GGE biplot analysis explained 93.6% (PC1:66.9%, PC2:26.7%) of the total variability. The GGE biplot for the environment view shows that E4 (Kalambaka) and E3 (Trikala) were close to the average environment. The GGE biplot for the genotype view shows that the most desirable genotypes were G3 (Livioletta) and G6 (Zt1), which showed relatively low stability. The G3 (Livioletta) genotype was placed nearly identical to the ideal genotype, which means that it has acceptable stability and performance for this trait.

Genotypic and Phenotypic Coefficients of Variation and Heritability
The traits of days to 50% flowering, the main stem length, main stem thickness, dry matter crude protein, and ash are components of fresh and dry forage yield. The traits of fresh forage yield and dry forage yield seem to have high variability, as described by the min and max of Table 5. The heritability estimate was 83.8% and 84.6%, respectively. These values are high [42]. The genetic variability and the GCV of these two traits are the highest among all other traits. This combination of high heritability (H 2 ) and high GCV is an indication that the variation among genotypes was largely due to the additive genetic part [43]. Abebe et al. [44] suggested that high heritability in these values, in a broad sense, indicate that the characters under study are less influenced by the environment in their expression. This means that the direct selection of the traits of fresh forage yield and dry forage yield could be effective. As far as the other traits, the heritability (H 2 ) was high, and the GCV was high to moderate, so the selection of these traits could be effective as well. The findings of this genetic analysis for the traits tested suggest that the selection of productive genotypes in order to create new varieties that are stable in all environments for conventional and organic cultivation is possible. The estimates of genetic parameters of forage dry matter crude protein content characterized by high heritability and high genetic variability, and GCV indicates that selecting for better quality, as described by the protein content, is possible.

Correlations between Traits
In our study, many correlations between traits displayed statistically significant results. Statistically significant correlations are useful for indirect breeding and selection of traits that show low stability through more stable traits that promote adaptation [7]. Positive correlations were also reported for other traits in common vetch and peas by Greveniotis et al. [2,45,46]. Georgieva et al. [47] reported significant correlations for many traits in field peas. We found positive relationships between the fresh forage yield and dry matter yield, which were expected, but also the main stem length and crude protein content.
Correlation studies are very important in the genetic improvement of cultivars [48,49]. Singh et al. [50] reported significant correlations between seed yield per plant and harvest index, as well as the biological yield per plant, plant height, number of seeds per pod, number of primary branches per plant, number of pods per plant, and 100-seed weight. Days to maturity and 100-seed weight and number of pods per plant showed a weak negative correlation with the seed yield per plant. In our results, the most interesting correlation was between the fresh forage yield and the stable characteristic, 'stem length', for indirect breeding purposes [46]. Linearity was not satisfactory in many cases due to low correlation coefficients. Cacan et al. [1] reported interesting yield performances for forage pea lines. They also reported statistically significant correlations between many traits studied. Kosev and Mikić [51] also reported high and significant correlations between many traits in peas and, most of all, with significant linearity.
Sayar and Han [41] used GGE biplot analysis in two growing seasons. Their results showed that two lines and cultivar Kirazli were superior for fresh forage yield, dry matter yield, plant height, and days to 50% flowering. PC2 scores of these genotypes were found near zero, making them stable genotypes. Bocianowski et al. [16] reported that AMMI analysis managed to depict certain cultivars for certain environments regarding seed yield. This was an encouraging result for practical farming.
Sayar and Han [41], as well as Yihunie and Gesesse [52], reported that the GGE biplot could be used as a tool for the discrimination of pea genotypes according to their productivity and stability and the selection of the most suitable genotype for cultivation.
Uzun et al. [53] assessed the dry matter performance for peas used for their forage yield. He reported that semi-leafless lines had significantly better standing ability than leafed peas. The leaf type had no effect on lodging scores at the seed-harvesting stage. Yihunie and Gesesse [52] used a GGE-biplot of field peas genotypes and defined the ideal genotype. Among the twelve environments used, three environments were the best for discrimination, while one genotype was found to be the most stable, the highest yielding, and it was recommended for wider cultivation in Northwestern Ethiopia and similar areas. Georgieva et al. [25] also reported the specific adaptation of certain genotypes in field peas. In our study, Vermio proved to be a stable genotype for forage yield performance in low-input farming in the Trikala and Kalambaka area, while Pisso was the best in Florina and Giannitsa areas and low-input farming systems. The two pea lines exhibited stable performance in Giannitsa and Florina areas, especially in low-input conditions. Livioletta was also a stable genotype.

Exploratory Data Analysis of Peas
To provide a certain classification for the studied pea genotypes and cultivation systems, a heat map (Figure 8) was carried out. Cluster analysis was previously used for classification purposes for various genetic materials (e.g., maize, sweet cherry), sometimes in combination with principal components analysis [54][55][56]. The available data were divided into groups of increasing dissimilarity. Based on these results, the peas were divided into two distinct clusters (C1 and C2), each one having two subclusters (SC1, SC2, SC3, and SC4, respectively). Grouping for each subcluster revealed differences among pea cultivations. More specifically, SC1 consisted of low-input cultivated genotypes, which were characterized by a lower forage yield and dry matter yield, as well as lower ash content. SC2 included various other subgroups, mainly containing conventionally cultivated genotypes, which exhibited mostly low forage and dry matter yields. SC3 contained genotypes with higher yields, with two distinct subgroups, one cultivated conventionally and the other cultivated with low input. Lastly, SC4 included genotypes cultivated conventionally during the first growing season of experimentation and exhibiting higher main stems. There were no identified specific clusters based on locality.

Crop Establishment and Experimental Procedures
Four different locations (Table 7) were employed for the field experiments, two of them in Northern Greece and another two in Central Greece, divergent regarding soil

Crop Establishment and Experimental Procedures
Four different locations (Table 7) were employed for the field experiments, two of them in Northern Greece and another two in Central Greece, divergent regarding soil type, altitude, and environmental conditions. Five cultivars (common in Greek cultivations) of peas, namely, cv. Olympos, cv. Pisso, cv. Livioletta, cv. Vermio, and cv. Dodoni, and two lines (Zt1, Zt2) were used.
Two types of cultivation approaches were selected: low-input and conventional farming systems. The plots cultivated under the conventional farming system were fertilized before sowing so that 40 kg ha −1 Nitrogen and 80 kg ha −1 P 2 O 5 were added into the soil. For low-input cultivation, no fertilizers or other agrochemicals were applied during the experiment, while prior to the establishment of the experiment in 2008, the fields had been in a two-year rotation consisting of bread wheat/legumes without nutritional supplementation or other agrochemical inputs. Weeds were fully controlled by hand.
All genotypes were sown in early November 2008 and 2009 according to a strip-plot design, with the seven genotypes randomized within each plot and a plot size of 8.75 m 2 . Replications were four for each plot. Each plot consisted of seven (7) rows 5 m long, spaced at 25 cm, and the number of plants per plot was around 1000 according to the sowing rate. The number of seeds was 120 per m 2 , and the depth of sowing was 4 cm.

Climatic Conditions
Experimentation lasted two growing seasons (2008-2009 and 2009-2010), and the mean monthly air temperatures (maximum, minimum, mean) and rainfall data during the study period are provided in Table 8 for each experimental area based on daily records.

Measurements
For each plot, the number of days from the sowing date to 50% of the flowering time was recorded. Ten random plants of each plot were selected at the flowering time and measured from the ground level to the top point with a ruler (1 mm sensitivity) after extending the plants upward. The arithmetic mean of the measurements (in cm) was accepted as the 'main stem length' for each plot. The main stem thickness (mm) was calculated by measuring the stem diameter at the top, middle, and bottom of each stem selected. These traits also served as correlation variables.
The chloromass (fresh forage) obtained from each plot right after harvesting in full flowering time was weighed, and the value was converted to a hectare basis in order to calculate the 'fresh forage yield (kg ha −1 )'. After, fresh forage samples (0.5 kg), harvested from each experimental plot, were placed in a drying oven at 70 • C for 48 h, left to cool, and weighed; the dry matter yield was determined for each plot, followed by the calculation on a hectare basis in order to obtain the 'dry matter yield (kg ha −1 )'.
In order to analyze the forage dry matter crude protein content (%) and ash content % of the dry matter, the forage dry matter was ground to pass through a 1 mm sieve and subsequently mixed for the analysis. Ash content was determined according to AOAC Official Method 942.05 [57], while total nitrogen was determined using AOAC Official Method 988.05 [57], followed by total protein content estimation.

Data Analysis
The experimental design was a combined analysis of seven genotypes in four replications over four locations for two cultivation systems and two years of experimentation. The formal ANOVA should include the interaction of years × locations or years × genotypes, etc., which were not the aims of our study and made no practical sense. To overcome such a problem, we created a simpler ANOVA with one degree of interaction less, and it did not affect the precision of the analysis for the genotypes in different environments, so we conducted an ANOVA as follows. In order for the ANOVA to be more informative, the combination of each year and location was assigned as an environment in the general meaning since locations and years contribute to the effect of the environment on the genotypes. In this way, we have fewer interactions in the ANOVA table and do not affect the variance of genotypes and the G × E (genotype × environment) interaction, which is crucial for proceeding in the stability analysis. Stability estimations were based on the stability index (x/s) 2 , where x and s are the entry mean trait and the standard deviation, respectively [8,58]. Trait correlations were examined using the Pearson coefficient according to Steel et al. [59], and the significance of all the statistics was checked at p < 0.05 using SPSS ver. 25. Stability analysis was performed using the free version of PB Tools v1.4 (International Rice Research Institute, Laguna, Philippines) over locations and years for each characteristic and the statistical tools were the AMMI1 and (GGE) biplot analysis.
The mean squared values of genotypes, genotype × environment, error, and replicates were used to estimate the variance components following the methods suggested by McIntosh [60], which were used for the estimation of genetic parameters for the tested traits as follows: Heritability in a broad sense (H 2 ) was calculated according to Johnson et al. [42] and Hanson et al. [61]: The genotypic coefficient of variation (GCV) and phenotypic coefficient of variation (PCV) were calculated for all tested traits according to Singh and Chaudhary [62]: where σ 2 g , σ 2 p , σ 2 gxe , σ 2 re , and x are the genotypic variance, phenotypic variance, genotype × environment variance, residual variance (error), and overall mean for each tested trait, respectively.
The mathematical processing of the data was performed by hierarchical cluster analysis (HCA) using Ward's method. HCA analysis was performed using JMP 14 (SAS Institute Inc., Cary, NC, USA). The results from the cluster analysis are presented in a dendrogram.

Conclusions
Correlations among various characteristics showed significant positive relationships between the forage yield along with the dry matter yield and forage dry matter crude protein content. Indirect forage yield stability improvement may be performed by improving the main stem length, which generally showed high stability indices.
Comparisons between conventional and low-input farming systems generally revealed genotypes that displayed highly stable performance, even in low-input farming systems. Stability index data could also serve to estimate the kind of heritability of various traits, either quantitative or qualitative.
AMMI analysis, and consequently, a GGE biplot, along with ANOVA data, showed that there is a strong interaction between genotypes and environments, as well as the farming system (conventional or low-input). Therefore, the necessity arises to propose certain genotypes of field peas for specific areas and farming systems so as to obtain the most stable performance. The Vermio cultivar proved to be a stable genotype for forage yield performance in low-input farming in Trikala and Kalambaka areas, while Pisso was the best in Florina and Giannitsa areas and low-input farming. The two pea lines displayed stable performance in Giannitsa and Florina areas, especially in low-input conditions. The stable behavior of some genotypes in low-input farming systems could be valuable for farmers that raise livestock in mountainous areas.
The genetic parameters showed that all traits were of high heritability and moderate to high GCV, and the direct selection for fresh forage yield and dry matter yield was expected to be effective.
Limitations of this study are related to the differences in environmental data through time (across years). Low rainfall may significantly affect the genotype behavior across different environments.