Maize inbred lines “18237” (recurrent parent and highly susceptible to GSR, P1) and “S72365” (donor parent and completely resistant to GSR, P2) were crossed to produce the F₁ hybrid, which was self-pollinated to generate the F₂ population, and backcrossed with “18237” to generate the BC₁F₁ population. These two populations were grown at the experiment farm of the Sichuan Academy of Agricultural Sciences (Chengdu, Sichuan, China). Each individual was artificially inoculated with F. graminearum.

Fusarium graminearum was cultured on potato dextrose agar in darkness at 25°C for 4–5 days. The maize kernel was prepared by first dipping in water for 20 h at 37°C in darkness, then in boiling water for 1 h. The kernels were then dried on a ventilated table, and autoclaved for 20 min at 121°C within plastic bags. Preparing for field evaluation, the sterilized kernels were inoculated with F. graminearum at 25°C in complete darkness for 15 days. Field inoculation of plants was conducted as described by Yang et al. (2010).

Plants were evaluated for stalk rot symptoms twice a week, beginning 1 month post-inoculation. Typical symptoms of stalk rot were observed, such as browning reactions in lower internodes, spongy stem, wilting, lodging, and plant death. Evaluating mycelial growth and root damage requires the stem to be cut. Incidents of stalk rot infection was scored using a disease assessment scale of 1–9. Scales 1–2 were regarded as resistant and 8–9 were regarded as susceptible. Plants with a score of 9 were dead and lodging with broken vascular tissue of the stem; a score of 8 was similar to 9, the plant lodging down but with an unbroken stem; plants with a score of 7 exhibited withered leaves and a soft stem, but no lodging; a score of 6 corresponded to symptoms of withered leaves, but with a harder stem than in plants with a score of 7; plants with a score of 5 exhibited withered leaves, and a slightly soft stem; a score of 4 was assigned when parts of leaves were withered, and a normal stem was observed; a score of 3 was given for the observation of only leaf chlorosis; a score of 2 indicated some yellow leaves; and a score of 1 indicated no obvious symptoms.

DNA was extracted from fresh young leaves of single plants using the standard CTAB protocol (Doyle and Doyle, 1990). For QTL-seq, two DNA pools, susceptible pool (S-pool) and resistant pool (R-pool) were constructed, respectively, by mixing an equal amount of DNA from 25 F₂ plants with the lowest disease scores and 25 F₂ plants with highest disease scores (Figure 1A). A DNA concentration of 2–5 μg from each of the P1, P2, R-pool, and S-pool were used to construct pair-end sequencing libraries (150 bp read length, which were sequenced using Illumina HiSeq 2500 (Illumina Inc., San Diego, CA, United States) by Gene Denovo Biotechnology Co. (Guangzhou, China). Raw reads with >10% unidentified nucleotides and with >50% bases having phred quality scores of <20 were filtered out to get high-quality clean reads. To identify SNPs, these clean reads were mapped and aligned to the maize reference genome (RefGen_V4¹) using the Burrows–Wheeler Aligner (BWA) software (Li and Durbin, 2009) with the settings as follow: mem 4 -k 32 -M. SNP-calling was performed for all samples using the SAM tools (Li and Durbin, 2009). The SNP positions with a read depth <6 and SNP-index <3 were filtered out. To confirm the physical positions of each SNP, the software ANNOVAR (Wang et al., 2010) was used to align and annotate SNPs.

In this study, the parameters of SNP-index and Δ (SNP-index) (Abe et al., 2012; Takagi et al., 2013a) were calculated to identify candidate regions for maize stalk rot resistant QTLs. The reference sequence for SNP-index calculation was developed by replacing the detected SNPs from one of the parental cultivars with those from the reference genome. The SNP-index represents frequencies of parental alleles in the population of pooled individuals. Slide window analyses with parameters “2 Mb windows size and 100 kb increment” was applied to SNP-index plots.

The Δ (SNP-index) was calculated based on subtraction of SNP-index between R-pool and S-pool. SNP-index is equal to “0” or “1” when entire reads contain genomic fragments from P1 or P2, respectively. The Δ (SNP-index) value will be significantly different from 0 in genomic regions with major QTL of the target gene (Takagi et al., 2013a). We calculated statistical confidence intervals of Δ (SNP-index) for all the SNP positions with given read depths under the null hypothesis of no QTLs, and plotted them along with Δ (SNP-index). With a 95% confidence interval in 10,000 bootstrap replicates, the Δ (SNP-index) was obtained for each read depth.

To verify the major QTL for GSR resistance from QTL-seq, polymorphic markers were developed in the predicted region of maize chromosome 8. SSR markers in the predicted region were searched using the MISA tool², and employed for polymorphism screening between two parental lines, and between the R-pool and S-pool. SNP markers that were in the predicted region of the QTL were converted to PCR-based markers, and the primers were designed using Primer 5³. Polymorphic markers were used to screen the F₂ population. The linkage analysis was performed using the software JoinMap 4.1 (Van Ooijen, 2006) and recombination values were converted to centiMorgan (cM) using the Kosambi mapping function (Kosambi, 1943). The genetic information together with phenotyping data was used for QTL analysis using the composite interval mapping (CIM) model in WinQTL cartographer 2.5 software (Wang et al., 2012).

The expression of candidate genes was investigated using real-time quantitative PCR (qPCR). Leaf samples were collected from P1, P2, F₁, F₂-S (susceptible to GSR) and F₂-R (resistant to GSR) individuals at the early stage of inoculated plants. Total RNA for all samples was extracted using Trizol Reagent (Invitrogen, Carlsbad, CA, United States) and digested with RNase-free DNase I (Takara Bio, Japan) for 30 min at 37°C. Reverse transcription was conducted by Super III Reverse Transcriptase (Invitrogen, San Diego, CA, United States). The qPCR primers for the candidate genes (Zm00001d011953 and Zm00001d011972) were 5′-CCAGCTGTACAGGAGCATGA-3′ (forward) and 5′-CCGGAACACGTCTTGGTAGT-3′ (reverse) for Zm00001d011953, 5′-AAAAGGCTTGTTGCTGGAGA-3′ (forward) and 5′-GGTGGAGGTGCATTTTGTCT-3′ for Zm00001d011972, respectively. qPCR was performed in a LightCycler^® 96 Real-time PCR Instrument (Roche, Swiss) with SYBR Green Real-time PCR Master Mix (Takara, Japan). The gene expression levels were determined using Ct value normalized with the formula 2^-ΔΔCt (Livak and Schmittgen, 2001). The maize Actin gene was employed as an endogenous control, with the following primers: 5′-GCCGGTTTCGCTGGTGATGATGCGCC-3′ (forward) and 5′-GTGATCTCCTTGCTCATACGATCGGC-3′. Three replicates were measured to calculate the average relative expression levels. A Student’s t-test was used to check the significant differences in expression levels among these five samples.

Frequency distribution of resistance to GSR is presented in Figures 1A–E. Following artificial inoculation with F. graminearum, the inbred line “18327” (P1) plants showed severe stalk rot symptoms (Figure 1B); “S72356” (P2) exhibited complete resistance to GSR and no symptoms were observed (Figure 1C). Most of the F₁ hybrids (85.3%) displayed high levels of resistance to GSR (Figure 1D), suggesting that the major GSR resistance allele might be dominant. The resistance to GSR in the F₂ population showed continuous variation (Figure 1A), and a skewed distribution of disease severity was observed in BC₁F₁ population (Figure 1E). Based on these results, it was suggested that the resistance to GSR in P2 was likely to be a quantitatively inherited trait.

Based on library construction and NGS-based high-throughput sequencing of P1, P2 and two DNA-pools, a total of 344 Gb of data was generated, including 2.93 billion of 150 bp high-quality clean reads, and 98.44–98.56% high-quality reads were mapped on the reference genome. The average sequence depths were 20-fold in parents and 30-fold in pools. The total number of variants was 16, 997, 640, including 15, 490, 449 SNPs and 1, 507, 191 indels. The Q20 ratio ranged from 95.66 to 96.21% (Table 1).

To identify the candidate genomic region conferring resistance to GSR, the SNP-index was calculated, based on each SNP identified. The average SNP-index was calculated with a sliding window of 2 Mb intervals with 100 kb increment for S-pool and R-pool to detect the candidate genomic regions. SNP-index graphs were generated for R-pool (Figure 2A) and S-pool (Figure 2B) by plotting the average SNP-index against the position of each sliding window in the P1 genome assembly. It was expected that the SNP-index graphs of the R-pool and S-pool would be identical for the genomic regions that are not relevant to phenotypic difference, whereas the genomic region harboring the GSR resistance QTL should exhibit unequal contribution from P1 and P2 parental genomes. In addition, the SNP-index of predicted regions for R-pool and S-pool would appear as mirror images (Takagi et al., 2013a). After calculating and combining the information of SNP-index in R-pool and S-pool, Δ (SNP-index) was calculated and plotted against the genome positions (Figure 2C).

In the present study, the region of chromosome 8 ranging from 161.001 to 170.6 Mb had an average SNP-index higher than 0.60 in R-pool with the highest equal to 0.80. Conversely, the average SNP-index in the region of S-pool was lower than 0.45 with the lowest equal to 0.3. The predicted genomic region harboring the resistance QTL to GSR identified by QTL-seq were determined by Δ (SNP-index) value. The threshold of the Δ (SNP-index) value was 0.25 at the 95% significance level. Results showed that the genomic region on the chromosome 8 from 161.001 to 170.6 Mb was significantly different from 0 (Figure 2C). These data demonstrated that in maize, a major QTL conferring GSR resistance was present in the 161.001–170.6 Mb region on chromosome 8. We named this region as Rgsr8.1.

Based on 565 non-synonymous variations, a total of 45 SNP markers (Supplementary Table S1) were developed, which distributed equally over the predicted region according to physical position. A total of 729 SSR markers were searched from the predicted region. Among these SSR markers, 165 SSR markers (Supplementary Table S2), equally distributed on the predicted region, were used to analyze the polymorphism.

The 45 SNP markers and 165 SSR markers were checked for polymorphisms between P1 and P2, R-pool, and S-pool. Of the SNP markers, 33 markers amplified well, and 12 SNP markers were found polymorphic between P1 and P2, R-pool, and S-pool (Table 2). Twenty-nine SSR markers were identified to be polymorphic between P1 and P2, R-pool, and S-pool (Table 3). In total, 12 SNP markers and 29 SSR markers were used for QTL analysis based on the F₂ populations. A major QTL for resistance to GSR, physically located in the region of 164.678–166.721 Mb on chromosome 8 (Figure 3), was flanked by two SNP markers (SNP-18 and SNP-25) with genetic distances of 4.57 and 6.62 cM. This result agreed with the QTL-seq analysis supporting a major GSR resistance QTL on chromosome 8. The LOD scores of the polymorphism markers within this region ranged from 0.26 to 45.23, and could explain 34.4% of the variance (Table 4). Additionally, we further narrowed down the Rgsr8.1 locus by using recombination test, based on 6 BC₁F₁ recombinants, which were recovered within the region on chromosome 8. To figure out the physical position where the recombinant events occurred, eight markers (Table 4) were used to analyze the P1, P2, and recombinants. The results showed that no recombinants were detected except for SSR-65 and SNP-25. Therefore, the mapping data narrowed the Rgsr8.1 locus down to a 2.04 Mb interval between the SSR-65 and SNP-25 (Figure 4). Furthermore, we used these eight markers to screen the F₂ and BC₁F₁ population. We estimated that the SSR-78 marker was linked with GSR resistance in “S72356” via phenotypic and genotypic identification.

Based on the Maize reference genome (RefGen_V4, see text footnote 1), 33 genes were located in the predicted region (Supplementary Table S3). Based on gene annotation of the region, two genes, Zm00001d011953 (Zm953) and Zm00001d011972 (Zm972), were chosen as candidate genes. Zm953 encodes an auxin response factor, and Zm972 encodes a disease resistance protein. The physical location of Zm953 (164, 991, 768) and Zm972 (165, 428, 843) were near to SSR-78, which had the highest LOD value (Table 4).

In addition, the expression levels of Zm953 and Zm972 were investigated in two parental lines, the F₁ plant, and the susceptible and resistant individuals in F₂ by using qPCR. The results showed that the expression level of Zm953 in the resistant plants from P2, F₁, and F₂-R was significantly higher than that in the susceptible plants from P1 and F₂-S (Figure 5). The expression level of Zm972 showed a similar result (Figure 5). Therefore, based on the gene annotation and results of expression analysis, we inferred that Zm953 and Zm972 may be the candidate genes for the major QTL conferring GSR resistance.