Predicting Yield and Stability Analysis of Wheat under Different Crop Management Systems across Agro-Ecosystems in India

The objectives of the study were as follows: 1) to evaluate the GxExM for wheat genotypes; 2) to predict yield performance and identify high stable wheat genotypes in different management practices; and 3) to make genotype-specific management and high performing genotype recommendations within and across agro-ecological regions. A diverse set of twenty-one genotypes was evaluated over three years (2012, 2013 and 2014) under two levels of crop management practices (CT and ZT) across three agro-ecological regions (BR, MP and PB) of India in replicated trials. Data were analyzed with SASGxE and RGxE programs using SAS and R programming languages, respectively. Across and within a location(s), the pattern of GxExM and GxMxY interactions (respectively) among univariate and multivariate stability statistics, grouping of genotypes in divisive clusters and estimates (with a prediction interval) of genotype varied in management practice CT and ZT. Across locations, the genotypes “Munal” and “HD-2967” were the best performers and high stable in CT and ZT, respectively. Genotypes “HD-2824” and “DPW-621-50”, and “Munal” may serve as diverse parents for developing high quality, climate smart, locally adapted genotypes for BR in CT and ZT, respectively. Genotypes “HD-2932”, “BAZ” and “JW-3288”, and “GW-322” and “HD-2967” are suitable for developing locally adapted stress tolerant genotypes for MP in management practices CT and ZT, respectively. Relatively small GxM and GxExM interactions in PB preclude in making definitive conclusions.


Introduction
Worldwide, wheat is the most important staple food as it provides 55% of the carbohydrates and 20% of the food calories and protein consumed globally [1].
Demand for wheat is projected to continue to grow over the coming decades, particularly in the developing world to feed an increasing population, and with wheat being a preferred food, continuing to account for a substantial share of human energy needs in 2050 [2].In India, wheat is the second most important cereal crop after rice.During the past century, India has achieved remarkable progress in increasing wheat production and productivity, primarily due to genetically improved cultivars, expansion of sown areas and better cultural methods.Improved agronomic practices played an important role in enhancing the dependability and sustainability of yields [3] [4].Yield is a complex quantitative trait and greatly influenced by external environment, which results in scale or rank shift in its performance [5] [6].This relative shift of genotype performance from one environment (location x year combination) to another is known as genotype x environment interaction (GxE) [7] [8].Substantial studies have been conducted to identify high yielding and consistent performing wheat genotypes (also known as stable genotypes).However, most of the high stable genotypes are less predictable across different crop management practices since plant breeders often perform analysis of two-way data (genotype x site or GxE) for several consecutive years to detect stable genotypes without taking crop management practices into account.Previous studies on tillage practices, a crop management practice, suggest that optimization of tillage practice alters the dynamic and complex soil in physical, chemical and biological properties [9]- [15].Proponents of conventional tillage (CT) and zero tillage (ZT) agree on the values of one over the other for crop yield variation across locations and cropping systems [16]- [21].Understanding interactions among genotype, environment and crop management are important in planning key farm decisions including genotype selection, sustainable farm management and economic planning.In this paper, we propose to identify high yielding and high stable wheat genotypes across environments and over two-crop management practices: CT and ZT.
The impenetrable interaction of a crop bio-system with the soil, the atmosphere, and the environment that a plant lives in introduces challenges when making breeding decisions because it complicates the demonstration of superiority of any genotype across environments.Genotype x environment interaction may result in low correlation between phenotypic and genotypic values, thereby reducing progress from selection.This reduction leads to bias in the estimation of heritability and in the prediction of genetic advance [22] [23] [24].
Several strategies have been proposed to deal with GxE, and the most powerful strategy is to exploit GxE either to develop locally adapted material or to use GxE to better characterize the genotypes [25] [26].Genotype x environment interaction can be characterized using statistical methods ranging from univariate to multivariate models.The univariate models include regression slope, deviation from regression, environmental variance, and Kang's yield-stability; multi-variate models include genotype main effects plus genotype by environment interaction (GGE) biplot, and additive main effects and multiplicative interaction (AMMI) [27]- [33].Multivariate models could be graphically represented through biplots where genotypes and environments are plotted in a single graph.Recently, hierarchical Bayesian and mixed models were introduced to model heterogeneous variance among environments and different correlation structures among environments [26] [34] [35] [36] [37].Mixed models allow more flexibility to model unbalanced data using restricted maximum likelihood estimates (REML).
Each statistical method reflects different aspects of the GxE, and no single method adequately explains genotype performance across environments [6].Stability statistics are best used in combination with trait performance (mean or BLUP: Best Linear Unbiased Predictor is an estimate of random effect) and have successfully been used in plant breeding.
In this study, we were interested in understanding differential sensitivity of certain wheat genotypes to different agro-ecological environments in India under different crop management practices for enhancing the selection of superior and stable genotypes.The objectives of the study were as follows: 1) to evaluate the genotype x environment x management interaction (GxExM) for wheat genotypes; 2) to predict yield performance and identify high stable wheat genotypes in different management practices; and 3) to make genotype-specific management and high performing genotype recommendations within and across agro-ecological regions.

Germplasm, Location and Management
Twenty-one genotypes of wheat were evaluated across three years (2012, 2103, and 2014) and three locations ranging from western to eastern Indo-Gangetic plains under two different tillage management practices.Locations were chosen to represent the major wheat growing agro-ecological conditions for major wheat production regions in India: Pusa (25˚57'08''N; 85˚40'13''E), Bihar (BR), Jabalpur (23˚10'7.6''N;79˚55'55''E), Madhya Pradesh (MP), and Ludhiana (30˚59'28''N; 75˚44'11''E), Punjab (PB) (Supplemental Figure S1 and Supplemental Figure S2).The soils of the experimental plots at three locations, PB, MP and BR, were sandy loam, clay loam and silty loam, respectively.Twenty-one genotypes were semi-dwarf spring wheat chosen to represent new vs. old release; stress resistant vs. susceptible; rust resistant vs. susceptible; eastern vs. northern adapted; tall vs. short stature; early vs. late maturity; low vs. high yield; and varied 1000 kernel weight, protein content, spikelet ear −1 and seeds spike −1 .These genotypes are frequently used advanced breeding lines or officially released cultivars.Hereafter the word "genotype" is used to indicate cultigen, cultivar, variety or genotype.Genotypes were considered to be random and representative samples of a wide range of genetic and phenotypic diversity in the wheat germplasm population (Supplemental Table S1).Wheat genotypes were evaluated at two levels of crop management, a conventional tillage (CT) and zero tillage (ZT).
The CT involved two plowings with disc harrow, the first plowing used a spring-tyne cultivator (about 0.15 m depth) followed by planking at optimum soil moisture content; and the second plowing involved seeding of different varieties using "Limit-Plot Planter".In ZT, different varieties were directly drilled without any preparatory tillage operations and in the presence of standing stubbles of rice residues (~15 cm) using "Limit-Plot Planter" with inverted "T" tyne openers [38] [39].A standard seed rate (100 kg•ha −1 ), seeding depth (4 cm), fertilizer nutrient application (120 kg N, 60 kg P 2 O 5 and 60 kg K 2 O ha −1 ), water management (irrigation at standard critical stages) was adopted for both CT and ZT.These tillage and other management practices correspond to the crop management system commonly used in commercial wheat production in India [40] [41].

Trials and Data Description
At each location, individual trials were established as a two-factor strip-plot (split-block) design with ten genotypes, three replications and two-crop management practices [42].Out of the ten genotypes, six genotypes were tried at all three locations for three years, and the remaining four genotypes were locationspecific and not necessarily duplicates.These four location-specific genotypes were evaluated at one or two location(s) for either one or all three years.Thus, total distinct genotypes were twenty-one (Supplemental Table S1).Within the blocks, the genotypes were arranged in sub-blocks, and the two-crop managements in the other sub-blocks were arranged perpendicularly to the sub-blocks with the genotypes.The experimental unit (harvest plot) size was 1.6 m × 10 m.
In each location and year, the ten wheat genotypes were evaluated for grain yield (t•ha −1 ).Wheat grains were harvested using the guide of15% moisture in grain, yellowing of spikelets and hard dough stage of grains [11].Data were missing at random from one replication, one location and one year for genotypes "BAZ", "CSW-16", "CSW16", "CSW-18", and "DPW-621-50".We excluded three genotypes ("GW-273", "GW-366" and "HD-2687") that were tested for a time period shorter than three years and only in one location to provide a sufficient representative sample of years and location as random and fixed factors, respectively.As a result of this data preparation, we obtained a four-way genotype x management x location x year interaction (GxMxLxY) dataset for the grain yield of eighteen wheat genotypes (Supplemental Table S1).A set of six genotypes evaluated in all locations were analyzed for across-location statistics.Similarly, a set of ten genotypes (six genotypes from all locations + four unique genotypes from each location) were analyzed for individual location statistics (Supplemental Table S1).

Data Analysis and Statistical Methods
The unbalanced grain yield GxMxLxY data were analyzed for genotype, environment, management, and genotype x environment x management interactions with the SASGxE [6] [43] [44] and RGxE [45] [46] programs using SAS and R programming language, respectively, in two steps.During the first stage, we im-puted the missing values using the mice() function of the mice (multivariate imputation by chained equations) package [47] of R. Parameters maximum iteration 50, predictive mean matching (pmm) method and random generator seed value 500 were used in the mice() function to generate five imputed datasets.Then, in the second stage, the five imputed datasets were combined across trials and years to obtain a balanced GxMxLxY mean data for all statistical analyses.
Years and genotypes, and locations and managements were analyzed as random and fixed effects, respectively.Estimates and significance of random effects were computed using RGxE.The random effect model was fit using the lmer() function of lme4 (linear mixed effects models) package [48].The F ratio (= MS between /MS within , where MS is mean square or variance estimate) and significance of fixed effects were computed using mixed() function of afex package [49].
The mixed() function computes type III such as p-values using the default method via the Kenward-Roger approximation for degrees of freedom.Similarly, the significance of random effects was computed using a likelihood ratio test to attain p-values.Likelihood is the probability of the data given a model.The logic of the likelihood ratio test is to compare the likelihood of two models with each other using restricted maximum likelihood (REML) methodology.The model without the factor of interest (the null model) is compared with the model with the factor of interest (the full model) using the anova() function.It gives a Chi-Square value, the associated degrees of freedom and p-value.According to Wilk's theorem, the negative two times the log likelihood ratio of two models approaches a Chi-Square distribution with k degrees of freedom, where k is number of random effects tested [50].
Univariate stability statistics [regression slope (b i ), deviation from regression ( 2 d S ), Shukla's stability variance ( 2 i σ ), and Kang's yield-stability statistics (YS i )], and BLUP for genotypes were computed using RGxE.Regression slope (b i ) and deviation from regression ( 2 d S ); Shukla's stability variance ( 2 i σ ) and Kang's yield-stability statistics (YS i ); and best linear unbiased predictor (BLUP) for genotypes were computed using lm() function of R [51] stability.par()function of the agricolae package [52] and ranef() function of lme4 package [48], respectively.Tests for significance were derived using a t-test for each b i and an F test for each 2 d S for statistical differences from one and zero, respectively, at 0.05, 0.01 and 0.001 levels of probability.
SASGxE provided R code that is ready to use in R statistical software [51] for the analysis of multivariate stability statistics (GGE biplot) [44].GGE biplot analysis was computed using the "GGEBiplotGUI" package [53], with the support in the helper application "RStudio" [54] in R statistical software.GGE biplot analysis was used to visually assess the presence of genotype x environment interaction and rank genotype based on stability and mean in each management practice [32] [55].For each management practice, input data of GGE biplot analysis consisted of genotype x environment matrix (2 × 2) of mean values.
Similarly performing genotypes across and within locations and years were clustered using PROC VARCLUS of SAS v9.4 [56].The VARLCUS procedure used user-defined second eigenvalue cutoff and underlying algorithm called di-visive clustering to split a given set of genotypes into two groups.Eigenvalues are the coefficients of principal component analysis.The value 1 of the second eigenvalue is a common choice for cut off because it represents the average size of the eigenvalues.However, we have used the smaller value of the second eigenvalue as 0.7 to account for sampling variability [57].PROC VARCLUS identified clusters and computed 1 − R 2 ratio ( ) tifies a cluster of genotypes that are highly correlated among themselves and not highly correlated with genotypes in other clusters.The graphical representation of divisive clustering, 1 − R 2 ratio, forest plot of BLUP along with prediction interval, and the stability statistics summary of each genotype in different management practices across and within locations were computed using SAS PROC TEMPLATE in conjugation with PROC SGRENDER while utilizing the graphical template language (GTL) of SAS v9.4.

Variance Analysis
The pooled analysis revealed statistically significant fixed effects (L, M, MxL) for grain yield (Table 1).The variance estimates of the yield for the year (Y) and the Table 1.Pooled variance analysis for wheat yield (t ha −1 ) of 6 genotypes tested in 3 years and 3 locations over 2 management practices; and location variance of 10 genotypes tested for 3 years over 2 management practices.interaction between year and location (YxL) were different (p < 0.01) from zero and accounted for 39% of the total variation.Except for GxY and MxY, estimates of the random effects were significant (ranged from 0.17% -24% of the total variance estimate).Large estimate of environment (LxY) variance dominated the expression of genotype effect (G) and interaction between genotype and management (GxM) (Table 1).Within location analysis, both fixed and random effects were significant in BR, MP and PB.In BR and MP, Y, G and GxM caused most of the variation in yield performance of wheat genotype (22% and 39%, 35% and 24%, and 26% and 21% of the total variance, respectively).In contrast, PB had large Y effect (53% of total variance), medium G effect (12% of total variation) and small GxM effect (3% of total variation) (Table 1).

Polygon View of GGE Biplot
The "polygon" (or which-won-where) view of the GGE biplot is drawn by joining extreme genotypes and perpendicular lines (rays) passing from the polygon sides divides the biplot into sector (Figure 1 and Figure 2).A crossover GxE pattern exists if environments fall into multiple sectors and, thus, different genotypes win for an environment or set of environments in a sector.Conversely, no GxE pattern exists if all environments fall into a single sector and, thus, a single genotype had the highest yield in all environments.In either case, the winning genotype is vertex genotype.If no environment is present in a sector, then vertex genotype of that sector is considered to be a poor performer in all test environments.Genotypes within the polygon were less responsive to location than the vertex genotypes.Across locations, the polygon view of the GGE biplot explained 96% and 97% of the genotype and genotype x environment variation for the management practice CT and ZT, respectively (Figure 1: Panel A and B).
For management practice CT, environments were grouped into two sectors with different winning (vertex) genotypes.This confirms the existence of GxE for CT.
Conversely, ZT had all the locations grouped into one sector, which suggests that GxE for ZT did not exist.However, the GGL biplots for individual year confirm that the location grouping varied across years for ZT.Results of GGL are presented in Supplemental Figure S3

Genotype BLUPs
BLUPs are the estimates of random effects.Across locations, estimates of genotype (random effect) for wheat yield ranged from 5.03 to 5.55 t•ha −1 and 5.01 to 5.64 t•ha −1 for management practice CT and ZT, respectively (Table 2).For both management practices (CT and ZT), the highest yield was estimated for genotype "HD-2967" and the lowest yield was estimated for genotype "CSW-18".
Within individual location under both management practices, the predicted yield ranged from 4.67 to 5.87 t•ha −1 and 4.64 to 6.09 t•ha −1 ; 4.84 to 5.17 t•ha −1 and 4.77 to 5.14 t•ha −1 ; and 5.06 to 5.50 t•ha −1 and 5.03 to 5.57 t•ha −1 for location BR, MP and PB, and management practice CT and ZT, respectively (Table 2).
Under both management practices, location BR had genotypes "HD-2967" and "CSW-16" with the highest and lowest predicted yield, respectively.For location MP, the highest yield was predicted for genotype "GW-322" and "MPO-1215" in management practice CT and ZT, respectively.Genotypes "JW-3288" and "HD-2932" had the lowest predicted yield for management practice CT and ZT, respectively.Location PB had the highest predicted yield for "Munal" and "HD- 2967" in CT and ZT, respectively, and the lowest predicted yield for "CSW-16" in both management practices.to unity under both management practices, except for "DPW-621-50" in ZT (Supplemental Table S2).Conversely, the 2 types was significantly different (p < 0.05) from zero, except for "BAZ" in CT.

Univariate Stability Statistics
Within individual location under both management practices b i and 2 d S values were close to unity (p > 0.001) and zero (p > 0.001), respectively, for all ten tested genotypes.However, the exceptions were "Munal" and "JW-3288" in location MP under management practice CT for b i and 2 d S values significantly different from one and zero, respectively, and "CSW-16" in location BR under management practice ZT for 2 d S value significantly different zero.For management practice CT, genotypes "HD-2967"; "HD-2932", "HD-2967" and "JW-3288"; and "DBW-88" had a negative b i value in location BR, MP, and PB, respectively.Similarly, genotype "MPO-1215" had a negative b i value in location MP under management practice ZT.

Mean vs. Stability, View of GGE Biplot
The "average environment coordinate" (AEC) view based on genotype-focused singular value partitioning (SVP = 1) can be referred to as the "mean vs. stability" view [58] of GGE biplot.That view facilitates genotype comparisons based on mean performance and stability across environments within a mega-environment.The "mean vs. stability" view of GGE biplot explained 96% and 99% of genotypic and genotype x environment variation across locations under management practices CT and ZT, respectively (Figure 1: Panel C and D).The arrow shown on the AEC abscissa points in the direction of higher yield performance of genotypes and ranks the genotypes with respect to yield performance.Thus, genotype "HD-2967" (G14) had the highest yield and "CSW-18" (G03) had the lowest under both management practices (Figure 1: Panel C and D).The stability of each genotype was explored by its projection onto the AEC vertical axis.The most stable genotype was located almost on the AEC abscissa (horizontal axis) and had a near-zero projection onto the AEC (vertical axis).Thus, "HD-2967" (G14) and "DPW-621-50" (G05), and "Munal" (G19) and "BAZ" (G01) were the most stable in management practices CT and ZT, respectively (Figure 1: Panel C and D).The ideal genotype has both the high trait mean and stable performance [46].An ideal genotype is represented by a circle on the head of arrow on the AEC abscissa (horizontal axis) (Figure 1: Panel C and D).Across locations, genotype "HD-2967" (G14) was the best in both management practices (Figure 1: Panel C and D).These findings are consistent across years (Supplemental Figure S4).

Divisive Clusters and 1 − R 2 ratio
Across locations under both management practices, similarly performing genotypes were grouped into two clusters (Figure 4: Panel A and Panel B).The most representative and distinct genotype within the cluster has high correlation with its own cluster and low correlation with other clusters [59].Thus, an ideal representative genotype has a low 1 − R 2 ratio ( ) ue. Across all the tested locations "Munal" and "HD-2967" were the representative genotype of cluster 1 and cluster 2, respectively, in CT (Figure 4: Panel A). and "−" indicate stable and non-stable according to Kang stability statistics (YS i ), respectively."HYHS", "HYLS", "MYHS", "MYLS", "LYHS" and "LYLS" indicate "high yield and high stable", "high yield and low stable", "medium yield and high stable", "medium yield and low stable", "low yield and high stable" and "low yield and low stable", respectively.High, medium and low yield is with reference to the center of "mean vs. stability" view of GGE biplot (0,0), which represents the mean yield.
Similarly, "BAZ" was the representative genotype of cluster 1 in ZT (Figure 4: Panel B).The zero value of 1 − R 2 ratio is the result of the presence of a single entity in the cluster and, thus, correlation within its own cluster is 1 [ ] ( ) For location BR, similarly performing genotypes across years were grouped into three clusters under both management practices (Figure 5: Panel A and Panel B).In CT, "HD-2824" and "DPW-621-50" were representative genotypes in cluster 1 and cluster 3, respectively (Figure 5: Panel A).In ZT, "Munal", "CSW-16" and "HD-2967" were representative in respective clusters (Figure 5: Panel B).Similarly, for location MP, 5 and 2 clusters were identified in CT and ZT, respectively (Figure 5: Panel C and Panel D).In CT, "HD-2932", "BAZ" and "JW-3288" were representative in cluster 1, cluster 3 and cluster 4, respectively (Figure 5: Panel C).In ZT, "GW-322" and "Munal" were representative in cluster 1 and cluster 2, respectively (Figure 5: Panel D).For location PB under both management practices, genotypes were grouped into 1 cluster, except for "DBW-88" in CT (Figure 6: Panel A and B).

Discussion
Across locations, the F ratio of the fixed effects (L, M, MxL) was statistically significant.The mean yield of genotypes for the crop management practices (M), CT and ZT, was different (results are not presented here) and these results are in agreement with our previous findings [11]  estimates of the random effects were significant (Table 1).Similar results were reported in the study by [62] [63] [64] [65] [66].The estimates of random effects varied from 0.17% -24% of the total variance.The large estimates of Y and YxL suggest that the agro-ecological conditions of the test locations were extremely different and accounted for most of the yield variation.Variations in meteorological data across locations confirm the diversity found in agro-ecological conditions of the test locations (Supplemental Figure S1 and Supplemental Figure S2).The small contribution of G in the total variance estimate is due to the fact that genotypes evaluated in this study were advanced breeding lines, elite cultivars or high yielding genotypes.These findings are confirmed with the common parents and number of selfing (inbreeding) found in their pedigree (Supplemental Table S1).The significant variance components for GxM, GxMxL and GxMxLxY  σ or Si) at 0.05, 0.01 and 0.001 levels of probability, respectively.NS is significantly not different." § " indicates negative slope."+" and "−" indicate stable and non-stable according to Kang stability statistics (YS i ), respectively."HYHS", "HYLS", "MYHS", "MYLS", "LYHS" and "LYLS" indicate "high yield and high stable", "high yield and low stable", "medium yield and high stable", "medium yield and low stable", "low yield and high stable" and "low yield and low stable", respectively.High, medium and low yield is with reference to the center of "mean vs. stability" view of GGE biplot (0,0), which represents the mean yield.0.05, 0.01 and 0.001 levels of probability, respectively.NS is significantly not different." § " indicates negative slope."+" and "-" indicate stable and non-stable according to Kang stability statistics (YS i ), respectively."HYHS", "HYLS", "MYHS", "MYLS", "LYHS" and "LYLS" indicate "high yield and high stable", "high yield and low stable", "medium yield and high stable", "medium yield and low stable", "low yield and high stable" and "low yield and low stable", respectively.High, medium and low yield is with reference to the center of "mean vs. stability" view of GGE biplot (0,0), which represents the mean yield.
could be due to positive effects of ZT on soils specially when the organic matter content of the soils of the study sites is low and is one of the major reasons for deteriorating soil health in this region [67] [68].In CT, intensive and shallow tillage leads to significant negative impact on soil health and quality parameter, whereas ZT facilitate buildup in soil organic matter, improve soil biological activity and soil structure due to maintenance of soil aggregates, and reduced oxidation of soil organic matter.These findings are in agreement with [69] and [70].Results of the long-term research by [19] suggests that ZT significantly reduced bulk density and penetration resistance; increased organic carbon, saturated hydraulic conductivity, water stable aggregate, microbial biomass carbon and soil enzyme activities compared to CT.Additionally, ZT facilitates advance planting (ZT drills works well on relatively high soil moisture which is not the case for CT), regulates (buffers) soil as well as canopy temperature that helps in developing better root system (mass, density, length) which in turn results into better nutrient uptake and plant withstanding against high velocity winds at grain filling/maturity (reduce lodging) [71].
Based on the mean predictive yield, the six genotypes tested across locations were grouped into high, medium, and low predictive yielders.In both management practices (CT and ZT), high, medium and low predictive yielding genotypes were "HD-2967" and "Munal"; "BAZ", "CSW-16" and "DPW-621-50"; and "CSW-18", respectively.In CT, high and medium predictive yielding genotypes were grouped into two separate divisive clusters, which confirm the differential performance of genotypes across locations (Figure 4: Panel A).Conversely, in ZT, all high and medium predictive yielding genotypes grouped into 1 cluster (Figure 4: Panel B).The different clustering pattern in CT and ZT suggests the existence of interaction of soil health dynamics due to CT with different genetic composition, plant architecture, and stress tolerant and disease resistance levels of high and medium predictive yielding genotypes (Supplemental Table S1).In ZT, the low predictive yielding genotype "CSW-18" was highly dissimilar among all the tested genotypes and, thus, formed a unique cluster (Figure 4: Panel B).In contrast, the performance profile of "CSW-18" across locations in CT was on par and shared the same cluster (cluster 1) with the medium and high predictive yielding genotypes ("DPW-621-50" and "Munal", respectively) (Figure 4: Panel A).Genotypes grouped into the same cluster with the duplicate 1 − R 2 ratio value provided the same information (redundancy exist) and can be used interchangeably to reduce testing costs and improve efficiency of breeding programs.However, in both CT and ZT, distinct values of 1 − R 2 ratio in each cluster indicated the existence of some extent of dissimilarity among genotypes within the cluster (Figure 4: Panel A and Panel B).This finding can be corroborated with the unique genetic makeup of or distinct parents being used in each genotype (Supplemental Table S1).A definitive conclusion on the most representative and distinct genotype within the cluster must be based on the low value of 1 − R 2 ratio.Thus, "Munal" and "HD-2967" were the most representative and distinct genotype in cluster 1 and cluster 2, respectively, in CT (Figure 4: Panel A).Similarly, "HD-2967" was the most representative and distinct genotype in ZT (Figure 4: Panel A).In CT and ZT, the most representative and distinct genotype(s) had the highest predictive yield (Figure 4: Panel A and Panel B).These findings suggest that genetic composition of "HD-2967" has the ability to outperform regardless of alteration made in soil health due to tillage practice.Therefore, growers must select these high performing genotypes as their preferred varieties according to tillage practice and cropping system is followed in their region.Similarly, wheat breeders can use them as a diverse parent for future breeding purposes while considering tillage practice into account so that actual effect of yield can be exploited.
Based on multiple stability measures, wheat genotypes were classified into three categories.Category 1 genotypes had high predictive yield and high stabil-ity, and are widely adopted across varied environmental conditions.These genotypes had b i equal to unity, low or non-nonsignificant 2 d S and 2 i σ , YS i higher than mean YS i , and close to ideal genotype and near zero projection onto the AEC (vertical axis) on "mean vs. stability" view of GGE biplot.Genotypes "HD-2967" (G14) and "Munal" (G19) qualified in category 1 for management practices CT and ZT, respectively (Figure 4: Panel A and Panel B).Significant was detected for category 1 genotypes [6].This contrast could be because the errors associated with the slopes of genotypes are not statistically independent [72].Category 2 genotypes exhibited high predictive yield but low stability.
These genotypes had b i greater than unity, high or non-nonsignificant 2 d S and 2 i σ , and high projections onto the AEC (vertical axis) on "mean vs. stability" view of GGE biplot.Thus, these genotypes are suited for specific locations.Results suggest that genotypes "HD-2967" (G14) and "Munal" (G19) were sensitive to environmental change and had greater specificity of adaptability to high yielding environment under management practices ZT and CT, respectively (Figure 1: Panel C and Panel D).Plant breeders can use these genotypes to develop high performers for specific location and management practices.However, according to YS i , category 2 genotypes were better than average and should be considered as stable.Unlike category 1 genotypes, projection onto AEC on "mean vs. stability" view of GGE biplot, category 2 genotypes had relatively opposite projections.This finding confirms the differential performance of the same genotypes in different management practices.Change of stability (high to low or vice versa) of "HD-2967" (G14) and "Munal" (G19) as tillage practices flipped (CT to ZT or vice versa) suggested that the expression of genes governing yield were dependent on management mediated soil health and related processes.Yield of "HD-2967" (G14) was less influenced when soil had higher bulk density, poor water retention capacity, higher hydraulic conductivity, lower soil organic carbon and lower biological activities (soil properties found under CT, [39].Thus, growers can use "HD-2967" for consistent performance across locations and in cropping system where CT is practiced.In contrast, "Munal" (G19) is recommended for high yield and high adaptability at specific location where CT is practiced or soil had higher bulk density, poor water retention capacity, higher hydraulic conductivity, lower soil organic carbon and lower biological activities.However, growers have to be careful in using location specific genotypes as they had greater specificity of adaptability to high yielding environments.Meaning small change in favorable environment (CT type soils) will result into large fluctuation in genotype performance.Similarly, genotype "Munal" (G19) and "HD-2967" can be used for consistent performance across and within a specific location(s), respectively, and in cropping system where ZT is practiced.Likewise, breeders have scope to further exploit the genes responsible for differential performance of these high yielding and high stable genotype as soil physical, chemical and biological properties are changed due to optimization of tillage practice.Category 3 genotypes had low to medium predictive yield and low stability.These genotypes are suitable for traits (other than yield) such as disease or lodging resistance, where low values are desired in high yielding environments.Category 3 genotypes include "CAW-16" and "CSW-18" in both management practices.
In BR, high predictive yielding genotypes tended to perform better in ZT over CT, whereas low predictive yielding genotypes performed better in CT.In both CT and ZT, different genotypes constituted the three clusters, which suggest the differential performance of genotypes across years and management practices within a location.In CT, genotypes "CSW-16", "HD-2824", "K-307", "DPW-621-50" and "HD-2967" had high predictive yield and stability across years based on multiple stability measures.These genotypes had greater adaptability in this location under management practice CT.Conversely, genotypes "DPW-621-50", "HD-2733", and "Munal" had high predictive yield and stability in ZT.
Wheat researchers can use high predictive yielding, high stable and low 1 − R 2 ratio genotypes from different clusters as a contrast parent for future breeding purposes.These genotypes include "HD-2824" and "DPW-621-50", and "Munal" in CT and ZT, respectively.Similar to BR, in MP, high and low predictive yielding genotypes tended to perform better in ZT and CT, respectively.Relatively, CT resulted in more diverse performances of ten tested genotypes, thus forming 5 clusters.Genotypes "HD-2967", "Munal", "BAZ" and "MPO-1215" had medium to high predictive yield and stable across years in CT.The significant b i value of "Munal" suggested that its performance improved in good years.In ZT, medium predictive yielding genotypes were stable.These genotypes include "BAZ" and "DPW-621-50".
High predictive yielding genotypes had large prediction intervals and non-nonsignificant 2 i σ .In CT, genotypes "HD-2932", "BAZ" and "JW-3288" with the lowest 1 − R 2 ratio had high stability.In contrast, in ZT, "GW-322" and "HD-2967" had the lowest 1 − R 2 ratio but with low stability.Pedigree of "GW-322" confirmed the presence of stem rust (Sr2, Sr11+), leaf rust (Lr1, Lr10, Lr13, Lr13+), yellow rust (Yr2+) and glutenin (Glu-A1b, Glu-B1b, Glu-D1a) genes (Supplemental Table S1).Glutenin (a type of glutelin) is the major protein of wheat flour.Similarly, "BAZ" and "HD-2932" have adult plant stem rust (Sr2+) and leaf rust (Lr13+) resistant genes, respectively.Additionally, previous studies have identified that "BAZ" and "GW-322" have capability to withstand heat and drought stress.Therefore, researchers can introgress genes from these genotypes into elite inbreds to make better hybrids with high stability.In PB, a small percentage of the total variation from GxM and GxMxY interactions (p < 0.01) suggests that the grain yield variation from these interactions is the least important.Furthermore, estimates of random effects of genotypes revealed mixed response of high and low predictive yielding genotypes in both management practices.Similarly, a divisive cluster did not identify dissimilarity among genotypes, except for genotype "DBW-88" in CT.Likewise, genotypes "DPW-621-50", "HD-   2. Significance value of regression coefficient, deviation from regression, Shukla's stability variance, and Kang's stability statistics for wheat yield (t•ha −1 ) of 6 genotypes tested in 3 years and 3 locations over 2 management practices; and for location BR, MP, and PB for wheat yield (t•ha −1 ) of 10 (6 same + 4 different) wheat genotypes tested for 3 years over 2 management practices.

Figure 1 .
Figure 1.The polygon (which-won-where) (Panel A and Panel B) and mean vs. stability (Panel C and Panel D) view of genotype main effects plus genotype x environment interaction effect (GGE) biplot of yield of 6 wheat genotypes tested in 3 years, 3 locations over 2 management practices.The biplots were based on Scaling = 0, Centering = 0, and SVP = 2. Key to the labels of genotype and location is presented in abbreviation section.

Figure 2 .
Figure 2. The polygon (which-won-where) view of genotype main effects plus genotype x year interaction effect (GGY) biplot of location BR (Panel A and Panel B), MP (Panel C and Panel D) and PB (Panel E and Panel F) for yield of 10 (6 same + 4 different) wheat genotypes tested in 3 years over 2 management practices.The biplots were based on Scaling = 0, Centering = 0, and SVP = 2. Key to the labels of genotype and location is presented in abbreviation section.
According toEberhart and Russell (1966), a regression coefficient (b i ) approximating unity, along with deviation from regression ( 2 d S ) near zero, indicates stability.Across locations, the b i value for all six genotypes was close (p > 0.001)

Figure 3 .
Figure 3.The mean vs. stability view of genotype main effects plus genotype x year interaction effect (GGY) biplot of location BR (Panel A and Panel B), MP (Panel C and Panel D) and PB (Panel E and Panel F) for yield of 10 (6 same + 4 different) wheat genotypes tested in 3 years over 2 management practices.The biplots were based on Scaling = 0, Centering = 0, and SVP = 2. Key to the labels of genotype and location is presented in abbreviation section.

Figure 4 . 2 dS
Figure 4. Divisive cluster, 1 − R 2 ratio, forest plot of BLUP along with 95% prediction interval, and stability statistics summary of yield of 6 wheat genotypes tested in 3 years, 3 locations over 2 management practices.† *, ** and *** significantly different from unity for the regression coefficients or slope (b i ) and from zero for the deviation from regression ( 2 d S ) and Shukla's stability variance ( 2 iσ or Si) at 0.05, 0.01 and 0.001 levels of probability, respectively.NS is significantly not different." § " indicates negative slope."+" and "−" indicate stable and non-stable according to Kang stability statistics (YS i ), respectively."HYHS", "HYLS", "MYHS", "MYLS", "LYHS" and "LYLS" indicate "high yield and high stable", "high yield and low stable", "medium yield and high stable", "medium yield and low stable", "low yield and high stable" and "low yield and low stable", respectively.High, medium and low yield is with reference to the center of "mean vs. stability" view of GGE biplot (0,0), which represents the mean yield.
[12] [18] [19][20] [21]  [60][61].The significant MxL suggests a variable response of genotype yield to the increased intensity of tillage practices across the test locations.Except for GxY and MxY, led to a different ranking of genotypes across environments under CT and ZT, justifying the development of a stable genotype that performs well over environments in different management practices.The ideal genotype should have a high mean and high stability.Across locations, the mean predicted (mean of BLUP) yield of ZT is more than CT (represented by the vertical line in Figure4: Panel A and Panel B).A similar trend of enhanced performance of ZT over CT across nine environments (three locations × three years) was noticed (Supplemental FigureS5).This

Figure 5 .
Figure 5. Divisive cluster, 1 − R 2 ratio, forest plot of BLUP along with 95% prediction interval, and stability statistics summary of locations BR (Panel A and Panel B) and MP (Panel C and Panel D) for yield of 10 (6 same + 4 different) wheat genotypes tested for 3 years over 2 management practices.† *, ** and *** significantly different from unity for the regression coefficients or slope (b i ) and from zero for the deviation from regression ( 2 d S ) and Shukla's stability variance(2i

Figure 6 .
Figure 6.Divisive cluster, 1 − R 2 ratio, forest plot of BLUP along with 95% prediction interval, and stability statistics summary of locations PB for yield of 10 (6 same + 4 different) wheat genotypes tested for 3 years over 2 management practices.† *, ** and *** significantly different from unity for the regression coefficients or slope (b i ) and from zero for the deviation from regression ( 2 d S ) and Shukla's stability variance ( 2 i σ or Si) at

2 dS
of both genotypes contrasts to the previous study where non-significant 2 d S 2967" and "PBW-550" had high predictive yield and high stability in both CT and ZT.The mixed response of high and low predictive yielding genotype could G21 or 21 = PBW-550 HYHS = High yield and high stable HYLS = High yield and low stable L = Location LYHS = Low yield and high stable LYLS = Low yield and low stable M = Management MET = Multi-environment trial MP = Jabalpur, Madhya Pradesh MYHS = Medium yield and high stable MYLS = Medium yield and low stable PB = Ludhiana, Punjab PC = Principal component RGxE = R language program for the analysis of genotype stability and location value REML = Restricted maximum likelihood SASGxE = SAS program for the analysis of genotype stability and location value SVP = Singular value partitioning Y = Year

2 iσ
or Si = Shukla's stability variance b i = Linear regression coefficient 2 d S = Deviation from regression YS i = Kang's stability statistic Supplementary Table