Identification of quantitative trait loci underlying resistance to southern root-knot and reniform nematodes in soybean accession PI 567516C

Soybean cyst nematode (SCN, Heterodera glycine Ichinohe), southern root-knot nematode [SRKN, Meloidogyne incognita (Kofoid and White) Chitwood] and reniform nematode (RN, Rotylenchulus reniformis Linford and Oliveira) are three important plant–parasitic pests in soybean. Previous study showed that plant introduction (PI) 567516C harbored novel quantitative trait loci (QTL) conferring SCN resistance to soybean. However, QTL underlying resistance to SRKN and RN in PI 567516C remain unknown. The objectives of this study were to identify QTL for resistance to SRKN and RN in PI 567516C. Two hundred and forty-seven F6:9 recombinant inbred lines, derived from a cross between cultivar Magellan and PI 567516C, were evaluated for resistance to SRKN and RN. Two hundred and thirty-eight simple sequence repeats and 687 single nucleotide polymorphism markers were used to construct a genetic linkage map. Three significant QTL associated with resistance to SRKN were mapped on chromosomes (Chrs.) 10, 13 and 17. Two significant QTL associated with resistance to RN were detected on Chrs. 11 and 18. Whole-genome resequencing revealed that there might be Peking-type Rhg1 in PI 567516C. Our study provides useful information to employ PI 567516C in soybean breeding in order to develop new cultivars with resistance to multiple nematodes.

Considerable efforts have been made to identify sources of resistance to SRKN or RN in soybean. Luzzi et al. (1987) evaluated 2370 soybean accessions in the USDA-ARS germplasm collection from maturity groups (MG) V, VI, VII and VIII for their resistances to SRKN. Four soybean genotypes, 'Amredo', PI 96354, PI 408088 and PI 417444, were found to be highly resistant to SRKN (Luzzi et al. 1987). Hussey et al. (1991) and Harris et al. (2003) evaluated 139 soybean cultivars and 608 soybean accessions for resistance to SRKN, respectively. Thirty-nine soybean cultivars and seven soybean accessions were found to be resistant to SRKN (Hussey et al. 1991;Harris et al. 2003). For sources of resistance to RN, Pickett and Dyer were the first two soybean cultivars which were reported to be resistant to RN (Rebois et al. 1968). After that, more sources of resistance to RN had been identified in soybean (Robbins et al. 1994a, b;Davis et al. 1996;Robbins and Rakes 1996). It was showed that Peking-derived SCN-resistant cultivars were resistant to RN and PI 88788-derived SCN-resistant cultivars were not (Rebois et al. 1970;Robbins et al. 1994a, b;Davis et al. 1996;Robbins and Rakes 1996). Peking and PI 88788 were two major sources of resistance to SCN that have been widely deployed in soybean breeding ). These two soybean accessions had different functional alleles of same QTL, Rhg1, for resistance to SCN Meksem et al. 2001;Brucker et al. 2005;Cook et al. 2012Cook et al. , 2014.
In comparison with about twenty-five QTL for resistance to SCN Wang et al. 2004;Kabelka et al. 2005;Guo et al. 2005Guo et al. , 2006Winter et al. 2007;Wu et al. 2009;Vuong et al. 2010), fewer QTL for resistance to SRKN or RN were mapped in soybean. Tamulonis et al. (1997) mapped two QTL with a major one on Chr. 10 and a minor one on Chr. 18 for resistance to SRKN in PI 96354. These two QTL were also confirmed in another study (Li et al. 2001). The major QTL on Chr. 10 was detected across 27 soybean cultivars with resistance to SRKN, indicating that this QTL had been widely deployed in soybean breeding (Ha et al. 2004). Recently, two independent studies reported the identification of candidate genes underlying this QTL from two different sources of resistance to SRKN, PI 96354 (Pham et al. 2013) and PI 438489B (Xu et al. 2013), respectively. Besides QTL on Chrs. 10 and 18, two other QTL were also reported to be associated with resistance to SRKN in soybean. One QTL was mapped on Chr. 7 and was responsible for resistance to SRKN race 2 in a soybean cultivar LS 5995 (Fourie et al. 2008). The other QTL was mapped on Chr. 6, which might account for a portion of genetic variation that the QTL on Chrs. 10 and 18 could not explain (Tamulonis et al. 1997;Li et al. 2001;Shearin et al. 2009).
Three major QTL responsible for resistance to RN were detected and genetically mapped on Chrs. 19, 11 and 18, respectively, in a soybean accession PI 437654 (Ha et al. 2007). Despite the fact that QTL on Chr. 19 was the largest one, combination of two QTL on Chrs. 11 and 18 was necessary and sufficient to confer high level of resistance to RN (Ha et al. 2007). These two QTLs were overlapped with the QTL for resistance to SCN (Wu et al. 2009), which provided evidences that there might be common genes responsible for resistance to SCN and RN in soybean.
Soybean accession PI 567516C was genetically unique from most sources of resistance to SCN (Arelli et al. 1997;Chen et al. 2006;Arelli et al. 2009;Vuong et al. 2010). QTL mapping for this accession identified two QTL, one on Chr. 10 and the other on Chr. 18, which were associated with resistance to SCN HG types 2.5.7, 0, 2.7, 1.3.5.6.7 and LY1 (Vuong et al. 2010). Thus, PI 567516C has great potential in breeding new soybean cultivars with resistance to SCN to fight against genetic shift of SCN populations due to long-term use of Rhg1 and Rhg4, two important genes for resistance to SCN. With the aim of identifying sources of resistance to multiple nematode species, we evaluated a subset of soybean accessions with resistance to SCN for resistance to SRKN and RN. The results showed that PI 567516C not only displayed resistance to SCN but also had high level of resistance to SRKN and RN, which further made this soybean line a good source of resistance to multiple nematode species in soybean breeding. However, QTL underlying resistance to SRKN and RN in PI 567516C are still unknown, hampering the utilization of this soybean accession in developing soybean cultivars with resistance to multiple nematode species through marker-assisted selection. The objectives of this study were to identify QTL responsible for resistance to SRKN and RN in PI 567516C.

Genetic population
A genetic population of 247 F 6:9 recombinant inbred lines (RILs) were advanced by single seed descent method from a cross between cultivar Magellan (Schapaugh et al. 1998) and PI 567516C (Arelli et al. 1997), in 2005 at the Bradford Research and Extension Center (BREC), University of Missouri-Columbia, MO, USA. Genomic DNA was extracted from a pooled sample of leaves from five plants of each RIL following a protocol as previously described (Vuong et al. 2010).

Southern root-knot and reniform nematode bioassays
Two hundred and forty-seven RILs and their parents, cv. Magellan and PI 567516C, were evaluated for resistance to southern root-knot nematode [SRKN, M. incognita (Kofoid and White) Chitwood] in the greenhouse facility at the University of Georgia, Athens, GA. A randomized complete block design (RCBD) with three replications was used for SRKN bioassays in this study. Each plant was inoculated with approximate 3000 SRKN eggs as described by Li et al. (2001). The SRKN galls that developed on each plant were counted at 36 days after inoculation.
The evaluation of the RIL mapping population for resistance to reniform nematode (RN, R. reniformis Linford and Oliveira) in the greenhouse follows a protocol as previously described (Robbins et al. 1999). A reproductive index (RI) was defined as the number of eggs ? vermiform nematodes at test termination (Pf)/initial infestation level (Pi). This value was calculated from five replicates for each RIL and used in statistical analysis and QTL mapping.

Statistical analysis
The SRKN gall numbers and RI among RILs were tested for normality using the PROC UNIVARIATE procedure of SAS 9.3 (SAS institute, Gary, NY, USA). A broad-sense heritability was calculated following the method described by Wu et al. (2009).

Linkage analysis and genetic mapping
The fluorescently labeled simple sequence repeats (SSR) markers and the universal soybean linkage panel 1.0 (the USLP 1.0) containing 1536 single nucleotide polymorphism (SNP) loci (Hyten et al. 2008) were utilized to genotype the mapping population using the Illumina GoldenGate assay (Fan et al. 2006;Hyten et al. 2010). The protocol for SSR marker analysis was described by Vuong et al. (2010). Briefly, polymerase chain reaction (PCR) was conducted with a final volume of 12.5 ll on the Eppendorf 96-well master cycler gradient (Eppendorf AG, Germany). Each PCR included 40-50 ng of template DNA, 0.13 lM of labeled forward primer (Applied Biosystems, Foster City, CA, USA) and 0.2 lM of reverse primer (IDT Inc., Coralville, IA, USA), 19 reaction buffer (20 mM Tris-HCl, pH 8.0, 50 mM KCl), 2.5 mM MgCl2, 0.2 mM of each of the dNTPs and 1 unit of Taq DNA polymerase (GenScript Corp., Piscataway, NJ, USA). PCR was performed at 95°C for 5 min, followed by 35 cycles of 94°C for 30 s, 47°C or 52°C for 45 s and 72°C for 1 min, with a final extension for 7 min at 72°C. After purified with the Whatman PCR cleanup procedure (Whatman Inc., Piscataway, NJ, USA), PCR products were analyzed with the ABI 3100 or 3730 DNA sequencer (Applied Biosystems, Foster City, CA, USA) and the GeneMapper 3.7 program (Applied Biosystems, Foster City, CA, USA). A genetic linkage map was constructed using JoinMap 4.0 (van Ooijen 2006). A likelihood of odds (LOD) threshold score of 3.0 and a maximum genetic distance of 50 cM were used for the initial linkage grouping of markers. The soybean genetic linkage groups (LGs) (Song et al. 2004) were replaced with the new assignments of corresponding chromosome numbers (Chr.) (Grant et al. 2010).
QTL analyses were performed using the multi-QTL method (MQM) with the program MapQTL 5.0 and the appropriate cofactor (van Ooijen 2004). The threshold of LOD score was chosen to be 3.4 to declare a putative significant QTL based on the results of permutation test with 1000 runs to determine the P = 0.05 genome-wide significance level. The proportion of the phenotypic variance explained by the QTL effects was estimated at the QTL peaks. Additive (A) effects of significant QTL were estimated from an output of the program MapQTL 5.0. The chromosomes with LOD plots were created using the MapChart 2.2 program (Voorrips 2002).
Whole-genome resequencing, SNP calling and copy number variation analysis Whole-genome sequencing of cv. Magellan, Peking, PI 437654, PI 424298, cv. Hutcheson, cv. Essex, PI 567516C, PI 567387, PI 209332, PI88788, PI 567305 and PI 404198B was conducted in Beijing Genome Institute (BGI), Shenzhen, China, using Illumina technology following a described protocol (Xu et al. 2013). The sequencing depth for each sample was about 159 coverage. SNP calling was conducted by the use of SAMtools (Li et al. 2009). Copy number variation (CNV) analysis was conducted using CNVseq software (Xie and Tammi 2009). Sequence alignment was visualized by IGV software (Thorvaldsdottir et al. 2013;Robinson et al. 2011).

Phenotypic variation and genetic map
Cultivar Magellan and PI 567516C showed significant difference in resistance to SRKN and RN ( Table 1). The frequency distribution of SRKN gall number as well as RN reproductive index (RI) among 247 F 6:9 RILs displayed large genetic variation, ranging from a low value of 0 to a high value of 156.3 and from 0.4 to 15.5, respectively (Table 1; Fig. 1). The normality test by the Shapiro-Wilk (w) statistic indicated that the gall number and RI were not normally distributed (Table 1; Fig. 1). The values of broad-sense heritability of gall number and RI were 0.75 and 0.32, respectively (Table 1).
Two hundred and thirty-eight polymorphic SSR markers as well as 1536 SNP array were utilized for genotyping 247 RILs. Among 1536 SNP loci, six hundred and eighty-seven were found to be polymorphic between two parents. Finally, two hundred and nineteen SSR and 681 SNP markers were used to construct a linkage map, spanning approximately 2816 cM with an average distance of 3.1 cM across 20 chromosomes (Supplemental Figure 1).

QTL associated with resistance to SRKN and RN
Three QTL were identified to be associated with resistance to SRKN in PI 567516C (Table 2; Fig. 2). The major one with a peak LOD value of 23.3 was  Fig. 2). The closest marker to the peak was Satt487 (Table 2). This QTL was located in a similar region to the ones in PI 96354 and PI 438489B (Tamulonis et al. 1997;Li et al. 2001;Xu et al. 2013), which indicated there might be the same QTL among these three resistant accessions. In PI 96354, four candidate genes, Glyma 10g02090 (Glyma.10g016600), Glyma02100 (Glyma. 10g016700), Glyma10g02140 (Glyma.10g017100) and Glyma10g02160 (Glyma.10g017200), were speculated to be responsible for resistance to SRKN (Pham et al. 2013). In PI 438489B, three candidate genes, Glyma10g02140 (Glyma.10g017100), Glyma10g02 150 and Glyma10g02160 (Glyma.10g017200), were identified in this QTL region (Xu et al. 2013 Table 1). These five accessions had not been reported to be resistant to SRKN before. Whole-genome resequencing was performed for these resistant accessions plus six susceptible ones (Supplemental Table 1). The results showed there were no CNVs of Glyma10g02140 (Glyma.10g 017100) and Glyma10g02160 (Glyma.10g017200) among all tested soybean lines. However, 44 SNP/ Indels in Glyma10g02140 were commonly found in all resistant PIs in comparison with Williams 82 (Supplemental Table 1). Among these SNP/Indels, six SNPs in the first exon were non-synonymous, one of which resulted in a stop codon (Supplemental Table 1). These six SNPs were identical to the ones in PI 96354 (Pham et al. 2013). For Glyma10g02160, there were 11 SNP/Indels common in all the resistant PIs (Supplemental Table 1). Among these, two were located in the exon. However, neither of them results in amino acid changes. The other two QTL were mapped on Chrs. 13 and 17, which explained 7.3 % and 5.1 % of total phenotypic variation, respectively (Table 2; Fig. 2). QTL on Chr. 13 was previously mapped in PI 438489B (Xu et al. 2013). QTL on Chr. 17 was a novel QTL for resistance to SRKN, which has not been mapped before.
Characterization of Rhg1 in PI 567516C through whole-genome resequencing method It has been reported that there are two different Rhg1 types, PI 88788-type and Peking-type Rhg1 (Brucker et al. 2005;Meksem et al. 2001). Recently, copy number and sequence variation for three genes, Gly-ma18g02580 (Glyma.18g022400), Glyma18g02590 (Glyma.18g022500) and Glyma18g02610 (Glyma.18g022700), between these two different types have been elucidated (Cook et al. 2012(Cook et al. , 2014. In our study, QTL for resistance to RN on Chr. 18 in PI positions of the markers are given in centimorgan (cM). SRKN southern root-knot nematode, RN reniform nematode 567516C was located in a similar region as the Rhg1 locus for resistance to SCN (Table 2; Fig. 2) . Interestingly, a previous study on resistance to SCN did not detect the Rhg1 locus in PI 567516C (Vuong et al. 2010). In order to determine whether PI 567516C has the Rhg1 locus, we analyzed copy number and sequence variation of the Rhg1 locus in PI 567516C using whole-genome resequencing approach. The results showed that PI 567516C had a similar level of CNV at the Rhg1 locus with Peking and PI 437654, two of which had Peking-type Rhg1 (Fig. 3a). These three soybean lines had a reduced copy number of Rhg1 in comparison with PI 88788 (Fig. 3a). In addition to CNVs, we found that PI 567516C also had identical non-synonymous SNP/ Insertions in the exon of Glyma18g02590 (Glyma.18g022500) with Peking and PI 437654 (Fig. 3b). For Glyma18g02580 (Glyma.18g022400), PI 567516C, PI 88788, Peking and PI 437654 had one identical SNP in the exon in comparison with Williams 82 and Hutcheson (Fig. 3b). For Glyma18g02610 (Glyma.18g022700), no SNP was found among these five soybean lines. All these results indicated that there might be Peking-type Rhg1 in PI 567516C.

Discussion
In our study, the value of broad-sense heritability of the RI after infestation with RN was 0.32, indicating that the error variance rather than genetic variance accounted for a large part of all phenotypic variances. This might be the reason why total phenotypic variation explained by two QTL for resistance to RN in PI 567516C was small, only 16.4 % ( Table 2). Large error variance might be attributed to environmental factor, which indicated that environmentcontrolled bioassays are warranted to reduce error variance.
Previous studies showed that the major QTL on Chr. 10 played an important role in conferring resistance to SRKN in soybean (Tamulonis et al. 1997;Li et al. 2001;Pham et al. 2013;Xu et al. 2013). In our study, this QTL was also detected and consistently mapped at the same genomic region on Chr. 10, which indicated that different sources of resistance to RKN might harbor the same gene for resistance. Glyma10g02140 (Glyma.10g017100) and Glyma10g02160 (Glyma.10g017200) encode pectin methylesterase inhibitors (PMEIs) that played a role in Fig. 3 Genomic characterization of the Rhg1 locus among PI 88788, Peking, PI 437654, PI 567516C, cv. Magellan and cv. Hutcheson using whole-genome resequencing approach. a Copy number variation at the Rhg1 locus. Hut Hutcheson. b Sequence variation in the exon of Glyma18g02590 (Glyma.18g022500).
Glyma18g02590 (Glyma.18g022500) was predicted to encode an a-soluble NSF attachment protein (a-SNAP). The differences were highlighted with specific colors. Hyphen represents that the position is identical to Williams 82 sequence. (Color figure online) regulating the activity of pectin methylesterases (PME) in plant (Giovane et al. 2004;Pham et al. 2013). Both of these two genes were reported to be upregulated in SRKN-susceptible reaction but suppressed in SRKN-resistant reaction (Pham et al. 2013). PME is an enzyme involved in the process of cell wall breakdown. In Arabidopsis, loss of function of PME3 led to higher resistance to sugar beet cyst nematode, while overexpression of PME3 increased susceptibility (Hewezi et al. 2008). We performed sequence analyses for these two genes among SRKN-resistant and SRKN-susceptible soybean lines. The results indicated that Glyma10g02140 (Glyma.10g017100) rather than Glyma10g02160 (Glyma.10g017200) was more likely to be the candidate gene for QTL on Chr. 10. Besides QTL on Chr. 10, one minor QTL on Chr. 13 was identified in PI 567516C (Table 2). This QTL was also mapped in PI 438489B (Xu et al. 2013). Confirmation of this QTL in our study suggested that it was conservative and could be potentially used in soybean breeding. QTL on Chr. 17 for resistance to RN was first reported by our study. This QTL explained similar ratio of total phenotypic variation as the QTL on Chr. 13 (Table 2), which indicated that it was also an important locus for resistance to SRKN. Identification of this QTL made PI 567516C to be unique from the other two sources of SRKN resistance, PI 96354 and PI 438489B, where only two SRKN-resistant QTL were mapped (Li et al. 2001;Xu et al. 2013).
Three QTL have been reported to be associated with resistance to RN in PI 437654 (Ha et al. 2007). Among them, two QTL, one on Chr. 18 and the other on Chr. 11, were found to play critical roles in resistance to RN (Ha et al. 2007). These two QTL were also identified to be responsible for resistance to SCN in PI 437654 (Wu et al. 2009), which provided evidences that resistance to RN and SCN might concurrently share common mechanisms in soybean. QTL on Chr. 18 for resistance to RN in PI 567516C was mapped in the Rhg1 region (Fig. 2), indicating that PI 567516C might have Rhg1. Through wholegenome resequencing method, we provided evidences that there might be Peking-type Rhg1 in PI 567516C. Without Rhg4 (Liu et al. 2012), Peking-type Rhg1 was not able to function as a gene for resistance to SCN (Meksem et al. 2001). This might be the reason why Rhg1 was not detected in PI 567516C in the previous study on resistance to SCN (Vuong et al. 2010). We speculated that Peking-type Rhg1 might be in much more soybean germplasm we expected.
In a previous study, the identification of a novel QTL on Chr. 10 for resistance to SCN in PI 567516C enabled this soybean accession to be a good candidate to be used in breeding cultivars with new resistance to SCN that are different from Rhg1-and Rhg4-containing source of resistance to SCN (Vuong et al. 2010). In our study, the identification and characterization of QTL underlying resistance to RKN and RN made this PI 567516C more valuable in soybean breeding programs to develop new cultivars with broad-based resistance to multiple nematode species, not only SCN but also SRKN and RN, through marker-assisted selection.