Crop Breeding and Applied Biotechnology - 21(1): e34402113, 2021 QTL mapping of spike fertility index in bread wheat

: Spike fertility index (SF) is a trait easily measured at maturity and strongly associated with the number of grains per unit area. In order to identify genomic regions involved in SF control, a biparental (Baguette 10 × Klein Chajá) population of 80 recombinant inbred lines (RIL) was used. Seven field trails were conducted to determine the SF BLUP value per RIL. RILs were genotyped using a commercial chip (Axiom® 35K SNP Wheat Breeder’s Array, Affimetrix). A linkage map was constructed with 857 SNP markers, and SF QTL mapping was performed. The narrow-sense heritability of SF was 0.89. Three genomic regions (QTL) associated with SF were found on chromosomes 2D, 4A, and 7A. The proportion of genetic variation explained by these three QTL was 32%, with no significant epistatic interaction between QTL.


INTRODUCTION
Bread wheat (Triticum aestivum L.) is one of the most important crops in the world. Given the current and future scenario of increased global demand for grains, breeding efforts must concentrate on improving grain yield (CIMMYT 2019). The identification of specific and efficient selection criteria, as well as advances in knowledge of the genetic and molecular basis of yield and yield components, will allow an increase in genetic gain.
The spike fertility index (SF), i.e., the number of grains per g of spike chaff [also termed "fruiting efficiency" (Ferrante et al. 2012)], has been widely proposed as a selection criterion in breeding programs (Slafer et al. 2015, Alonso et al. 2018, Fischer and Rebetzke 2018, Valvo et al. 2018, Gerard et al. 2019) mainly due to its association with number of grains per unit area (Abbate et al. 1998, Foulkes et al. 2015, Ferrante et al. 2017).
However, the method of reference for determination of the spike fertility index, first described by Fischer (1984), uses spike dry weight at anthesis, making it a destructive method that is highly sensitive to the exact phenological stage in which measurement is carried out (Fischer and Rebetzke 2018). Abbate et al. (2013) proposed the alternative of the spike fertility index measured at maturity (i.e., calculated with spike chaff weight at maturity), as a selection criterion in breeding programs. This trait has been then shown to have good association with NG (Alonso et al. 2018, Pradhan et al. 2019) and moderate to MP Alonso et al. high heritability (Martino et al. 2015, Alonso et al. 2018, Pretini et al. 2020a, as well as transgressive segregation and a low genotype × environment interaction (Martino et al. 2015, Mirabella et al. 2016, Alonso et al. 2018).
In addition, because it is a simpler, non-destructive method, it is ideal for use as a high-throughput measurement method for selection in early generations. Although there are no conclusive results of how accurate SF at maturity is in estimating SF at anthesis, Slafer et al. (2015) suggest that there are indications of a small overestimation of the SF at maturity, due to the fact that spike dry matter may increase from anthesis to maturity. In this respect, Pretini et al. (2020a) observed instability in this estimator in different environments. However, Abbate et al. (2013) found high association between these two indices (r > 0.7), and Alonso et al. (2018), Fischer and Rebetzke (2018), and Pradhan et al. (2019) detected a positive association between SF at maturity and NG. Furthermore, Terrile et al. (2017) and Alonso et al. (2018) found that selection for high SF resulted in stable and positive genetic gain in grain yield.
Most important agronomic traits in cereals are quantitatively inherited and the genes underlying their variation have been difficult to detect (Neumann et al. 2011). Several studies carried out during the last ~20 years have pursued identification of quantitative trait loci (QTL) of yield and yield-related traits in mapping populations (Börner et al. 2002, McCartney et al. 2005, Kumar et al. 2007, Neumann et al. 2011, Hussain et al. 2017. However, there is little information in the literature about mapping QTL for SF. The first evidence in the literature revealed QTL associated with SF at anthesis based on genome wide association studies (GWAS) in chromosomes 2A, 2D, 4D, 5A, and 7A (Basile et al. 2019, Gerard et al. 2019. As for SF at maturity, Guo et al. (2017) found a QTL on chromosome 2A, with a relatively small effect. Ramirez et al. (2018) reported significant yet small effects of Ppd-B1 and Ppd-D1 on SF at maturity that were independent from those on the flowering date. Recently, a GWAS conducted by Pradhan et al. (2019) allowed identification of 15 markertrait associations (MTAs) with SF in chromosomes 1B, 3B, 4A, 6B, and 7D. Fourteen of these MTAs were located in four regions apparently involved in abiotic and biotic stress response pathways. Recently, Pretini et al. (2020b) validated two QTL for SF in a doubled haploid mapping population, derived from Baguette19 and BioINTA2002. Those QTL were located in chromosomes 3A and 5A. In summary, evidence of the existence of genomic regions significantly associated with SF at maturity is only fragmentary and it has been scarcely addressed through the use of biparental populations.
Nevertheless, extensive work has been carried out in a recombinant inbred line (RIL) population derived from two Argentinian cultivars (Baguette 10 and Klein Chajá) contrasting for SF and other yield-related traits. In this population, the mode of inheritance of SF (Martino et al. 2015, Mirabella et al. 2016, Alonso et al. 2018, the possibility of using SF as a selection criterion in breeding for grain yield (Alonso et al. 2018), and the existence of significant molecular marker-SF associations [in a preliminary study by Panelo et al. (2019)] have been established. Therefore, the objective of this study was to identify QTL associated with SF in this population, evaluated in Balcarce and Marcos Juarez, Argentina.

Plant material
Quantitative trait loci mapping was conducted using a population of 80 RIL from the cross between 'Baguette 10' and 'Klein Chajá', both Argentinian semi-dwarf hard spring wheat varieties released in 2000. 'Baguette 10' (B10; pedigree ARCHE/GENIAL) has been classified as having high SF, whereas 'Klein Chajá' (KCJ; pedigree NINJING/3/BUC'S'//H697/ DKBL) has been classified as having low SF (Martino et al. 2015). Both cultivars have similar intermediate growth cycles and are well adapted to Argentinian wheat-growing areas. They show several differences in spike architecture: B10 has a compact, short, dense spike, with very thin glumes and rachis, whereas KCJ has a longer, laxer spike with a large number of spikelets and grains and a heavy chaff structure (Martino et al. 2015, Alonso et al. 2018).

Phenotypic evaluation
Six field experiments were carried out at the Balcarce (BCE) Experimental Station (lat 37º 45' S; long 55º 18' W, alt 130 m asl) of the Instituto Nacional de Tecnología Agropecuaria (INTA), Buenos Aires Province, Argentina, and one experiment was carried out at the Marcos Juárez (MJ) Experimental Station (lat 32° 43' S; long 62° 06' W, alt 112 m asl), INTA, Córdoba Province, Argentina. General crop management of experiments 1-3 has been described in Alonso et al. (2018); the remaining experiments were conducted similarly except for the fact that no irrigation was applied. Specific information on each experiment is presented in Table 1.
The heading, anthesis, and physiological maturity dates of each plot were registered in the field when 50% of the spikes reached those stages, using the Zadoks scale (Zadoks et al. 1974). Physiological maturity was determined as loss of green from the peduncle. As expected from previous data (Martino et al. 2015, Alonso et al. 2018, each of the phenological stages was concentrated within one week for >90% of the RILs and within 10 days for the entire population (data not shown). Weather conditions were recorded daily with a standard meteorological station located at each experimental station.
At maturity, 15-20 spikes were drawn at random from each plot. When seven-row plots were used, only the five central rows were sampled. Spikes were cut at the lowest spikelet level, counted, weighed, and threshed. Spike chaff dry weight (g) was calculated as the difference between total spike dry weight (i.e., before threshing) and total grain weight. Grains were counted using an electronic counter. Then, SF (grains g -1 ) was calculated as the quotient between number of grains and spike chaff dry weight (Abbate et al. 2013).

Linkage map construction
Eighty RIL and the parents of the population were genotyped with the Axiom® 35K SNP Wheat Breeder's Array (Affimetrix) (Allen et al. 2017). For genotyping, genomic DNA was extracted from a single seedling leaf tissue according to Haymes (1996). Only those polymorphic SNPs showing less than 10% missing data and segregation distortion under 20% were considered for construction of the linkage map. In addition, independent genotyping of the population with two functional markers for the Vrn-A1 (Yan et al. 2004) and Rht-D1 (Ellis et al. 2002) genes was added to the analysis. Those markers that had identical or compatible segregation in the whole population were grouped using the merger. pysoftware (https://github.com/juancrescente/lmap) in order to minimize the redundant information. Once the markers were grouped, the linkage map was constructed using the R/qtl software (Broman et al. 2003). Genetic distances between the markers were calculated based on the Haldane mapping function. Throughout the 21 wheat chromosomes, several SNPs were anchored to the Ref Seq 1.0 genome assembly (Appels et al. 2018). Relationships between the genetic and physical positions of the SNPs were then established.

Statistical analysis and QTL detection
Linear mixed models were fitted for SF using the lme function from the nlme package (Pinheiro et al. 2013) of the R software (R Core Team 2015). The models included replications within environments (years and/or sowing dates within a year) and environments as fixed factors, and genotypes and the genotype × environment interaction as random factors. The critical level of significance used was 0.05. Variance components and narrow-sense heritability were estimated according to Alonso et al. (2018).
Quantitative trait loci analysis was conducted with Composite Interval Mapping (CIM) using QTL Cartographer software (Wang et al. 2012a). Best linear unbiased predictors (BLUP) for each RIL, obtained from the mixed model, were used in QTL analysis. Threshold was calculated with 500 permutations and a 0.05 critical level of significance. Up to ten markers showing the highest F value after the forward-backward stepwise regression analysis were added as cofactors in the CIM step [model 6, using a moving window size of 10 centiMorgan (cM) and a walking speed of 1 cM]. The most likely position of the QTL was determined as the point with the maximum logarithm of the odds (LOD) score. The confidence interval (CI) of each QTL was defined as the map interval corresponding to a LOD-2 decrease to each side of the LOD The physical position of QTL was considered to be that of the marker nearest the peak LOD score. Using the flanking sequence for each SNP marker, provided by the chip manufacturer (Axiom, Affimetrix), a local alignment was performed using the BLAST algorithm (Altschul et al. 1990) against the reference sequence IWGSC RefSeq v1.0 of the bread wheat genome (Appels et al. 2018) to verify its position.
Haplotype was constructed for each RIL using the marker associated with the peak of the maximum LOD-score of each QTL identified. Haplotype proportions were tested with the Chi-square test (p value = 0.05, df = 7). Due to heteroscedasticity between haplotype groups, Welch's ANOVA test was carried out using userfriendlyscience package (Verboon et al. 2018) of the R software, and haplotype differences were estimated using the Games-Howell nonparametric test.

RESULTS AND DISCUSSION
Mean SF in the RIL population was 98.3, 91.9, 89.5, 92.4, 89.5, 90.7, and 101.7 grains g -1 for experiments 1 to 7, respectively. In each environment, SF values showed a symmetric, bell-shaped distribution with a wide range of variation from ~55 to ~135 grains g -1 [partial data is available in Alonso et al. (2018) and Figure S1], whereas the parents showed SF values as expected from previously published data (Alonso et al. 2018 and Figure S1). The relative weight of variance components and heritability were also coincident with previous reports, as genetic, genetic × environmental, and residual variances were 54.1, 9.3, and 67.2, respectively, and heritability was 0.89. These results are in line with previous studies, which showed that SF fit an oligogenic model comprising a few loci with relatively large effects and high heritability (Martino et al. 2015, Mirabella et al. 2016, Alonso et al. 2018; i.e., a high probability of a finite number of markers explaining a considerable proportion of the phenotypic variance. In this study we included data from seven environments, and SF was accurately determined and predicted with BLUPs, which increases the reliability of the QTL significance (Piepho et al. 2008, Segura et al. 2009, Sadok et al. 2013).
This study was carried out with genotypic data from 80 individuals. The genetic linkage map consisted of 368 loci on the 21 chromosomes of bread wheat and spanned 3674 cM (Table S1). These loci included 857 SNP and two functional markers corresponding to the Rht-D1 and Vrn-A1 genes. Linkage maps per chromosome are shown in Figure S2. The fact that the functional markers Rht-D1 and Vrn-A1 were correctly mapped in the genome validates the mapping procedure. The use of mapping populations of reduced size, as was the case in the present study, may allow detection of QTL with major effects, but limits the detection of additional, small, yet real QTL (Beavis 1998, Vales et al. 2005, Cavalcanti et al. 2012, Wang et al. 2012b. Schön et al. (2004) recommend the use of a conservative threshold, such as the one used here, if the aim of a study is to identify a few large QTL controlling a limited proportion of the genetic variance. On the other hand, the number of SNP markers on each chromosome does not allow restriction of QTL to smaller map distances (Li et al. 2010).
The elucidation of the molecular and genetic basis of yield-related traits can contribute not only to understanding how yield is determined, but also to the development of technologies for speeding up the selection process, leading to high yielding cultivars. In this study, we detected three QTL associated with SF, a trait closely linked to number of grains per unit area, which is the main yield component in bread wheat. Three QTL for SF were detected in chromosomes 2D, 4A, and 7A (Table 2); the latter two have not been previously reported. All three QTL showed high stability across seven environments, which spanned differences in temperature, water regime, solar radiation, and photoperiod. The proportion of genetic variation explained by these three QTL was 32%, with no significant epistatic interaction between QTL.
The QTL on chromosome 2D (referred to as Q sf.bfe.2DL ) had an additive effect of 7.44 grains g -1 , expressed as the difference from the population mean. The donor of the SF-increasing allele was B10 (i.e., the parent with higher SF). The confidence interval of this QTL was 9 cM, and the physical position of the SNP marker associated with the peak of the maximum LOD score was ~648 Mbp (Table 2).
On chromosome 4A, the QTL Q sf.bfe.4AL showed an additive effect of 2.44 grains g -1 , with B10 as the donor parent. Located at ~472 Mbp, the confidence interval of Q sf.bfe.4AL comprised 5.4 cM. A third QTL, Q sf.bfe.7A , was detected in chromosome 7A donated by B10. Q sf.bfe.7A showed an additive effect of 1.74 grains g -1 , and a greater confidence interval (15.1 cM) in comparison with the remaining two QTL. Regarding its physical position, the SNP marker associated with the peak of maximum LOD was located in ~119 Mbp. Basile et al. (2019) also detected SF QTL in the latter two regions. One of the QTL was located in the vicinity, but not in the same position, of the newly reported Q sf.bfe.4AL . Indeed, the QTL reported by Basile et al. (2019)  MP Alonso et al. to colocalize with Q sf.bfe.7A , detected in this study. Guo et al. (2017) proposed two candidate genes for SF and other associated traits, named CONSTANS4 (CO4) and Six-rowed spike 1 (Vrs1), located in chromosome 2AL. Vrs1, also known as GNI-A1, has recently been associated with the number of fertile flowers and grains per spikelet traits (Sakuma et al. 2019. However, considering the genetic positions of markers associated with CO4, this QTL would not be homoeologous to Q sf.bce.2DL . The genetic position of Vrs1 has not yet been defined in the IWGSC RefSeq v1.0 reference sequence of the bread wheat genome (Appels et al. 2018). As for the QTL mapped by Gerard et al. (2019) on chromosome 2D, the physical position of the associated marker (Kukri_rep_c68068_95) is 641.1 Mb, which is very close to AX-94501170 and AX-95232269 (648 Mb), the SNP marker associated with the peak of maximum LOD of Q sf.bce.2DL detected in the present study. Gerard et al. (2019) measured SF (fruiting efficiency -FE) at anthesis, not at maturity, as it was done in this study. Abbate et al. (2013) have shown the existence of differences between FE at anthesis and at maturity. Therefore, the QTL detected by Gerard et al. (2019) may not colocalize with Q sf.bce.2DL . None of the QTL found by Pretini et al. (2020b) in chromosomes 3A (685.12 Mbp) and 5A (461.49 Mbp) were present in the B10xKCJ RIL population.
The effect of the Ppd genes, reported by Ramirez et al. (2018) as being associated with SF, could not be evaluated in this population because their parents, B10 and KCJ, are monomorphic for these genes (Vanzetti et al. 2013). a Position of the LOD score peak in centiMorgan (cM). b LOD score value at the LOD score peak. c LOD-2 confidence intervals in cM. d AE: Additive Effect: Difference from the population mean. e SNP marker associated with the LOD score peak.  The results of haplotype analysis are shown in Table 3. The haplotype proportions corresponded with the expected values (χ 2 = 0.71). The best performing group, BBB (i.e., the one with B10 alleles at all three QTL), showed an average increase of 4.9 grains g -1 in SF with respect to the overall mean, whereas group AAB (i.e., with two KCJ alleles and one B10 allele) decreased average SF by 4.4 grains g -1 compared with the overall mean. The haplotype group BBB, as well as all three groups carrying at least two B alleles, showed significantly higher SF increases than did group AAB. The AAA group was significantly different only from BBB. These results show the positive effect of B10 alleles on SF compared to the presence of KCJ alleles.
The present study provides information on genomic regions controlling SF in bread wheat. Peak SNP markers (Table  2) may be used to develop fine mapping populations in order to detect candidate genes which control the trait and to design strategies of marker-assisted selection for a complex quantitative trait such as grain yield.. Table 3. Average increase in spike fertility index (SF, grains g -1 ) with respect to the population mean as analyzed by the haplotype group at three QTL