Forage Performance and Detection of Marker Trait Associations with Potential for Napier Grass ( Cenchrus purpureus ) Improvement

: The evaluation of forage crops for adaptability and performance across production systems and environments is one of the main strategies used to improve forage production. To enhance the genetic resource base and identify traits responsible for increased feed potential of Napier grass, forty-ﬁve genotypes from Empresa Brasileira de Pesquisa Agropecu á ria (EMBRAPA), Brazil, were evaluated for forage biomass yield and feed nutritional quality in a replicated trial under wet and dry season conditions in Ethiopia. The results revealed signiﬁcant variation in forage yield and feed nutritional qualities among the genotypes and between the wet and dry seasons. Feed ﬁber components were lower in the dry season, while crude protein, in vitro organic matter digestibility, and metabolizable energy were higher. Based on the cumulative biomass and metabolizable energy yield, top performing genotypes were identiﬁed that are candidates for future forage improvement studies. Furthermore, the marker-trait association study identiﬁed diagnostic single nucleotide polymorphisms (SNP) and SilicoDArT markers and potential candidate genes that could di ﬀ erentiate high biomass yielding and high metabolizable energy genotypes in the collection.


Introduction
Forage provides an essential component of the feed resource for small scale dairy-producing farmers in developing countries [1]. Production by the farmers, however, is affected by seasonal fluctuations of forage performance, both in terms of quantity and quality. The productivity of forage crops is also hampered by abiotic and biotic stresses [2,3]. Hence, several national research institutes are aiming to improve forage productivity and to provide new forage cultivars by the identification of adaptable, higher quality forage accessions, mainly through the introduction and evaluation of different options in their forage development programs [4,5].
Napier grass, (Cenchrus purpureus Schumach., syn. Pennisetum purpureum Schumach), also called elephant grass, is one of the most popular tropical forage crops [6]. Napier grass is a C4 perennial grass species of the Poaceae family, native to Sub-Saharan Africa, and it is widely distributed across the tropical and subtropical parts of the world [7]. Napier grass is valued for its high biomass production, pest resistance, feed quality, and year-round availability due to its perennial nature [8,9]. The grass is used in the 'push-pull' strategy for the management of insect pests that can cause significant damage to maize production, and it can grow in degraded and marginalized lands; hence, it can play a role in preventing soil erosion [10]. Due to its superior biomass yield, Napier grass is also a candidate species for biofuel production [11]. To exploit the potential of this grass, various conservation centers maintain and develop superior cultivars for forage production. Among these centers, the forage genebanks of the International Livestock Research Institute (ILRI) in Ethiopia and the Brazilian Agricultural Research Corporation, Empresa Brasileira de Pesquisa Agropecuária (EMBRAPA), maintain a large number of Napier grass germplasm that has been acquired from a range of different countries.
Napier grass has limited genomic information available, and, until recently, only a handful of molecular markers have been applied, mainly in studies targeting the assessment of genetic diversity. Although the recent advances in the genotyping by sequencing (GBS) approach have enabled the generation of a large number of high-density genome-wide markers [12][13][14], this species still lacks a reference genome. Hence, the reference genome sequence of pearl millet [15], a closely related species to Napier grass, has been used for mapping the genome-wide markers.
Several studies have been directed towards the characterization and evaluation for forage yield, feed quality, and resistance to the major diseases affecting Napier grass [16]. From the ILRI genebank collection, some genotypes shown to be promising for high biomass yield and nutritional quality have been identified and utilized by national programs, including in Ethiopia and Kenya [17][18][19]. Due to the impacts of climate change in East and Central Africa (ECA) and the potential risks from the emergence and spread of pests and diseases [20], the selection and development of new forage varieties with desirable traits that fit with the climatic conditions and production systems of ECA needs continuous assessment of a large number of genotypes. The success of such improvement programs will depend on the available gene pool and the genetic diversity contained within. In line with these goals, a collaborative study was conducted between ILRI and EMBRAPA to identify unique genotypes, including some elite lines, from the collections of both centers [21]. As a result of this targeted introduction, 53 Napier grass genotypes were transferred to the ILRI gene bank from EMBRAPA in 2013, an activity that greatly enhanced the genetic diversity contained in the ILRI collection and provided potential new options for future utilization.
Forage quality and the production potential of an individual genotype can be location-and season-specific; hence, the assessment of plant performance and adaptability to a new environment are major activities employed when new genotypes are introduced into a collection. Consequently, the introduced Napier grass genotypes were evaluated for their agronomic performance and feed quality under wet and dry season conditions in the Ethiopian environment. Here, we report on the performance of 45 Napier grass genotypes for forage biomass yield (total fresh weight (TFW) and total dry weight (TDW)]) and feed quality traits in the wet and dry seasons of Bishoftu, Ethiopia, using field phenotypic evaluation and molecular analyses. We also report on the identification of single nucleotide polymorphisms (SNPs), SilicoDArT markers, and potential candidate genes associated with high TDW and metabolizable energy. In general, the results demonstrate that there is significant variation in growth and performance among the genotypes and between the wet and dry seasons. Marker-trait associations identified discriminatory SNP and SilicoDArT markers that were able to differentiate high biomass yielding genotypes in the collection. One of these markers aligned with markers previously identified to be associated with improved agronomic performance in pearl millet (Pennisetum glaucum), a closely related species to Napier grass. This finding may be used in the future to improve the selection efficiency in Napier grass through the application of marker assisted selection (MAS).

Description of the Experimental Area
The study was conducted in Bishoftu, Ethiopia (008 • 47 20 N and 038 • 59 15 E), over a 21-month period, from May 2016 to January 2018. The altitude of Bishoftu is 1890 meters above sea level with an Alfisol soil type. The monthly precipitation and temperature for the harvest period of 2016 to 2017 are provided in Figure 1. During the dry periods, plants were supplemented with water by sprinkler irrigation.

Description of the Experimental Area
The study was conducted in Bishoftu, Ethiopia (008°47′20′′N and 038°59′15′′E), over a 21-month period, from May 2016 to January 2018. The altitude of Bishoftu is 1890 meters above sea level with an Alfisol soil type. The monthly precipitation and temperature for the harvest period of 2016 to 2017 are provided in Figure 1. During the dry periods, plants were supplemented with water by sprinkler irrigation.

Planting Materials and Trial Establishment
Forty-five genotypes of Napier grass from EMBRAPA were used for the study (Table 1). Cuttings from established Napier grass plants were planted in Bishoftu using a randomized complete block design, with three replications, at the start of the main rainy season on 6 May 2016. Each genotype was planted in a single 5 m long row, with a spacing of 1 m between rows and 0.5 m between plants within a row. Four months later, on 8 September 2016, the established plants were uniformly cut at 5 cm above ground level. Land preparation, planting, weeding, harvesting, and related management practices were uniformly applied to all plots. Inorganic fertilizer, urea, was applied during the rainy season in 2017 at a rate of 6.2 g/plant.

Planting Materials and Trial Establishment
Forty-five genotypes of Napier grass from EMBRAPA were used for the study (Table 1). Cuttings from established Napier grass plants were planted in Bishoftu using a randomized complete block design, with three replications, at the start of the main rainy season on 6 May 2016. Each genotype was planted in a single 5 m long row, with a spacing of 1 m between rows and 0.5 m between plants within a row. Four months later, on 8 September 2016, the established plants were uniformly cut at 5 cm above ground level. Land preparation, planting, weeding, harvesting, and related management practices were uniformly applied to all plots. Inorganic fertilizer, urea, was applied during the rainy season in 2017 at a rate of 6.2 g/plant.

Data Collection
Plants were harvested after every 8 weeks of regrowth for eight consecutive harvests between November 2016 and January 2018. During each harvest, biomass data were collected from the middle 8 plants in each row. To determine the forage biomass yield, i.e., total fresh weight (TFW) and total dry weight (TDW), forage harvesting was undertaken by chopping with a sickle, leaving a residual stubble height of 5 cm. TFW was measured immediately after each harvest. A sub-sample of 300 g of whole plant material was then chopped and placed in a labeled paper bag. An additional 600 g sample was separated into leaf and stem, weighed immediately, chopped, and placed into separate labelled paper bags for leaf and stem. The labeled paper bags were oven dried at 60 • C for 48 h. The sub-samples from the bags were used to calculate the TDW, and the leaf to stem ratio (LSR) was calculated by dividing the leaf dry weight by the stem dry weight of each genotype.
For genetic analysis, leaf sample collection, DNA extraction, and genotyping were conducted, as described previously by Muktar et al [12].

Feed Quality Analysis
The oven-dried tissue samples from the whole plant, leaf, and stem components were ground separately until fine enough to pass through a 1 mm sieve. The ground samples were scanned using Near Infrared Spectroscopy (NIRS) (FOSS Forage Analyzer 5000 with software package WinISI II) to predict feed quality traits, as indicated by Choudhary et al. [22]. The traits analyzed were: acid detergent fiber (ADF), neutral detergent fiber (NDF), acidic detergent lignin (ADL), organic matter (OM), dDry matter (DM), inorganic matter (ash), total nitrogen (N), crude protein (CP), in vitro organic matter digestibility (IVOMD), and metabolizable energy (ME), i.e., values were predicted using previously developed equations calibrated against conventional wet chemistry analyses for Napier grass. The predicted values were corrected by DM, and the percentage values were used for statistical analyses. The ME yield was calculated by multiplying ME by TDW.

Data Analysis
In order to improve the precision and to reliably estimate genetic values for genotypes, a spatial analysis was conducted using TFW, TDW, and LSR. As described by Velazco et al. [23], the mixed model features of spatial analysis were used to fit a smooth surface to account for all sources of environmental variation. We, therefore, used the first-order autoregressive (AR1) model [24] in GenStat software version 19 [25], where blocks (columns) and plots (rows) in each block were used as random terms and genotypes as fixed terms. The plot values were then fitted based on the spatial variations. The fitted data were subjected to an analysis of variance (ANOVA) using GenStat software version 19 [25] to determine the significance of the main effects and the interactions using the model set out below.
where Y ijk is the response (Biomass yield and chemical composition) during each season, µ = overall mean, Gi = effect of Napier grass genotypes, S j = effect of season of harvest (j = dry and wet), G * S ij = the interaction of ith genotype and jth season, and ε ijk = the residual error. The least significant difference (LSD), for comparison of mean values, was employed to compare genotypes for traits with significant p values. Genetic parameters, genotypic coefficient of variation (GCV), and phenotypic coefficient of variation (PCV) were estimated according to Singh and Chaudary [26], using the formulae: where GCV = genotypic coefficient of variation, PCV = phenotypic coefficient of variation, σ 2 g = genotypic variance, σ 2 p = phenotypic variance, and X = grand mean.
To evaluate Napier genotypes based on annual TDW and ME, cumulative TDW and ME yield were computed. The cumulative annual TDW was obtained by adding the mean TDW of wet and dry season harvests of each genotype. The cumulative ME yield was obtained by adding the mean ME of wet and dry season harvests and multiplying by the respective cumulative TDW.
Pattern analyses, i.e., cluster analysis and principal component analysis (PCA), were applied to visualize the pattern of genotypes using DeltaGen 3_1 [27]. Pearson correlation analysis among forage biomass and feed quality traits were conducted using GenStat software (version 19, VSN international LTD, UK) [25].
To examine scaling relationships between biomass yield (TDW) with feed quality traits and leaf to stem ratio, linear regression on log-transformed trait values were conducted using the GenStat platform [25]. The scaling relationship was analyzed based on the allometric equation of Gilles and Gilles [28].
where Y is the response variable, i.e., feed quality or leaf to stem ratio traits, and X is the fixed variable, i.e., total dry weight; a is the y-intercept, and b is the slope of the scaling function, representing the allometric exponent.
To represent the relationship between total dry weight with feed quality traits in wet and dry seasons a scatter plot was drawn using observed and predicted values of linear regression.

Marker-Trait Association Analysis
Marker-trait associations were detected with R statistical software (https://www.r-project.org) using the R package GWASpoly [29], which facilitates taking kinship into account in a mixed model to correct against spurious associations due to the level of relatedness between genotypes. In addition, the simple ANOVA model, which assumes that genotypes are independent, was employed to test the association of each marker with the trait. The collection was genotyped using the DArTseq platform (http://www.diversityarrays.com) to generate genome wide (SNP and SilicoDArT) markers, as described previously [12]. The genome-wide marker data, together with the phenotype data for annual dry weight production and metabolizable energy, were used in the marker-trait association analysis.
The DArTseq sequence reads corresponding to the SNP and SilicoDArT markers were aligned to Napier grass genomic [30] and transcriptome [31] sequences using the bwa mem sequence alignment method [32], and the sequences were annotated using the GeneBank NCBI blastx tool (https://blast. ncbi.nlm.nih.gov/Blast.cgi).

Forage Biomass Productivity of Napier Grass Genotypes
In order to assess the adaptation and performance of the Napier grass genotypes sourced from EMBRAPA, a statistical analysis was conducted using analysis of variance (ANOVA). Highly significant (p < 0.001) variation for total fresh weight (TFW), total dry weight (TDW), and leaf to stem ratio (LSR) was observed between genotypes and harvest seasons ( Table 2). The results show that, with the exception of LSR, seasonal variation was more significant compared to the genotypes contribution. As expected, forage yield was lower during the dry season, with the mean total dry weight for wet season harvests almost three-fold higher compared to the dry season harvests (Table 3). Similarly, LSR was over two-fold higher during the dry season compared to the wet season, presumably reflecting the fact that genotypes produce a higher proportion of stem associated with the increased yield during eight weeks of regrowth in the wet season (Table 3). In line with this finding, the top yielding genotypes during the wet season harvests, BAGCE 94, BAGCE 86, and BAGCE 100, had a low to medium LSR reflecting the contribution of stem elongation to the high yield. The highest LSRs during the dry season were obtained from genotypes BAGCE 7, BAGCE 43, and BAGCE 25, while the lowest LSRs were observed from BAGCE 80, BAGCE 56, and CNPGL-93-06-1 (Table A3; Appendix B). Total dry weight (TDW), total fresh weight (TFW), leaf to stem ratio (LSR), mean square (MS), coefficient of determination (R 2 ); genotype by season interaction (G * S); coefficient of variation (CV); * p < 0.05 (5%), ** p < 0.01 (1%), and *** p < 0.001 (0.1%). Total dry weight (TDW), total fresh weight (TFW), leaf to stem ratio (LSR), genotypic coefficient of variation (GCV), and phenotypic coefficient of variation (PCV).
The phenotypic coefficient of variation (PCV) and the genotypic coefficient of variation (GCV) were calculated to assess the contribution of these factors to the forage biomass yield traits ( Table 3).
The GCV values of both TFW and TDW for the wet season harvests were higher than for the dry season harvests and, conversely, PCV were higher in the dry season than the wet season harvests. This indicates that the environment had a greater influence on forage production in the dry season than in the wet season. However, overall the GCV values of both traits were close to their respective PCV values, indicating that genetic factors contributed a major effect to the variation observed in the trial. On the other hand, PCV values for LSR were more than two-fold higher than the GCV values, indicating that the environment plays a significant role in the variation observed among genotypes for this trait ( Table 3). Some of the top performing genotypes identified for TFW and TDW during the wet season were also top performers during the dry season. For example, the genotypes BAGCE 94 and BAGCE 100 were among the highest yielding in both seasons, while genotypes CNPGL 92-190-1, BAGCE 75, and CNPGL 96-21-1 were amongst those that produced the lowest (Table A3; Appendix B).

Feed Qualities of Napier Grass Genotypes
The results of the feed quality analysis are shown in Table 4. Significant differences were observed among the genotypes in the traits: dry matter (DM); organic matter (OM); inorganic matter (ash); acid detergent fiber (ADF); acid detergent lignin (ADL); neutral detergent fiber (NDF); crude protein (CP); total nitrogen (N); in vitro organic matter digestibility (IVOMD); and metabolizable energy (ME) in both the leaf and stem (L + S) and the leaf (L) tissue samples taken across harvests (Table 4). In addition, highly significant (p < 0.001) seasonal variation for most feed quality parameters, with the exception of DM, OM, and ash, were observed, but the genotype-by-season interaction was not significant. The data on fiber fractions showed that, in general, the genotypes had a higher percentage of NDF, ADF, and ADL in the wet season than in the dry season, which had significantly lower (p < 0.05) mean values (Table 5). Conversely, the mean values for CP, IVOMD, and ME were significantly higher (p < 0.05) during the dry season (Table 5). During the wet season, PCV values for most feed quality parameters were almost two-fold higher than their respective GCV values (Table 5). Similarly, during the dry season, the PCV values were also approximately two-fold higher than the respective GCV values (Table 5). Taken together, these results indicate that the environmental effect was equally as important as the genotypic effect for the observed variation in feed nutrition quality traits. Generally, genotypes that showed highest mean performance for NDF had the lowest mean performance for CP, IVOMD, and ME in both wet and dry seasons (Table A4; Appendix B). Feed quality traits are expressed on a dry matter basis (%DM); dry matter (DM); organic matter (OM); minerals (ash); neutral detergent fiber (NDF); acid detergent fiber (ADF); acidic detergent lignin (ADL); total nitrogen (N); crude protein (CP); in vitro organic matter digestibility (IVOMD); metabolizable energy (ME); mean square (MS), coefficient of determination (R 2 ); genotype by season interaction (G * S); coefficient of variation (CV); ns = non-significant, * p < 0.05 (5%), ** p < 0.01 (1%), and *** p < 0.001 (0.1%). Feed quality traits are expressed on a dry matter basis (%DM); genotypic coefficient of variation (GCV); phenotypic coefficient of variation (PCV); dry matter (DM); organic matter (OM); minerals (ash); neutral detergent fiber (NDF); acid detergent fiber (ADF); acidic detergent lignin (ADL); total nitrogen (N); crude protein (CP); in vitro organic matter digestibility (IVOMD); metabolizable energy (ME).

Phenotypic Variability and Association Between Biomass Yield and Feed Quality Traits
In order to understand the variation between genotypes and the correlation between traits, a pattern analysis was performed to depict patterns of genotypes using forage biomass and feed quality traits across the wet and dry season harvests ( Figure 2; Table 6 and Figure A1, Appendix B). The biplot from the pattern analysis revealed that the first two principal components explain 70.8% and 68.3% of the total variation during the wet and dry season harvests, respectively ( Figure 2). The first component (PC1) accounted for 53.0% and 48.4% of the total variation in the wet and dry season harvests, respectively. Of the different traits, CP and N were the main contributors during the wet season, while NDF and ADF were the main contributors during the dry season. The second component (PC2), mainly associated with DM, IVOMD, and ME, accounted for 17.8% and 19.9% of the total variation in the wet and dry season harvests, respectively. The biplot analyses revealed three cluster groups where genotypes in the first group are high in CP, IVOMD, and ME, while genotypes in the second group are high in TFW, TDW, NDF, and ADF, indicating that these groups of genotypes differ in their feed qualities and biomass yield. The genotypes in the third cluster were distributed between the two groups. For example, genotypes BAGCE 94, BAGCE 100, and CNPGL 94-13-1 are clustered in a group that performed best for forage biomass yield and were also high for fiber components, while CNPGL 92-133-3, CNPGL 00-1-1, and CNPGL 94-07-2 were clustered in a group that showed the lowest in IVOMD, ME, and CP percentages in both seasons, although this was associated with a lower biomass yield ( Figure 2).

Phenotypic Variability and Association Between Biomass Yield and Feed Quality Traits
In order to understand the variation between genotypes and the correlation between traits, a pattern analysis was performed to depict patterns of genotypes using forage biomass and feed quality traits across the wet and dry season harvests ( Figure 2 ; Table 6 and Figure A1, Appendix B). The biplot from the pattern analysis revealed that the first two principal components explain 70.8% and 68.3% of the total variation during the wet and dry season harvests, respectively ( Figure 2). The first component (PC1) accounted for 53.0% and 48.4% of the total variation in the wet and dry season harvests, respectively. Of the different traits, CP and N were the main contributors during the wet season, while NDF and ADF were the main contributors during the dry season. The second component (PC2), mainly associated with DM, IVOMD, and ME, accounted for 17.8% and 19.9% of the total variation in the wet and dry season harvests, respectively. The biplot analyses revealed three cluster groups where genotypes in the first group are high in CP, IVOMD, and ME, while genotypes in the second group are high in TFW, TDW, NDF, and ADF, indicating that these groups of genotypes differ in their feed qualities and biomass yield. The genotypes in the third cluster were distributed between the two groups. For example, genotypes BAGCE 94, BAGCE 100, and CNPGL 94-13-1 are clustered in a group that performed best for forage biomass yield and were also high for fiber components, while CNPGL 92-133-3, CNPGL 00-1-1, and CNPGL 94-07-2 were clustered in a group that showed the lowest in IVOMD, ME, and CP percentages in both seasons, although this was associated with a lower biomass yield ( Figure 2). The degree of correlation among forage biomass yield (TFW and TDW) and feed quality traits are presented in Table 6. Forage biomass yield was positively correlated with NDF, ADF, ADL, and OM and negatively correlated with LSR, ash, CP, and N, in both wet and dry season harvests. Strong positive correlations were observed between the fiber traits, while the correlations of fiber traits with ash, CP, N, IVOMD, and ME were negative. The scaling relationship between TDW and feed quality traits (NDF, CP, IVOMD, and ME) is presented in Table 7. The slopes of linear regression between TDW and feed quality traits were statistically significant, except for LSR in the wet season, indicating dry season conditions. Total fresh weight (TFW); total dry weight (TDW); dry matter (DM); organic matter (OM); inorganic matter (ash); neutral detergent fiber (NDF); acid detergent fiber (ADF); acidic detergent lignin (ADL); crude protein (CP); in vitro organic matter digestibility (IVOMD); metabolizable energy (ME). Feed quality traits are expressed on dry matter bases (%DM) and forage yields are expressed in (t/ha). The blue, green, and pink colors depict the three cluster groups of Napier genotypes.
The degree of correlation among forage biomass yield (TFW and TDW) and feed quality traits are presented in Table 6. Forage biomass yield was positively correlated with NDF, ADF, ADL, and OM and negatively correlated with LSR, ash, CP, and N, in both wet and dry season harvests. Strong positive correlations were observed between the fiber traits, while the correlations of fiber traits with ash, CP, N, IVOMD, and ME were negative. The scaling relationship between TDW and feed quality traits (NDF, CP, IVOMD, and ME) is presented in Table 7. The slopes of linear regression between TDW and feed quality traits were statistically significant, except for LSR in the wet season, indicating that feed quality traits are allometrically associated with TDW ( Table 7). The scaling slope increased slightly in the dry season compared to the wet season for the relationship between TDW and CP, IVOMD, and ME, while the scaling relationship between TDW and NDF was similar in both seasons. The scatter plot of the observed and predicted values indicated a decrease of CP, IVOMD, and ME as TDW increased, and NDF increased as TDW increased ( Figure 3). Consequently, the results from both the pattern analysis and correlation analysis indicate that selection of genotypes for increased biomass yield could also increase the amount of fiber and simultaneously reduce the digestibility and protein content of the forage. Taken together, these results indicate the need for careful consideration when selecting for improved performance of Napier grass to support livestock productivity.
Agronomy 2020, 10, x FOR PEER REVIEW 2 of 29 that feed quality traits are allometrically associated with TDW ( Table 7). The scaling slope increased slightly in the dry season compared to the wet season for the relationship between TDW and CP, IVOMD, and ME, while the scaling relationship between TDW and NDF was similar in both seasons. The scatter plot of the observed and predicted values indicated a decrease of CP, IVOMD, and ME as TDW increased, and NDF increased as TDW increased ( Figure 3). Consequently, the results from both the pattern analysis and correlation analysis indicate that selection of genotypes for increased biomass yield could also increase the amount of fiber and simultaneously reduce the digestibility and protein content of the forage. Taken together, these results indicate the need for careful consideration when selecting for improved performance of Napier grass to support livestock productivity.
(a) (b) Figure 3. Relationship between neutral detergent fibre (NDF) ,crude protein (CP), invitro organic matter digestibility (IVOMD) and metabolizable energy (ME) with total dry weight (TDW) of 45 Napier grass genotypes during wet season (purple) and dry season (blue) harvests of observed data (a) and fitted data using allometric relationship (b) under Bishoftu, Ethiopia. Feed quality traits are expressed on dry matter bases (%DM) and TDW is expressed in (t/ha).    All data were log10-transformed before analysis; Leaf to stem ratio (LSR), neutral detergent fiber (NDF), crude protein (CP) in vitro organic matter digestibility (IVOMD), metabolizable energy (ME), total dry weight (TDW), R 2 , determination coefficient, standard error (S.E.). Test1 (p-value) indicates the statistical significance of the slope value being different from 1.

Cumulative Dry Weight and Metabolizable Energy Production of Genotypes
In order to evaluate forage performance, genotypes were assessed for the combination of TDW and ME as these traits are major components underpinning feed value. Genotypes BAGCE 94, BAGCE 100, BAGCE 53, and CNPGL 94-13-1 were highest for annual TDW. There was no significant differences in ME; in terms of MJ/Kg DM/year, the genotypes that produced the highest annual TDW also produced the highest ME yield (Table 8).

Marker-Trait Associations
For the association study, annual dry weight yields and the amount of metabolizable energy quantified for the 45 EMBRAPA Napier grass genotypes were used. Normality analysis showed the normal distribution of the data (Figure 4), which was also supported by a non-significant (Table 9) Shapiro-Wilk normality test [33], demonstrating the suitability of the data for analysis of variance (ANOVA) and genome-wide association study (GWAS) analyses. The broad-sense heritability values, which show the amount of the phenotypic variation in a population, attributed to the genetic differences between individuals, were calculated from the phenotypic variance across replications, which revealed broad-sense heritability values of above 60%, with 64% for annual dry weight yield and 65% for metabolizable energy (Table 9). Hence, the estimated broad-sense heritability values indicated that the variation observed in the traits was mainly due to genetic factors rather than environmental factors, further demonstrating the suitability of the traits in the collection for marker-trait association studies. Most feed quality component traits showed low heritability values and hence were not suitable for marker-trait association analysis.

Marker-Trait Associations
For the association study, annual dry weight yields and the amount of metabolizable energy quantified for the 45 EMBRAPA Napier grass genotypes were used. Normality analysis showed the normal distribution of the data (Figure 4), which was also supported by a non-significant (Table 9) Shapiro-Wilk normality test [33], demonstrating the suitability of the data for analysis of variance (ANOVA) and genome-wide association study (GWAS) analyses. The broad-sense heritability values, which show the amount of the phenotypic variation in a population, attributed to the genetic differences between individuals, were calculated from the phenotypic variance across replications, which revealed broad-sense heritability values of above 60%, with 64% for annual dry weight yield and 65% for metabolizable energy (Table 9). Hence, the estimated broad-sense heritability values indicated that the variation observed in the traits was mainly due to genetic factors rather than environmental factors, further demonstrating the suitability of the traits in the collection for markertrait association studies. Most feed quality component traits showed low heritability values and hence were not suitable for marker-trait association analysis.  Table 9. Estimated generalized heritability (H 2 g), ANOVA and normal distribution of the phenotypic traits used in marker-trait association analysis.   Table 9. Estimated generalized heritability (H 2 g ), ANOVA and normal distribution of the phenotypic traits used in marker-trait association analysis. Molecular markers used in the analysis were as previously described [16]. For the genome-wide association analysis, a total of 87,289 markers, comprising of 25,773 SNP and 61,516 silicoDArT markers, were used. The markers minor allele frequency (MAF) and missing values were greater than 5% and less than 10%, respectively. Association analysis using two statistical models, the mixed model using kinship matrix as a covariate (MM) and ANOVA, identified 64 SNP and 24 silicoDArT markers that were associated with dry weight yield at a −log10 p value threshold of > 2.5, using both models (Table A5; Appendix B). Similarly, 60 SNP and 23 silicoDArT markers were associated (−log10 p > 2.5) with the metabolizable energy trait. From the identified 88 markers, 82 markers were found to be associated with both traits, which would be anticipated as the dry weight is used in the quantification of cumulative metabolizable energy. A number of the markers were used to detect most of the high performing genotypes identified in the agronomic performance analysis. Although the genomic location of the majority of the associated markers was not identified in the pearl millet reference genome, at least one associated marker was detected on each of the seven chromosomes. Eleven of the associated markers were mapped on chromosome 3, potentially indicating the importance of the chromosome for future association analyses.

Analysis of Nucleotide Sequences Corresponding to the Associated Markers
Nucleotide sequences of the associated markers were retrieved from the available Napier grass genomic [30] and transcriptomic [31] sequence data by aligning the DArTseq short sequence reads (approximately 55 nucleotides (nt)) corresponding to each SilicoDArT marker (Supplementary Table S1). Seven of the SilicoDArT markers were aligned with the assembled genomic sequences [30]. The genomic sequences ranged from 278 to 6143 nt in length, which is suitable for the design of PCR primers, amplicon sequencing and in silico functional analysis of the associated markers. A total of 36 markers were able to be aligned with the transcriptome sequences (approximately 150 nt in length), indicating that many of the associated markers are found within the gene coding regions. For 23 of the markers, sequence annotation information was retrieved from the GeneBank NCBI databases using the blastx tool. Most of the annotated sequences were detected in the transcriptome sequences. Out of the eleven associated markers located on chromosome 3, five were annotated as a cytochrome P450, glutamate decarboxylase 2, mitochondrial fission protein ELM1 isoform X2, inositol transporter 2, and U-box domain-containing protein 35-like.

Seasonal Variation Determines Forage Yield of Napier Grass Genotypes
The performance of introduced EMBRAPA Napier grass genotypes from Brazil was evaluated at Bishoftu, Ethiopia, for forage biomass yield and feed quality traits under sequential wet and dry season harvests, where the wet season is between May and September, while dry season extends from October to April in both 2016 and 2017 (Figure 1). Generally, the study revealed marked differences among genotypes for forage biomass yield and feed quality traits, which indicates the existence of genotypic variability amongst the genotypes. Similarly, season and genotype-by-season interaction significantly affected most of the traits, but their effect was more pronounced in forage biomass yield than feed nutrition quality traits. The seasonal change was reflected by the distinct distribution of rainfall (Figure 1) in the wet and dry seasons, which could greatly influence Napier grass plants growth and development in the respective seasons. For example, Napier grass genotypes in the trial had a lower LSR in the wet season than in the dry season, indicating that plant growth had been stimulated by the rains resulting in increased stem development during the wet season. Hence, for efficient exploitation of Napier grass, evaluation of forage performance relies on the evaluation of genotypes across years, particularly during the dry season, where major feed shortages affect livestock production and productivity.

Forage Yield of Napier Grass Genotypes
Forage biomass yield is one of the major targets in forage crop improvement, where fast-growing and high LSR types are preferred [34]. The current study revealed that TFW and TDW were higher in the wet season than in the dry season, indicating that moisture availability was a determinant factor for growth and development of Napier grass. For example, top productive genotypes BAGCE 94, BAGCE 100, and BAGCE 30 had mean TDW of 5.7 t/ha, 5.5 t/ha, and 4.9 t/ha, respectively, during the wet season, while in the dry season, their respective mean TDW were 2.5 t/ha, 1.64 t/ha, and 1.46 t/ha, respectively (Table A3; Appendix B). These results indicate that the potential performance for high levels of forage production is associated with adequate moisture levels in the soil. With regards to TDW in the wet season, TDW in the current study was much lower than previously reported [19,35] for some Napier grass genotypes from the ILRI collection grown under Ethiopian conditions. While these variations might be attributed to differences in genotype performance, the influences of the prevailing environmental conditions and agronomic management practices among these studies could not be ruled out.
In relation to LSR, the proportion of leaf to stem was higher in the dry season than the wet season, which indicated slower stem development during the dry season. In the wet season, the highest mean LSR was obtained from genotypes CNPGL 93-08-1, BAGCE 63, and BAGCE 80, with mean values of 4.9%, 4.8%, and 4.1%, respectively, while mean LSR values for top biomass producing genotypes had low to medium LSR values. Thus, top yielding genotypes can be classified as fast growing with a higher stem proportion, which might have implications in feed nutritive quality because leaves have been shown to have a high level of feed quality compared to the stem fraction [36].
Significant variation for TFW, TDW, and LSR was observed among genotypes, regardless of the season (wet or dry), indicating the existence of genetic variability among the Napier grass genotypes. As Napier grass is a perennial fodder, the identification of consistently productive cultivars across the seasons and years is an important parameter. Results from the current study showed that top biomass yielding genotypes in the dry season were also top producing during the wet season. Genotypes that show consistently high TDW throughout the wet and dry seasons were BAGCE 94, BAGCE 100, and BAGCE 53. These genotypes performed better than some of the current ILRI 'best bet' genotypes planted in a replicated trial in adjacent plots (Tables A1 and A2, Appendix A), indicating the potential of these genotypes to be used for improved biomass yield in the tested environment. Under the Brazilian environment, however, BAGCE 53 was an early flowering and low biomass producing genotype [37], which indicates differential performance of this genotype in the respective environments.
The cumulative annual yield also reflects that these top yielding genotypes produced the highest annual forage yield and can therefore be selected for stable forage production across environments (Table 8). Furthermore, the observation of relatively high and considerable GCV and PCV values for TFW and TDW in the current study indicates the importance of the genotypic effect in the expression of these traits.

Genotypic Performance for Feed Nutritional Quality
Nutritional quality of forage crops depends mainly on the digestibility and amount of essential nutrients [38]. Results from the feed quality analyses revealed the presence of genotypic variability both in wet and dry seasons, but the GCV and PCV values were low for feed nutritional quality traits in this study, similar to a previous report [39] that showed low genotypic and phenotypic coefficient of variation for quality traits. In the present study, the mean NDF value of studied genotypes was higher than the maximum expected mean NDF value for forage grasses [40,41]; however, the observed mean value for IVOMD and CP contents were higher than the minimum requirement for maintaining ruminants [42,43]. In livestock production, energy is one of the limiting factors in animal performance [44]. ME is the commonly used trait for evaluating energy content of feed [45]. Napier grass genotypes were within the range of ME content for tropical grasses [44]. Genotypes that had the highest ME content were BAGCE 86, CNPGL 96-24-1, and CNPGL 93-08-1, both in wet and dry seasons. Overall, the observed nutritional performance of genotypes indicated that these genotypes could be an important resource for improving feed quality.
Generally, top biomass yielding genotypes had the highest fiber components, i.e., NDF and ADF, while the respective, CP, IVOMD, and ME values were low. In contrast, low forage biomass yielding genotypes had high CP, IVOMD, and ME, with low fiber components. It is also interesting to note that purple Napier grass genotypes, such as CNPGL 93-18-2, CNPGL 93-08-1, and CNPGL 92-133-3, had high CP and ME contents. These findings were consistent with the observed negative correlation between fiber components and CP, IVOMD, and ME. A previous report [46] also indicated that an increase in fiber components reduces cellular nutrients, such as crude protein and digestibility. Furthermore, a PCA partitioned the genotypes into three cluster groups; for example, genotypes in cluster one showed high CP, IVOMD, and ME values, while genotypes in cluster two showed high fiber components and forage biomass yield in both wet and dry seasons. The observed values of feed nutritional quality traits were highly dependent on the season for all Napier grass genotypes in this study. For example, fiber components NDF and ADF were higher in the wet season, while IVOMD, CP, and ME were lower in the wet season. A decline in organic matter digestibility and metabolizable energy during the wet season is associated with increased phenological development and forage production, which stimulates stem production, resulting in a higher stem proportion [47,48]. Functionally, this decline in digestibility is attributed to an increase in structural components and cell maturation [28]. Therefore, the increased concentration of CP, IVOMD, and ME during the dry season could potentially compensate for reduced forage biomass as reflected by relatively lower TDW and higher LSR in the dry season. However, it should be noted that our findings are not in agreement with the results of Reference [49], who reported a decline in CP, IVOMD, and ME values and an increase in NDF, ADF, and ADL in the dry season. These differences might be attributed to the phenological differences of plants at time of harvest.
In line with the scaling relationship between TDW and feed quality traits, the slopes were significantly different from 1 and tended to increase in the dry season compared to the wet season, which might show a flexible relationship between TDW and feed quality traits. The difference in CP, IVOMD, and ME between wet and dry seasons for the same TDW might indicate: (1) an increased structural development, such as stem elongation and reduced leaf area, which would negatively affect these traits; and (2) accurate comparison of genotypes for feed quality traits can be done by considering the effect of plant size [28]. In addition, the rapid decline in the dry season for CP, IVOMD, and ME as compared to the wet season might be attributed to an abiotic stress response due to increased temperature and reduced soil moisture in the dry season that could trigger fiber deposition within the cell wall [50,51]. The allometry based on TDW explained about 30% of the variation; hence, the residual variation would be an important factor to test the significant variation of feed quality traits among genotypes (Table 7). Furthermore, the observed low R 2 value would indicate that prediction based on biomass yield alone might not explain the changes in feed quality traits (Table 7).
Top forage biomass yielding genotypes produced the highest cumulative annual ME yield, indicating increased annual forage biomass yield in high yielding genotypes complements the observed low ME content per kilogram of dry matter. This emphasizes that energy production, per se, coupled with forage biomass production is crucial for characterizing and selecting Napier grass genotypes for livestock/dairy production. Therefore, the efficiency of selection for improved feed quality performance is influenced by how traits are associated. For example, selection only for higher forage biomass yield could compromise CP, IVOMD, and ME. Since no genotypes were entirely high in TDW and ME, exploitation of plant breeding and marker assisted selection (MAS) could be an alternative strategy to develop improved Napier grass for both forage biomass yield and feed quality traits.

Markers Associated with Dry Weight Yields and Metabolizable Energy
Eighty-two markers (SNPs and SilicoDArTs) associated with annual dry weight yield and metabolizable energy were identified, using both the mixed model (MM) and ANOVA statistical model, in marker-trait association analysis. The MM was corrected with pairwise kinship matrix data, while the ANOVA model did not include any correction. As compared to MM, the ANOVA detected a higher number of markers associated with the traits. However, we report here only on markers detected by both models for robust selection. Many of the identified markers detected most of the high performing genotypes identified in the agronomic performance trial. Ten markers (most of them SNPs) were highly discriminant between high and low yielding genotypes. For example, the heterozygous form of one SNP marker (ID number 23617359) identified the top performing accessions for TDW and ME. However, four additional genotypes (BAGCE_93, BAGCE_1, BAGCE_97, and BAGCE_56, all high yielding) also showed the same heterozygous form of the maker. By using a combination of two SNP markers (ID numbers 23617359 and 30283369), two of the accessions (BAGCE_93 and BAGCE_1) were excluded, showing that the use of a combination of markers could improve the diagnostic ability of the markers.
The DArTseq genotyping platform [52] produces short sequence reads, averaging about 55 nt in length, corresponding to each of the SilicoDArT and SNP markers [12]. In this study, the short sequence reads corresponding to the associated markers were compared with previously reported Napier grass genomic [30] and transcriptome [31] sequences. Very few of the DArTseq sequence reads aligned with the assembled genomic sequences, as only assembled sequences longer than 200 nt were used in the comparison. However, many of the DArTseq reads aligned with the transcriptome sequences, indicating that many of the associated markers are found within the gene coding regions. DNA variation in the coding regions is likely to provide a significant contribution to adaptation and productivity [53]. One of these sequences (corresponding with marker 23644354) was annotated as an 'enhanced ethylene response protein', which is one of the genes involved in source-sink communication and sucrose-mediated regulation of starch synthesis [54].
A greater number of associations were detected on chromosome 3, which might indicate the position of quantitative trait locus (QTL) governing biomass yield and metabolizable energy. One of these associated markers was annotated as cytochrome P450. In a previous GWAS analysis in pearl millet [15], cytochrome P450 was associated with plant population, grain number, panicle number, and tiller numbers. Plant cytochrome P450s are involved in a wide range of biosynthetic pathways and play critical roles in the synthesis of lignins and various fatty acid conjugates, hormones, pigments, defense compounds, and signaling molecules [55].
The identified markers could be useful in the implementation of marker assisted selection in Napier grass to improve the efficiency and precision of selecting genotypes for higher dry weight and metabolizable energy. When compared to the field phenotyping and evaluation, the use of the markers is simpler and can save time, resources, and effort, as the selection can be carried out at the seedling stage. In addition, the markers can be important in the identification and mapping of QTLs controlling the traits [56]. In most of the associated markers, the minor allele frequency was associated with higher dry weight yield and metabolizable energy, hence, increasing the frequency of the minor alleles in the population by breeding could improve the productivity of Napier grass with crucial implications for livestock productivity.

Conclusions
Napier grass genotypes showed significant variation for forage biomass yield and feed nutritional quality traits under both wet and dry seasons. Furthermore, genotypic performance was highly influenced by season as reflected by higher forage biomass yield (TFW and TDW) in the wet season than in the dry season. It is interesting to note that top forage biomass producing genotypes under wet season conditions were also top yielding in the dry season. This signifies the inherent genotypic potential of Napier grass genotypes for high biomass production under different environments. Thus, the identification of genotypes with stable high yields across environments is essential for small scale livestock production systems. With regard to cumulative annual forage biomass production, genotypes BAGCE 94, BAGCE 100, and BAGCE 53 are identified as top forage biomass producers for the Bishoftu conditions. Feed nutritional quality traits, such as fiber components, were lower in the dry season, while CP, IVOMD, and ME were higher in the dry season due, at least in part, to a higher LSR. Overall, Napier grass genotypes with higher feed qualities in terms of CP and digestibility were identified suggesting the importance of these genotypes for future utilization. However, the results from the scaling relationship between feed quality traits with forage biomass suggested the need for considering structural changes of the genotypes for accurate feed quality evaluation. In terms of ME, the top forage biomass producing genotypes also had the highest ME yield per year. However, the respective ME per unit dry matter of the high biomass yielding genotypes were low to medium. GWAS analysis identified SNP and SilicoDArT markers that discriminate high biomass yielding Napier grass genotypes. The observation of some identified markers aligning with agronomic performance in pearl millet [47] might be important in implementing MAS in Napier grass to improve the efficiency and precision of selection for high dry matter yields and metabolizable energy. Acknowledgments: The authors would like to thank Tesfaye, Yirsaw Wubete and Fetene Argaw for their contribution to field data collection. The authors would also like to thank EMBRAPA for making their germplasm and breeding lines available for the study through the Africa Brazil Agricultural Innovation Marketplace Project. The authors also would like to thank Linkai Huang for providing genomic sequence of purple Napier grass genotype.

Conflicts of Interest:
The authors declare no conflict of interest.