Multi-Trait Selection of Quinoa Ideotypes at Different Levels of Cutting and Spacing

: Climate change has affected the food supply chain and raised serious food concerns for humans and animals worldwide. The present investigation aimed to assess the effect of environmental factors along with three different levels of cutting (i


Introduction
Improvements in agricultural productivity are directly associated with upgrading universal living standards and sustaining crop production under rising challenges to the environment and the agriculture sector [1].Crop production has increased as a result of investments in agricultural research sectors, but this boost in productivity has not been distributed evenly across the world [2], and there are indications that it is decelerating in particular regions [3].Simultaneously, several human activities such as industrialization and globalization have raised the global temperature to 2 • C during the last half-century.This rise in temperature ultimately affects the global climate pattern, which is directly affecting crop yield [4]. Climate change has now risen as an additional burden, particularly in developing countries like Pakistan, in accomplishing goals related to food security.Several studies have reported the negative impacts of climate change on several plant species including their seed yield, biomass production, and quality as fodder for animals [5,6].Moreover, this change in climate patterns is also adversely affecting soil microbiota by suppressing and reducing the population of beneficial microbes.These beneficial microbes increase soil productivity (increasing the concentration of nutrients that are necessary for plant growth and development) and suppress soil-borne pathogens.Keeping in view the demands of fodder crops as feed for animal consumption, it is important to introduce some new crop species as fodder that have the potential to meet both the qualitative and quantitative criteria of fodder production.Fodder production greatly varies from country to country and region to region depending on weather conditions and other factors such as soil nutrients and fodder crop species.Ensuring high forage production by introducing exotic germplasm under new climate conditions is very crucial for sustainable agriculture.Crop breeding and genetics and agronomic practices together play a major role in ensuring maximum fodder production worldwide.Overall, ensuring long-term food security and resource conservation, minimizing environmental impacts, and fostering resilience in the face of changing climatic circumstances all depend on the combination of crop genetics and agronomic techniques in sustainable agricultural systems [7,8].
Crop selection based on targeted traits (such as yield-related and quality traits in fodder) is another challenging factor for plant breeders.Experienced breeders often target several plant traits, and if introduced into a new genotype, they would result in high crop performance.Hence, this new genotype is treated as an ideotype, a concept introduced by Donald [9].These multi-trait selection indexes are often very challenging and difficult, using statistical tools with several limitations and errors.For example, The Smith-Hazel (SH) index is a widely used phenotypic selection index based on phenotypic and genotypic covariance [10,11].The SH index requires inverting a phenotypic covariance matrix (i.e., multicollinearity), which often appears when assessing several traits and may thus lead to biased index coefficients and poorly conditioned matrices.Furthermore, in addition to multicollinearity, breeders also face difficulties in explaining the real economic value of particular traits by converting them into realistic economic weightings.The MGIDI was therefore introduced to overcome these restrictions [12].The MGIDI is designed to select genotypes based on multi-trait information, overcoming the abovementioned limitations.So far, the MGIDI has been successfully used in selecting superior genotypes [13][14][15][16].
Quinoa (Chenopodium quinoa Willd.) is well known for its highly nutritious profile and is currently grown in several regions of the world.It was first cultivated by pre-Columbian cultures about 7000 years ago and called the "mother grain" of the Incan Empire [17].Quinoa has already developed resistance to different abiotic stresses on Andean Altiplano (>3500 m above sea level) [18].Global demand for quinoa has risen rapidly due to its highly nutritious seeds, which are free of gluten, have a low glycemic index [19], and contain an excellent amount of minerals (zinc, sodium, potassium, phosphorus, iron, manganese, magnesium, copper, and calcium), vitamins (B, E, and C), carbohydrates, lipids, fiber, and essential amino acids [20].Quinoa is also known to possess more proteins than other cereals like wheat, rice, maize, and barley.In addition to quinoa's nutritious value, it also has the potential to overcome challenges associated with food security by cultivating it on marginal lands that are not fit for other paramount crops.This potential of quinoa was acknowledged in 2013 when the United Nations declared the year as "Quinoa International Year".Till now, only three crops (rice in 2004 and potato in 2008) have received this award.Despite the agronomic potential of quinoa, there are not many active breeding initiatives for it, making it an underutilized crop [21,22].Plant breeders across the world have recognized and introduced quinoa as a multipurpose crop from human food to animal feed.All breeding programs of quinoa are relatively focused on developing high-yielding seed quinoa genotypes despite the complex interactions between environmental factors and other agronomic traits.The biomass production of quinoa is also as important as seed production because fresh leaves and stems can also be used as forage for animal consumption.Different investigations have reported that quinoa leaves contain more nutrients than quinoa seeds, and therefore, quinoa leaves and stems (together, green biomass) can be an important source of fodder for animals.In addition to quinoa, oat (Avena sativa L.) is another important fodder crop that is grown worldwide to feed animals [23].It is a multi-cut crop and well placed in farming systems for quantity and quality of fodder supply during the shortage period (December to April) to fill up the gap [24].In Pakistan, oat is mainly cultivated for animal feeding in both irrigated and rainfed environments at high altitudes (1000-2300 m) [25].However, oat biomass production is low at higher altitudes, i.e., above 2300 m, and still cannot fulfill the fodder feed gap.Moreover, climate change has adversely affected oat biomass worldwide, but, on the other hand, quinoa has naturally developed a resistance to a wide range of abiotic and biotic stresses [18,26,27].
The fresh biomass of quinoa, especially harvested chuff and leaves, are highly preferred by animals like buffalo, cow, ox, goat, and sheep.So far, very few experiments have been conducted on different ways to use quinoa in animal diets, such as for silage and forage.Quinoa dry matter yield is acceptable for animal consumption because it is easily digestible and contains an excellent amount of crude fat (CF), crude protein (CP), neutral detergent fiber (NDF), and low acid detergent fiber (ADF) content, which make quinoa an excellent quality forage [28,29].Forage quality is directly associated with low levels of anti-nutritional quality traits (e.g., ADF and NDF), which increase its digestibility.ADF and NDF are actually cell wall components; NDF is composed of hemicelluloses, and ADF consists of lignin and cellulose [30,31].Furthermore, some other factors such as management practices, genotype, and plant growth stages also determine forage quality.The concentrations of CP, CF, ADF, and NDF in leaves and stems also depend on different plant-to-plant spacing levels along with other agronomic practices [32][33][34].Different studies have recommended various plant-to-plant spacing levels and have examined the effect of spacing levels on quinoa biomass production [35,36].
Quinoa has recently been introduced as a new crop in Pakistan.Therefore, this investigation aimed to assess (i) the effect of the local climatic conditions of Faisalabad on the growth and development of quinoa; (ii) the effect of different levels of cutting and spacing on morphological and quality traits of quinoa; (iii) the selection of suitable genotypes based on multiple recorded parameters using a selection index (SI); and (iv) the comparison of selected quinoa genotypes with oat accessions.

Experimental Site, Design, and Agronomic Managements
The current experiment was conducted for two consecutive years (i.e., 2018 to 2019 and 2019 to 2020) under the local climatic conditions of Faisalabad (31 • 25 7.3740 N and 73 • 4 44.7924E; 186 m above sea level) at the Agronomy research area of the University of Agriculture, Faisalabad, Pakistan.The climate of Faisalabad is semi-arid and classified as very hot and humid in summer with cool and dry winters.Faisalabad receives an average annual rainfall of 375 mm and almost half of this falls during the monsoon season (i.e., July to August).A total of 10 accessions (R11, R126, R15, R24, R3, R30, R4, R9, and UAF7 originated from New Mexico, United States, and R45 originated from Chile) of Quinoa seeds (Figure 1) were collected from the department of crop physiology, University of Agriculture Faisalabad.Before sowing, the land was irrigated and prepared by two ploughings to a depth of 15 cm followed by planking to crush hard clods to smoothen the soil surface and preserve moisture in order to ensure maximum germination/emergence.The quinoa accessions for two growing seasons (i.e., year 1 and year 2) were sown on November 15th and November 20th, respectively.The seeds were sown manually using a hand dibbler (2-3 cm depth) with a 30 cm row-to-row distance.Plant-to-plant distance/withinrow spacing was kept different, i.e., 21, 23, and 26 cm, and the data for different traits were recorded at three different cutting levels (i.e., cutting 1, cutting 2, and cutting 3) at different growth stages (vegetative growth stage, budding stage, and flowering stage) (Figure 2).The experimental layout was a split-split plot design under the randomized complete block (RCBD) arrangements with three replications, where levels of cutting were assigned in the main plot, levels of spacing were assigned in the subplot, and genotypes were assigned in the sub-subplots (Figure 1).The soil at the experimental site was a sandy loam and contained potassium (169 mg kg −1 ), nitrogen (0.4 g kg −1 ), phosphorus (8.1 mg kg −1 ), a pH of 7.82, and organic matter (0.74%).All agricultural practices such as irrigation (2 irrigations in first year and 1 in the second year due to heavy rains), fertilizer application (N:P:K: 75:50:50 kg ha −1 ), and weeding were conducted rigorously in a timely fashion.During the growing seasons, no major weather hazards such as insect/pest attacks occurred.) represent plant residues after harvesting.Note: For leaves' quality traits, the prefix "L" was used with each trait, in other words, LCP, LCF, LADF, and LNDF.Similarly, for stems, the prefix "S" was used with each trait, in other words, SCP, SCF, SADF and SNDF.

Measurements of Morphological and Quality Parameters
Data for different quality and morphological traits were recorded in three different levels of cutting (harvestings) at different growth stages.The first cutting was performed after the plants reached a height of 30 cm followed by cutting 2 and cutting 3-with 15-day intervals between each cutting.Five plants from each plot were selected randomly, and data were recorded for plant height (PH cm), leaf area (LA cm 2 ), number of branches/plant (Br), number of leaves/plant (NoL), intermodal distance (ID cm), dry weight (DW g plant −1 ), and fresh weight (FW g plant −1 ).Fresh weight was immediately measured after each cutting using a weighing balance (Mhand, DJ-V-B).To measure dry weight, plants were transferred into paper bags, dried in the sun for 2 days, and then transferred to the oven for further drying at 65 ± 2 °C until a constant weight was obtained.After dry weight measurements, stem and leaves were separated and ground in a grinder machine, and data for fodder quality traits such as crude fat (CF %), crude protein (CP %), neutral detergent fiber (NDF %), acid detergent fiber (ADF %) [37], and moisture (Mois %) were recorded separately for both stems and leaves using a near-infrared (NIR) instrument (AgriNIR™) [38].  ) represent plant residues after harvesting.Note: For leaves' quality traits, the prefix "L" was used with each trait, in other words, LCP, LCF, LADF, and LNDF.Similarly, for stems, the prefix "S" was used with each trait, in other words, SCP, SCF, SADF and SNDF.

Measurements of Morphological and Quality Parameters
Data for different quality and morphological traits were recorded in three different levels of cutting (harvestings) at different growth stages.The first cutting was performed after the plants reached a height of 30 cm followed by cutting 2 and cutting 3-with 15-day intervals between each cutting.Five plants from each plot were selected randomly, and data were recorded for plant height (PH cm), leaf area (LA cm 2 ), number of branches/plant (Br), number of leaves/plant (NoL), intermodal distance (ID cm), dry weight (DW g plant −1 ), and fresh weight (FW g plant −1 ).Fresh weight was immediately measured after each cutting using a weighing balance (Mhand, DJ-V-B).To measure dry weight, plants were transferred into paper bags, dried in the sun for 2 days, and then transferred to the oven for further drying at 65 ± 2 °C until a constant weight was obtained.After dry weight measurements, stem and leaves were separated and ground in a grinder machine, and data for fodder quality traits such as crude fat (CF %), crude protein (CP %), neutral detergent fiber (NDF %), acid detergent fiber (ADF %) [37], and moisture (Mois %) were recorded separately for both stems and leaves using a near-infrared (NIR) instrument (AgriNIR™) [38].
) represent plant residues after harvesting.Note: For leaves' quality traits, the prefix "L" was used with each trait, in other words, LCP, LCF, LADF, and LNDF.Similarly, for stems, the prefix "S" was used with each trait, in other words, SCP, SCF, SADF and SNDF.

Measurements of Morphological and Quality Parameters
Data for different quality and morphological traits were recorded in three different levels of cutting (harvestings) at different growth stages.The first cutting was performed after the plants reached a height of 30 cm followed by cutting 2 and cutting 3-with 15-day intervals between each cutting.Five plants from each plot were selected randomly, and data were recorded for plant height (PH cm), leaf area (LA cm 2 ), number of branches/plant (Br), number of leaves/plant (NoL), intermodal distance (ID cm), dry weight (DW g plant −1 ), and fresh weight (FW g plant −1 ).Fresh weight was immediately measured after each cutting using a weighing balance (Mhand, DJ-V-B).To measure dry weight, plants were transferred into paper bags, dried in the sun for 2 days, and then transferred to the oven for further drying at 65 ± 2 • C until a constant weight was obtained.After dry weight measurements, stem and leaves were separated and ground in a grinder machine, and data for fodder quality traits such as crude fat (CF %), crude protein (CP %), neutral detergent fiber (NDF %), acid detergent fiber (ADF %) [37], and moisture (Mois %) were recorded separately for both stems and leaves using a near-infrared (NIR) instrument (AgriNIR™) [38].

Statistical Analysis
To test the normality and homogeneity of variance, Shapiro-Wilk and Bartlett's tests were used, respectively.Data were also checked for any outliers and multicollinearity among the variables before performing ANOVA.Data of all traits were subjected to splitsplit plot ANOVA to measure the main and interaction effects in each year separately (Supplementary Table S1).To analyze the overall effect and interactions of cutting (cut), spacing (space), genotype (gen), and year for all traits, a 4-way multivariate analysis of variance (MANOVA) was used to obtain a cumulative significance value based on Pillai's trace test [39] using the model mentioned below: where Y represents a response variable matrix of n × p, n represents total plot numbers, i.e., the combinations of different levels for cutting (c = 1, 2, 3), spacing (s = 1, 2, 3), genotypes (g = 1, 2, 10), and years (y = 1, 2), and p represents the total number of traits.Furthermore, X, b, and e represent n × m model matrix (as m = c × s × g), m × p coefficients model matrix, and n × p residuals model matrix, respectively.The estimated means for all treatments and traits were organized into rows and columns in a two-way table using the projected values from significant terms.Tukey's multiple comparison test (p < 0.05) was undertaken to measure differences using means of genotypes and treatments.In order to identify the best combination of cut × space treatment to obtain the maximum benefits and the best performing genotypes, the MGIDI [40] was used with the help of the metan package: where MGIDIi represents multi-trait genotype-ideotype distance index for the ith row/ genotype/treatment; ij represents the ith score row/genotype/treatment in the j th factor (i = 1; 2;. ..; g; j =1; 2;. ..; f), with g and f being the number of rows/genotypes/treatments and factors, respectively; and j represents the ideotype j th score.The combinations that depict the lowest MGIDI will be adjacent to the ideotype and, thus, exhibit preferred values for all measured traits.Selection intensity (SI) of ~20% for main effects and 10% for interactions was considered for selection differential for all traits.Based on the ideotype concept, all the traits were assigned a value between 0-100, where 0 represents a leastvalued trait, and 100 represents a most-desired trait.In this study, LADF, LNDF, SADF, and SNDF were selected as traits in which a decreased value was more desirable, while for the rest of the traits, an increased value was desirable.Moreover, weights were also added to traits to provide differential importance to each trait.The morphological traits were given a value of 0.5, whereas the quality parameters were assigned a value of 1 to indicate that, during selection of genotypes/cutting/spacing, these traits will be considered first, followed by the morphological traits.
All the data manipulation and analyses were carried out using R software (version 4.2.2).Correlation matrix between different morphological and quality traits along with all treatment combinations were generated and plotted using qgraph package, bar-graphs, and a climograph for weather data (Figure 3), and other figures were drawn using the ggplot2, dplyr, patchwork, plotly, and metan packages in R software.

Results
MANOVA results revealed that all the main effects (genotypes, cut, space, and years) showed significant (p < 0.05 & 0.01) differences.All the possible two-way interactions studied between the main effects were found significant except space × year, indicating that the response of the genotypes was influenced by cutting, spacing, and year effects.As for the three-way interactions, only cut × space × gen and cut × gen × year were found significant, while the four-way interaction was found non-significant (Supplementary Table S2).

Climatic Conditions
The maximum and minimum temperature during the first growing season was recorded 19 and 11 °C in November and January, respectively.In the second growing season, the maximum (19 °C) and minimum (10 °C) temperature was recorded in November and January, respectively (Figure 3).The average rainfall per year (last 55 years) in Faisalabad in November, December, January, and February has recorded 2.9, 3.1, 16.9, and 26.3 mm.During the first growing season, 0, 3.3, 20.9, and 49.4 mm of average rainfall was recorded in November, December, January, and February, respectively.While in the second growing season, the average rainfall in all months (except February) remained much higher than in the first growing season.During the second growing season, 16.1, 14.6, 60.3, and 15.7 mm of rainfall was recorded in November, December, January, and February, respectively.Moreover, the same pattern was observed for relative humidity, where, during the first growing season, it was 23% in November, 25% in December, 39%

Results
MANOVA results revealed that all the main effects (genotypes, cut, space, and years) showed significant (p < 0.05 & 0.01) differences.All the possible two-way interactions studied between the main effects were found significant except space × year, indicating that the response of the genotypes was influenced by cutting, spacing, and year effects.As for the three-way interactions, only cut × space × gen and cut × gen × year were found significant, while the four-way interaction was found non-significant (Supplementary Table S2).

Climatic Conditions
The maximum and minimum temperature during the first growing season was recorded 19 and 11 • C in November and January, respectively.In the second growing season, the maximum (19 • C) and minimum (10 • C) temperature was recorded in November and January, respectively (Figure 3).The average rainfall per year (last 55 years) in Faisalabad in November, December, January, and February has recorded 2.9, 3.1, 16.9, and 26.3 mm.During the first growing season, 0, 3.3, 20.9, and 49.4 mm of average rainfall was recorded in November, December, January, and February, respectively.While in the second growing season, the average rainfall in all months (except February) remained much higher than in the first growing season.During the second growing season, 16.1, 14.6, 60.3, and 15.7 mm of rainfall was recorded in November, December, January, and February, respectively.Moreover, the same pattern was observed for relative humidity, where, during the first growing season, it was 23% in November, 25% in December, 39% in January, and 52% in February, while in the second growing season, humidity remained at 37% in November, 30% in December, 54% in January, and 39% in February.

Morphological and Forage Quality Traits of Quinoa Genotypes
Bar graphs (mean values) of quality and morphological traits are presented along with Tukey's mean comparison test results (Supplementary Figure S1).The year effect for most morphological traits was significant.PH, NoL, Br, FW, and DW showed significantly higher values at cutting level 3 for all the genotypes than at C1 and C2 levels, while for LA and ID, no significant effects of cutting, spacing, and year were recorded.S2 (23 cm) was more effective for morphological traits than S1 (21 cm) and S3 (26 cm) (Supplementary Figure S1).
Among leaf quality parameters, all the quality traits were stable across years except LNDF and LCF.The concentrations of all quality parameters increased with an increase in growing time, and hence, concentrations at cutting level 1 were lower than at cutting levels 2 and 3. Furthermore, S2 was more effective for all quality parameters (except LNDF) compared with S1 and S3 (Supplementary Figure S2).Stem quality traits also showed a similar pattern, where significant differences were observed in spacing and year for all traits except SADF (Supplementary Figure S3).However, SCP increased with an increase in growth time, i.e., from cutting 1 to cutting 2 but decreased at cutting 3 (Supplementary Figure S2A).SCF, SADF, SNDF, and SMois slightly increased at cutting level 3 compared to cutting level 2 (Supplementary Figure S2B-E).

MGIDI Index Selection
Keeping in view the significance and importance of factors, the MGIDI was calculated for (i) cut × space to identify the best possible combination for quinoa and (ii) genotype selection and to identify the best genotypes that performed ideally at all levels of cutting and spacing in both years.Based on the studied pattern, traits were grouped into the same factor.

Selection of Best Application combination (Cut × Space) to Obtain the Maximum Desired Gains
The MGIDI was computed to rank the treatment combination (cut × space) based on the desired values of each trait.The factor analysis grouped all the studied traits into single factors (FA1) based on their similar pattern as follows (Table 1).Combination ranking according to the MGIDI showed that cutting level 3 combined with spacing level 2 (23 cm) was more effective for achieving the desired results, followed by C3-S3 (Supplementary Figure S1).The computed MGIDI displayed the desired selection differentials (SDs) for 13 out of 17 studied traits, which means that our goal for the selected 13 traits was achieved (Figure 4).The SDs quantify the population's mean trait value alterations between preselection and post-selection.The SD for the traits with less-desired values ranged from 6.73% (SADF) to 11.7% (LADF).The SD for the traits with more highly desired values ranged from 2.19% (SCP) to 72.6% (DW).Overall, the computed MGIDI index displayed lower total gain, i.e., 28.8% and 8.57% for the traits that tended to increase and decrease, respectively.The positive SD observed for productive traits was mainly due to C3-S2 and C3-S3, indicating that a combination of C3 with S2 would be more fruitful for obtaining the maximum desired gain for the quinoa genotypes.

Selection of Best Performing Genotypes
In the significant factor analysis of genotype, 17 traits were grouped into five factors (FA1-5) (Table 2).The MGIDI provided a total SD of 5.31% for traits that tended to increase and −3.43% for traits that tended to decrease.The SDs calculated with these genotypes were desired for a total of 15 traits out of 17.The SD ranged from −8.83% (LADF) to 1.70% (SNDF) for traits desired with low values and from −1.61% (ID) to 15.9% (FW) for traits desired with high values.Overall, all the genotypes performed better in the second growing season, with R9 and R3 ranked as the most suitable genotypes.Both genotypes showed maximum contribution in FA1 and FA4, mainly relying on quality-related traits (Figure 5).The selected genotypes were then used to compute the SD.The strengths and weaknesses plot and PCA (Figure 5b,c) showed that the selected genotype's strengths were mostly associated with FA1 and FA4, i.e., higher values for most of the morphological and quality traits compared to the other genotypes.Factors that contributed the most and least were drawn at the center and edges of the plot, respectively.

Correlation Analysis
Forage biomass and quality depend on the relationship among the traits that directly contribute to quality and yield.In this study, correlation analysis revealed moderate (p < 0.05) and high (p < 0.01) significant positive correlations among the studied traits in both years under all levels of cutting and spacing, but some traits also showed negative correlations (Supplementary Tables S3-S8).
An increase in biomass resulted in the degradation of anti-nutritional quality traits (ADF and NDF).Biomass-related traits like PH, NoL, LA, Br, FW, and DW showed highly and moderately significant negative correlations with LADF and LNDF, ranging from r = −0.501 to −0.935 in both years under all levels of cutting and spacing (Supplementary Figure S4A-C subplot a-b-c).

Comparison with Oat
The quality traits of selected quinoa genotypes (R3 and R9) were compared with oat genotypes (G3 and G7).Quinoa genotypes at cutting level 2 and 3 under spacing level 2 (23 cm) were compared with oat quality parameters.Our results revealed that the CP and CF contents of quinoa in leaves and stems were significantly higher than those of oat (Figure 6A,B).In quinoa leaves, NDF (LNDF) was significantly lower at both cuttings than in oat accessions, while ADF (LADF) was slightly higher at both cuttings than in oat (Figure 6A).In stem, however, quinoa showed desirable results as the levels of both SADF and SNDF were significantly lower than oat accessions under both cutting levels (Figure 6B).Comparing the forage quality characteristics of quinoa and oat revealed low fiber (ADF and NDF) and high crude protein and fat (CP and CF) levels in quinoa, highlighting the importance of further research into this species.
(Figure 6A).In stem, however, quinoa showed desirable results as the levels of both SADF and SNDF were significantly lower than oat accessions under both cutting levels (Figure 6B).Comparing the forage quality characteristics of quinoa and oat revealed low fiber (ADF and NDF) and high crude protein and fat (CP and CF) levels in quinoa, highlighting the importance of further research into this species.

Environmental Factors Markedly Influenced Morphological and Forage Quality Traits of Quinoa in Both Years
Climatic factors like air temperature, rainfall, and humidity are known to have significant effects on quinoa growth and biomass production.Quinoa has already devel-  .The optimal temperature for better quinoa growth has been reported as 20 • C [42].In present investigation, we also perceived the effects of the abovementioned climatic factors: in the second growing season, the quinoa gained more biomass in terms of PH, NoL, and Br compared to the first growing season.The average temperature in the second growing season in all months was closer to the quinoa optimal temperature compared to the first growing season and therefore resulted in more FW and DW (Supplementary Figure S1F,G).Quinoa can grow relatively well in humidity ranging from 40% to 88% [43].In the second growing season, the relative humidity in all growing months was more stable and favorable for quinoa growth compared with the first growing season (Figure 3).This stability and increase in relative humidity in the second growing season were mainly due to higher rainfall across all the growing months compared with the first growing season (Figure 3).Thus, the effect of these climatic factors on quinoa growth can be clearly seen from the outcomes of the present investigation where the year effect was highly significant (p < 0.05) on most morphological traits except for LA and ID.In addition to biomass, forage quality is also directly associated with plant growth, nutrients in the tissues, community composition, and forage digestibility.In our investigation, we observed a significant (p < 0.05) year effect on forage quality parameters in both leaves and stems.During the second growing season, LCP, SCP, and SCF levels increased in stems and leaves compared with the first growing season (Supplementary Figure S2A).Moreover, LCP, LCF, SCP, and SCF levels increased with the increase in forage biomass, and these findings were consistent with those of a two-year field experiment conducted in [44].

Effect of Different Treatments on Morphological and Quality Traits of Quinoa
In addition to climatic factors, various other practices such as different levels of cuttings and plant-to-plant spacing also affect quinoa biomass and quality.Cuttings at different growth stages significantly affect crop (oat, alfalfa, canola, and millet) biomass and forage quality [45][46][47][48].However, very limited research data are available on the effect of cuttings at different growth stages of quinoa.Yilmaz et al. [49] evaluated the effect of cuttings at three different stages to determine the quality and forage yield of quinoa.Their results revealed that different cutting levels had significant effects on dry forage yield and other morphological and forage quality parameters.In our study, we also observed the same effects of cuttings on quinoa morphological and quality traits.Different cutting levels significantly affected PH, NoL, Br, FW, DW, LCP, SCP, and SCF, and maximum biomass was recorded at cutting level 3 in both years (Supplementary Figure S1A-G).Furthermore, LA and ID showed mixed responses in which a few genotypes like R4, R45, and UAF performed better at cutting level 2 rather than cutting level 3.In the case of quality, the concentrations of quality traits such as LCP, LCF, SCP, SCF, LADF, LNDF, SADF, SNDF, LMois, and SMois remained lower at cutting level 1, increased with the increase in growth time, and slightly varied between cutting levels 2 and 3 (Supplementary Figure S2A-E).Our results were in agreement with multi-field trials under two cutting levels, in which maximum biomass and quality were produced at cutting level 2 rather than cutting level 1 [50].
In addition to cutting levels, plant-to-plant distance (spacing level) is another important factor affecting quinoa biomass and quality [35,51].Finding the appropriate plant-toplant distance for quinoa cultivation is a major challenge as a slight change in plant-to-plant distance can result in large variations in PH, NoL, Br, CP, CF, ADF, and NDF.Zulkadir [36] planted quinoa variety Q-52 using three different spacing levels, i.e., 20, 40, and 60 cm to determine the most appropriate plant-to-plant distance to gain maximum biomass production.Zulkadir [36] identified 20 cm as the most appropriate spacing level to gain maximum grain yield and biological yield followed by 40 and 60 cm, respectively.Similarly, Reddy et al. [52] also planted quinoa genotypes with three different spacing levels (20 × 10 cm, 25 × 10 cm, and 30 × 10 cm) along with nitrogen (N) levels to evaluate the effects of the given treatments on yield.They reported 30 × 10 cm + 25% N as the most appropriate treatment to gain maximum stover yield and seed yield, while 20 × 10 cm + 25% N was reported as the most appropriate to gain the maximum B:C ratio, net returns, and gross returns.In the present study, three different spacing levels, viz., 21, 23, and 26 cm were used to determine the appropriate spacing level to achieve maximum biomass and forage quality.Just like cutting, different levels of spacing also showed significant effects on quinoa biomass and quality.The maximum biomass in terms of PH, NoL, Br, FW, and DW was observed under spacing level two (S2), i.e., 23 cm followed by S1 (21 cm) and S3 (26 cm) (Supplementary Figure S1A-G).However, no significant effects of spacing levels were observed on LA and ID.Similarly, S2 was also more effective for quality parameters like CP, CF, NDF, and ADF in both leaves and stems except for NDF in leaves and ADF in stems (Supplementary Figure S1A-E).These findings were in agreement with the findings reported by Asher et al. [35], Temel and Yolcu [53], TEMEL and KESKIN [54], and Sief et al. [55].

Treatment Combination and Genotypes Selection Using the MGIDI
The selection of suitable genotypes has always been a challenge for plant breeders.Most plant breeders use classical stability indexes like deviation from regression, regression, and mean to select the most stable genotype [56,57].However, these statistical tools are insufficient to detect the weaknesses and strengths of genotypes and to select individuals with desired average performance and stability [58].The MGIDI is a sophisticated quantitative genetic technique used to exploit desired accessions in all crop species [59] and is free from multi-collinearity problems [60].In the current study, the MGIDI was computed for the best treatment combination (cut × space) and for the selection of suitable genotypes by grouping the studied traits into the same factor.According to the MGIDI results, cutting level 3 in combination with spacing level 2 (23 cm) was more effective in achieving the maximum biomass and quality of quinoa than other treatments (Figure 4a).Moreover, among all the genotypes, R9 and R3 gathered the desired mean performance and stability of several traits such as PH, NoL, Br, FW, DW, LCP, LCF, LADF, LMois, SCP, SCF, SADF, and SMois (Figures 4c and 5a-d).Moreover, the selection of appropriate quinoa genotypes aligned with sustainable agricultural principles has several potential benefits.Firstly, it can contribute to food security by improving crop yields and nutritional quality.Secondly, the identification of genotypes with resource-efficient traits can promote sustainable resource management, reducing water, fertilizer, and pesticide use.Additionally, by selecting genotypes with resistance to pests and diseases, farmers can reduce the need for chemical interventions, thus minimizing environmental impacts.Furthermore, the adaptability of genotypes to changing climatic conditions can enhance the resilience of agricultural systems.
Furthermore, correlation analysis revealed a highly significant and moderately significant positive correlation among the studied traits under all levels of cutting and spacing in both years.Among all productive traits, PH and NoL showed a highly significant positive correlation with most of the traits that directly contribute to biomass such as LA, Br, ID, FW, and DW.An increase in PH is directly associated with an increase in NoL [61,62], LA [63,64], FW [65,66], and DW [67,68].Similar results were seen for quality traits in which an increase in the productive traits such as PH, LA, ID, and DW led to increasing LCF, LCP, and SCF (Supplementary Figure S4A,B subplot a-b).Forage digestibility (energy) depends on the levels of anti-quality traits such as ADF and NDF [69].ADF is a fiber portion that consists of lignin and cellulose.ADF is mainly associated with the digestibility of forage and is, therefore, necessary to compute forage net energy (NE) or total digestible nutrients (TDN) for silage, haylage, and hay [70].On the other hand, NDF is the most common measure of fiber used for animal feed analysis but does not represent a unique class of chemical compounds.NDF measures most of the structural components of plant cells such as cellulose, hemicellulose, and lignin but not pectin [71].In general, both ADF and NDF values increased as the forage matured.Forages with low ADF and NDF levels are usually higher in energy [72].LADF had a highly positive significant correlation with LNDF (Supplementary Figure S4A-C subplot a-b), indicating that an increase in NDF level is directly associated with an increase in ADF.Furthermore, increases in PH, NoL, LA, FW, and DW decreased ADF and NDF levels in both leaves and stems (Supplementary Figure S4A-C subplot a-b-c).These findings were consistent with other findings, in which an increase in forage yield resulted in increases in CP and CF and decreases in NDF and ADF [43,[73][74][75].
Oat (Avena sativa L.) is commonly used as an annual forage crop in irrigated and rainfed areas worldwide.Oat has great biomass production with high forage quality [76,77].Therefore, to further analyze the fitness of quinoa as potential forage, we compared the quality traits of two selected genotypes (R3 and R9) with two oat accessions (G3 and G7).Quinoa had a greater amount of crude protein and crude fat (CP and CF) in both leaves and stems compared to oat (Figure 6A).The concentrations of anti-quality traits such as ADF and NDF were lower in quinoa leaves (except for LADF) and stems compared to oat (Figure 6B).Furthermore, this comparison revealed the latent potential of quinoa species as an alternative source of high-quality forage.

Conclusions
Our findings reveal that quinoa biomass production and quality traits are greatly influenced by environmental factors.Optimum air temperature and relative humidity are essential for the growth and development of quinoa to obtain maximum biomass and high-quality forage.In the second growing season, the environmental conditions were more favorable for quinoa than the first growing season and therefore showed remarkable growth in the second season.In addition to environmental factors, other practices such as different cutting and spacing levels also affected quinoa growth in terms of PH, NoL, Br, FW, and DW.Under all cuttings and spacing levels, cutting level 3 in combination with spacing level 2 (23 cm) was more effective for quinoa biomass and quality traits across both years.In the current study, the selected quinoa genotypes (R3 and R9) using MGIDI, had shown significant improvement in both morphological and quality traits.Furthermore, morphological traits such as PH, NoL, LA, FW, and DW showed moderately significant negative correlations with anti-quality traits, indicating that an increase in biomass can reduce concentrations of anti-nutritional quality traits to enhance the forage digestibility.The selected quinoa genotypes were then compared to the quality parameters of oat accessions.The comparison revealed that quinoa has more CP and CF and low NDF in both leaves and stems than oat.Overall, the findings indicated that quinoa has latent potential as an alternative forage to fulfill the forage demand gap in Pakistan.

Supplementary Materials:
The following supporting information can be downloaded at https://www.mdpi.com/article/10.3390/su151411446/s1.  S2: MANOVA results of 17 traits f 10 Quinoa genotypes evaluated at cutting and spacing levels for two years.Pillai's trace test was used to get a cumulative significance value for the studied main effects and interactions; Table S3: Correlation Matrix for the cutting1 × spacing 1, 2, and 3 during first growing season (year 1) for all the 17 traits of 10 Quinoa genotypes; Table S4: Correlation Matrix for the cutting 2 × spacing 1, 2, and 3 during first growing season (year 1) for all the 17 traits of 10 Quinoa genotypes; Table S5: Correlation Matrix for the cutting 3 × spacing 1, 2, and 3 during first growing season (year 1) for all the 17 traits of 10 Quinoa genotypes; Table S6: Correlation Matrix for the cutting 1 × spacing 1, 2, and 3 during second growing season (year 2) for all the 17 traits of 10 Quinoa genotypes; Table S7: Correlation Matrix for the cutting 2 × spacing 1, 2, and 3 during second growing season (year 2) for all the 17 traits of 10 Quinoa genotypes; Table S8: Correlation Matrix for the cutting 3 × spacing 1, 2, and 3 during second growing season (year 2) for all the 17 traits of 10 Quinoa genotypes.

Figure 2 .
Figure 2. Illustration for different levels of cutting at different growth stages.Each plot was divided into three equal portions for each cutting.() represent plant residues after harvesting.Note: For leaves' quality traits, the prefix "L" was used with each trait, in other words, LCP, LCF, LADF, and LNDF.Similarly, for stems, the prefix "S" was used with each trait, in other words, SCP, SCF, SADF and SNDF.

Figure 2 .
Figure 2. Illustration for different levels of cutting at different growth stages.Each plot was divided into three equal portions for each cutting.(

Figure 2 .
Figure 2. Illustration for different levels of cutting at different growth stages.Each plot was divided into three equal portions for each cutting.() represent plant residues after harvesting.Note: For leaves' quality traits, the prefix "L" was used with each trait, in other words, LCP, LCF, LADF, and LNDF.Similarly, for stems, the prefix "S" was used with each trait, in other words, SCP, SCF, SADF and SNDF.

Figure 3 .
Figure 3. Climograph of the experimental site based on the air temperature (T_min and T_max), rainfall (mm), and humidity (%) distribution during the two experimental years (four months/year).(A) represents the climatic conditions of the first growing season (year 2018-2019) i.e., from November to February while (B) represents the climatic conditions of the second growing season (year 2019-2020) i.e., from November to February.

Figure 3 .
Figure 3. Climograph of the experimental site based on the air temperature (T_min and T_max), rainfall (mm), and humidity (%) distribution during the two experimental years (four months/year).(A) represents the climatic conditions of the first growing season (year 2018-2019) i.e., from November to February while (B) represents the climatic conditions of the second growing season (year 2019-2020) i.e., from November to February.

Figure 4 .
Figure 4. Cutting × spacing interaction ranking based on the MGIDI: (a) the selected combination (highlighted in red) based on their index value, considering 20% selection intensity (red line); (b) PCA plot shows the distribution of the studied traits and treatments; (c) selection differentials (%)

Figure 4 .
Figure 4. Cutting × spacing interaction ranking based on the MGIDI: (a) the selected combination (highlighted in red) based on their index value, considering 20% selection intensity (red line); (b) PCA plot shows the distribution of the studied traits and treatments; (c) selection differentials (%) for morphological, stem, and leaf quality traits were obtained from the selected interactions.Facets grouped the traits based on desired SDs.

Figure 5 .
Figure 5.The MGIDI ranks genotypes in descending order.(a) The non-selected and selected genotypes are shown in gray and red circles, respectively, and depending on the intensity (10%), the red circle depicts the cut point.(b) The proportion of every factor on the computed MGIDI represents the strengths and weaknesses of the selected genotypes.The smaller the ratio defined by a factor (near the outer edge), the closer the traits within that factor are to the ideotype.Dashed line represents the theoretical value if all the factors contributed equally.(c) Principal component analysis (PCA) biplot performed with the studied traits for the substrate factor.(d) Selection differentials (%) for morphological and stem and leaf quality traits were obtained from the selected genotypes.The traits are grouped by facets according to the desired selection differentials.

Figure 5 .
Figure 5.The MGIDI ranks genotypes in descending order.(a) The non-selected and selected genotypes are shown in gray and red circles, respectively, and depending on the intensity (10%), the red circle depicts the cut point.(b) The proportion of every factor on the computed MGIDI represents the strengths and weaknesses of the selected genotypes.The smaller the ratio defined by a factor (near the outer edge), the closer the traits within that factor are to the ideotype.Dashed line represents the theoretical value if all the factors contributed equally.(c) Principal component analysis (PCA) biplot performed with the studied traits for the substrate factor.(d) Selection differentials (%) for morphological and stem and leaf quality traits were obtained from the selected genotypes.The traits are grouped by facets according to the desired selection differentials.

Figure 6 .
Figure 6.Comparison of leaf and stem qualities of selected quinoa genotypes with oat varieties being used as standard for forage.Oat data were recorded at the grain-filling stage, while quinoa genotype data at both cutting levels (C2 and C3; S: Spacing) were utilized for comparing the outcomes.(A) Leaf quality parameters: LCP: Leaf crude protein; LCF: leaf crude fat; LADF: leaf acid detergent fiber; LNDF: leaf neutral detergent fiber; LMois: leaf moisture.(B) Stem quality parameters: SCP: Stem crude protein; SCF: stem crude fat; SADF: stem acid detergent fiber; SNDF: stem neutral detergent fiber; SMois: stem moisture.The bars represent the standard error of mean.

Figure 6 .
Figure 6.Comparison of leaf and stem qualities of selected quinoa genotypes with oat varieties being used as standard for forage.Oat data were recorded at the grain-filling stage, while quinoa genotype data at both cutting levels (C2 and C3; S: Spacing) were utilized for comparing the outcomes.(A) Leaf quality parameters: LCP: Leaf crude protein; LCF: leaf crude fat; LADF: leaf acid detergent fiber; LNDF: leaf neutral detergent fiber; LMois: leaf moisture.(B) Stem quality parameters: SCP: Stem crude protein; SCF: stem crude fat; SADF: stem acid detergent fiber; SNDF: stem neutral detergent fiber; SMois: stem moisture.The bars represent the standard error of mean.

1 .
Environmental Factors Markedly Influenced Morphological and Forage Quality Traits of Quinoa in Both Years Climatic factors like air temperature, rainfall, and humidity are known to have significant effects on quinoa growth and biomass production.Quinoa has already developed resistance to various abiotic factors and can resist a broad range of air temperature ranging from −8 to 35 • C [41] Figure S1: Morphological Traits of 10 Quinoa genotypes studied during two years at different cutting and spacing levels; Figure S2: Leaf Quality Traits of 10 Quinoa genotypes studied during two years at different cutting and spacing levels; Figure S3: Stem Quality Traits of 10 Quinoa genotypes studied during two years at different cutting and spacing levels; Figure S4: Correlations among matrixes of various traits under different cuts and spacing levels; Table S1: Split-Split Plot ANOVA of both Years (2018-2019 & 2019-2020) for all Traits of 10 Quinoa Genotypes tested at different Cuttings and Spacing Levels; Table

Table 1 .
Predicted selection differentials, goals, and trait distributions of factors for the studied cutting × spacing interaction.

Table 2 .
Predicted selection differentials, goals, and trait distributions of factors for the studied genotypes.