QTL for Yield Traits and Their Association with Functional Genes in Response to Phosphorus Deficiency in Brassica napus

Background Oilseed rape (Brassica napus L.) is one of the most important oil crops. A primary limitation to the cultivation of this crop is the lack of available phosphorus (P) in soils. To elucidate the genetic control of P deficiency tolerance in Brassica napus, quantitative trait locus (QTL) for seed yield and yield related-traits in response to P deficiency were identified using a double haploid mapping population (TN DH) derived from a cross between a P-efficient cultivar, Ningyou 7 and a P-inefficient cultivar, Tapidor. Results Three field trials were conducted to determine seed yield (SY), plant height (PH), number of primary branches (BN), height to the first primary branch (FBH), relative first primary branch height (RBH), pod number per plant (PN), seed number per pod (SN) and seed weight of 1,000 seeds (SW) in 188 lines of TN DH population exposed to low P (LP) and optimal P (OP) conditions. P deficiency decreased PH, BN, SN, PN and SY, and increased FBH and RBH with no effect on SW. Three reproducible LP-specific QTL regions were identified on chromosomes A2, A3 and A5 that controlled SN, PN and SW respectively. In addition, six reproducible constitutive regions were also mapped with two each for SY-LP on A2, and FBH-LP on C6 and one each for PH-LP and SW-LP on A3. About 30 markers derived from 19 orthologous genes involved in Arabidopsis P homeostasis were mapped on 24 QTL regions by comparative mapping between Arabidopsis and Brassica napus. Among these genes, GPT1, MGD2 and SIZ1 were associated with two major loci regulating SY-LP and other yield-related traits on A2 between 77.1 and 95.0 cM. Conclusion The stable QTLs detected under LP conditions and their candidate genes may provide useful information for marker-assisted selection in breeding high-yield B. napus varieties with improved P efficiency.


Introduction
Phosphorus (P) is one of the essential macronutrients for plants and is a major limitation for plant development worldwide due to its low availability and inaccessibility in soil [1][2][3]. To mitigate this, P fertilizers are applied to improve P availability in soils. However, cost and low recovery rates of P fertilizers not only increase the cost of crop production but also lead to serious environmental pollution and exhaust non-renewable phosphate resources. Therefore, a sustainable strategy for crop production would be to breed crops with high efficiency in acquiring P from native soil reserves or fertilizer sources.
Nutrient-efficient plants produce higher yields per unit of nutrient applied or absorbed compared with other plants grown under similar agroecological conditions [4]. One way to determine the genetic basis of tolerance to nutrient deficiency is to map quantitative trait loci (QTL) based on molecular markers. A number of QTLs associated with P efficiency traits have been detected in rice [5][6][7], maize [8], wheat [9,10], common bean [11,12] and soybean [13,14]. Using a soybean F 9 recombinant inbred lines (RILs), 13 QTLs associated with root traits and 18 associated with P-efficiency were identified under low P and high P conditions, and three QTL clusters were found to be associated with traits for root and P-efficiency at low P levels [14]. In wheat, each QTL detected for P use efficiency (PUE) was found to be linked to QTLs regulating agronomic traits [10]. Pup1 is a major locus that confers P deficiency tolerance in rice [5,7,15,16]. Overexpression of a Pup1-specific protein kinase gene (PSTOL1) in phosphorus-starvation-intolerant varieties significantly enhanced grain yield in phosphorus-deficient soil [17]. These results provide strong evidence for the hypothesis that enhancing P efficiency would improve agronomic performance of crops.
Brassica napus is one of the most important oilseed crops worldwide and is extremely sensitive to P availability in soil. Since high seed yield is one of the desired targets in B. napus breeding, various studies have reported the importance of identifying Pefficient germplasm in B. napus [18]. QTLs for seed yield relatedtraits in response to P deficiency were first identified in B. napus using an F 10 RIL population in two crop seasons. In these experiments, a number of low-P-specific QTLs were detected across seasons. Four markers developed from genes involved in Arabidopsis P homeostasis were mapped to the confidence intervals of QTLs associated with a P efficiency coefficient (the ratio of seed yield at low P level to that at adequate P level) or low-P-specific QTLs [19]. Identification of QTLs associated with root traits and P-efficiency parameters partially revealed the genetic basis of tolerance to P deficiency and also confirmed the importance of root morphology traits in adaption to low P levels [20,21]. High colinearity and synteny between the Brassica species and Arabidopsis genomes have allowed in silico mapping of multiple Arabidopsis thaliana genes involved in P metabolic pathway as well as genes controlling agronomical traits in B. napus QTLs. These include QTLs for shoot mineral concentrations [22], seed mineral concentrations [23], root morphological traits [20,21], flowering time [24] and seed yield traits [19,25]. The sequence of B. rapa genome [26] and B. oleracea genomes (http://www.ocri-genomics. org/) has further facilitated identification of conserved genome segments and prediction of homologous genes for target traits in B. napus, in addition to confirming QTLs detected previously for shoot P content and P uptake efficiency related traits in Brassica species [21][22][23][27][28][29][30].
In the present study, a double haploid (DH) population derived from a cross between P-efficient cv. Ningyou 7 and P-inefficient cv. Tapidor (TN DH population) was used to determine seed yield and yield-related traits under low P (LP) and optimal P (OP) conditions in three trials conducted during four years. In addition, QTLs controlling seed yield and yield-related traits conferring Pdeficiency tolerance were also determined. Comparative mapping between the linkage groups of TN DH map (B. napus), and the B. rapa and B. oleracea genomes was also performed to associate QTLs with orthologous genes involved in tolerance to P deficiency. QTLs that confer P efficiency and their candidate genes would provide useful information to better understand P efficiency and to develop molecular markers to improve seed yield in B. napus cultivated in environments with low P.

Phenotypic analysis of the tested traits in TN DH population
Phenotypic differences in the tested traits between the two parental lines under low P (LP) and optimal P (OP) conditions are presented in Figure 1. Compared to the P-inefficient parent 'Tapidor', the P-efficient parent 'Ningyou 7' had higher number of primary branches (BN), seed weight of 1,000 seeds (SW), seed number per pod (SN) and seed yield (SY), and lower height to the first primary branch (FBH) and relative first primary branch height (the ratio of FBH to plant height (PH), RBH) under both P conditions. However, average pod number per plant (PN) in Ningyou 7 was lower than Tapidor under LP conditions. There was no significant difference in PH between the TN DH parents under two P conditions in Tri.1 (field trial conducted from Sept. 2008 to May 2009) and Tri.2 (field trial conducted from Sept. 2009 to May 2010), while PH in Ningyou 7 was higher than Tapidor in Tri.3 (field trial conducted from Sept. 2010 to May 2011). Each parent experienced lower PH, BN, SN, PN and SY in P deficient than optimal P conditions. In contrast, each parent had greater FBH and RBH in low than optimal conditions and experienced no different in SW across P levels.
TN DH population exposed to P deficiency also showed lower PH, BN, SN, PN and SY, and higher FBH and RBH with no effect on SW (Table 1 and Table S1; Figure 2 . Genotypic variation in each tested trait also existed in TN DH population under both P conditions. Moreover, variation coefficients for BN, SN, PN and SY in LP conditions were larger than in OP conditions in all three trials (Table 1). In addition to the strong effects of environment, genotype, envir-onment6P level and genotype6P level interactions for most traits, there were also strong effects of P level on all traits ( Table 2). Transgressive segregation and continuous distribution was observed for each trait in TN DH population exposed to both P conditions in all three trials ( Figure 2).
Pearson's correlation coefficient between traits was also calculated ( Table 3). The phenotypic correlation between SY and yield-related traits under LP was similar to OP conditions in the TN DH population. SY was positively correlated to PH, BN, SN and PN, and negatively correlated to RBH. However, higher correlation coefficients between SY and PH, RBH, SN, and lower correlation coefficients between SY and PN were observed under LP conditions than under OP conditions (Table 3).

QTL for seed yield and yield-related traits
A total of 155 significant QTLs were identified for the eight tested traits from all three trials, including 79 QTLs detected under LP conditions and 76 under OP conditions ( Figure 3 and Table S2). Intervals of several QTLs overlapped with each other resulting in clustering of the 155 QTLs into 81 regions on 15 chromosomes. A majority (80.0%) of these regions was distributed on A2, A3, A5, A9, C6 and C9 (Figure 3). Of the 81 QTL regions, 29 were specifically detected under LP condition, and 29 specifically under OP condition. The remaining 23 QTL regions were detected under both P conditions ( Figure 4 and Table S2).
The 29 OP-specific QTL regions were integrated from 40 QTLs detected under OP conditions, with individual R 2 ranging from 5.9-17.6% ( Figure 4 and Table S2). Among the 29 OPspecific regions, two regions were consistently mapped between 59.8-63.4 cM on A3 and 36.9-43.7 cM on A9 that regulate SW and PH respectively in two of the three field trials. In region 56.9-64.5 cM on C1, qSY-OP1-C1 coincided with qSN-OP2-C1 and the QTLs for FBH and RBH were co-located in five regions: 0.0-5.7 cM on A3, 18.2-28.1 cM and 128.2-141.2 cM on A9, 28.9-46.4 cM on C6 and 54.1-70.3 cM on C9. In addition, QTL for FBH also overlapped with QTL for PH in two regions: 28.4-34.8 cM on A3 and 28.7-31.0 cM on A9.
The 23 constitutive QTL regions detected under both P conditions were integrated from 73 QTLs, with individual R 2 varying from 5.7 to 20.2% ( Figure 4 and Table S2). Six regions were considered as major loci that respond to LP deficiency because these QTLs were consistently associated with the same trait under LP conditions in at least two trials. QTLs for SY-LP were detected in two of the three trials and mapped in the 77.1-91.3 cM and 91.3-95.0 cM regions on A2 accounting for 10.8-13.0% and 8.4-14.2% of SY-LP variation, respectively, and both also coincided with QTLs for SN-OP. The reproducible locus for PH-LP was located between 50.2-59.9 cM on A3, accounting for 7.9-9.2% of PH-LP variation. The reproducible locus for SW-LP was located in 70.2-77.0 cM region on A3, explaining 7.8-11.5% of SW-LP variation. The two major regions affecting FBH-LP were detected in all three field trials and mapped in the intervals of 46.0-50.8 cM and 51.6-61.5 cM on C6, explaining 9.2-17.8% and 12.4-18.5% of FBH-LP variation respectively.

Discussion
The purpose of this study was to detect P efficiency related QTLs regulating seed yield and yield-related traits using the TN DH population of B. napus derived from a cross between Pinefficient cv. Tapidor and P-efficient cv. Ningyou 7, and to associate QTLs with functional genes involved in P homeostasis.

Effect of P deficiency on seed yield and yield-related traits
Both parents, Tapidor and Ningyou 7, showed significant differences in the tested traits under both low P (LP) and optimal P (OP) conditions ( Figure 1). Ningyou 7 had higher number of primary branches (BN), seed weight of 1,000 seeds (SW), seed number per pod (SN) and seed yield (SY) and lower height to the first primary branch (FBH), relative first primary branch height (RBH) and pod number per plant (PN) than Tapidor under both P conditions. QTL mapping revealed that majority of the Ningyou 7 alleles among the identified QTLs were associated with increased BN, SW, SN and SY but decreased PN under both LP and OP conditions ( Table 4), indicating that the contribution of Ningyou 7 to seed yield was higher than Tapidor under both P conditions. P deficiency significantly decreased PN, BN, SN and plant height (PH), increased FBH and RBH, and had no effect on SW in the double haploid mapping population (TN DH) derived from a Table 1. Trait differences in a double haploid mapping population (TN DH) derived from a cross between a P-efficient cultivar, Ningyou 7 and a P-inefficient cultivar, Tapidor under low (LP) and optimal phosphorus (OP) conditions in three field trials. cross between Ningyou 7 and Tapidor (Tables 1 and 2) consistent with previous studies [19,25] indicating that SW was genetically more tightly controlled and less influenced by environmental factors than other yield-related traits. Of the seven yield-related traits, PN, BN and RBH were most severely affected by Pdeficiency (Table 1) in addition to having high correlation with SY (Table 3). Therefore, higher levels of PN and BN, and lower levels of RBH should be highly considered in the breeding of B. napus cultivars tolerant of environments with low P.
Dissection of stable QTLs associated with seed yield and yield-related traits conferring P-deficiency tolerance In this study, a total of 155 putative QTLs associated with seed yield and seven yield-related traits were identified under LP and OP conditions in all three field trials (Table 4). These QTLs were clustered and integrated into 29 LP-specific, 29 OP-specific and 23 constitutive QTL regions ( Figure 4 and Table S2). The genetic architecture of seed yield in B. napus was previously analyzed by QTL mapping of seed yield-related traits [25,[31][32][33][34][35]. About 85 QTLs associated with seed yield were identified along with 785 QTLs for eight yield-associated traits using TN DH population from 10 natural environments with optimal P [25]. In the present study, a majority of these OP-specific QTLs were identified in only one of the three trials and 10 OP-specific QTLs, i.e. all QTLs for PH-OP mapped on A9, were located in the QTL associated with traits detected in semi-winter environments by Shi et al. [25]. These results indicated that quantitative traits were influenced both by heredity and growth conditions.
In recent studies on B. napus LP-specific QTLs were reported for yield-related traits [19], root morphological traits [20] and shoot traits at the vegetative stage [21]. Further analysis of these QTLs suggested that the presence of specific genes or differential expression of some genes regulated the response to P deficiency. In this study, nine QTL regions were detected for the same traits under LP conditions in at least two field trials, including three LPspecific regions and six constitutive QTL regions ( Figure 4; Table  S2). These stable QTL regions associated with the response to P deficiency have been compared with QTLs for P efficiency identified in a BE RIL population by comparative mapping [19][20][21]23].  Table S2), the regions for SN-LP and PN-LP were co-located with the QTLs Table 2. Significance of three-way ANOVA analysis for the eight traits among a double haploid mapping population (TN DH) derived from a cross between a P-efficient cultivar, Ningyou 7 and a P-inefficient cultivar, Tapidor under low and optimal phosphorus conditions in three field trials. for SN and PN identified in the BE RIL population exposed to LP, respectively [19]. This region also coincided with QTLs for SN and PN in BE population exposed to LP [19]. The two major loci associated with FBH-LP that were repeatedly detected in all three trials was located on C6 between 46.0-50.8 and 51.6-61.5 cM, accounting for 9.2-17.8% and 12.4-18.5% of FBH-LP variation, respectively. The stable region for FBH-LP (51.6-61.5 cM on C6) associated with PH-LP and PH-OP was co-located with the reproducible QTL for PH in BE RIL population exposed to LP. The reproducible locus for SW-LP located between 70.2 and 77.0 cM on A3 also coincided with the stable loci for SW in BE RIL population exposed to LP [19]. These results indicated that common mechanisms may be involved in the response of different B. napus cultivars to P deficiency. In addition, stable QTL regions, 50.2-59.9 cM on A3 regulating PH-LP and 46.0-50.8 cM on C6 associated with FBH-LP, colocalized with QTL regions for root morphology traits conferring P-efficiency in BE RIL population [20,21]. This indicated that selection of early vigor plants with more developed root system at the seedling stage could enhance P acquisition, improve seed yieldrelated traits and therefore increase the seed yield at low P.  Table 3. Pearson's correlation coefficients among traits in a double haploid mapping population (TN DH) derived from a cross between a P-efficient cultivar, Ningyou 7 and a P-inefficient cultivar, Tapidor under low (below diagonal) and optimal phosphorus (above diagonal) conditions.    Detection of these reproducible loci in different populations that respond to P deficiency would provide useful information to conduct marker-assisted selection to breed B. napus with improved P-efficiency.

Gene prediction in Brassica by in silico comparative mapping
Unlike anonymous molecular markers which were widely used to construct B. napus genetic maps such as SSR and AFLP [24,36], GBMs developed from specific genes that are closely related to phenotypic traits could be useful in associating target traits with function [37,38]. For example, Yang et al. [20] found that the two mapped GBMs, BnIPS2 and BnGPT1, were linked to the QTL controlling root morphology in response to P starvation. The RI lines with favorable alleles related to BnIPS2 and BnGPT1 could develop larger root system, produce higher dry weight, acquire more P under LP conditions and produce more seeds [19,20]. In this study, a total of 30 GBMs derived from 19 genes involved in P homeostasis in Arabidopsis were found to be linked to 24 QTL regions associated with yield traits (Figure 4 and Table 4). BnWRKY75-A3 and BnPHT3-A3 designed from the phosphatestarvation-response genes AtWRKY75 and AtPHT3;1 respectively, were associated with the constitutive QTL region regulating FBH in LP and OP conditions, and RBH and SW in OP condition and located between 12.2-18.6 cM on A3 (Table 4). They are also associated with primary root length (PRL) in LP and OP conditions, and total root length (TRL) in OP condition (data not shown). The association between functional genes and QTLs may accelerate identification of candidate genes regulating target traits.
Comparative genome analysis is commonly performed between Brassica species and Arabidopsis genome to align conserved genomic blocks and associate genes with QTLs identified in Brassica species [24,[39][40][41]. For example, 40 synteny blocks and 84 islands were aligned between A. thaliana pseudochromosomes and the TN linkage groups and 153 orthologs of the flowering time genes in Arabidopsis were mapped in the QTL regions by in silico mapping [24]. QTL clusters on chromosome C3 associated with root morphology under P starvation of B. napus and QTLs for shoot dry weight, shoot P content and P use efficiency under LP conditions in B. oleracea were found to be located in the same block which also has a number of orthologous genes for root development, auxin transport and P metabolism in Arabidopsis [20,29]. However, many chromosomal regions from the present linkage maps could not be aligned with the Arabidopsis genome because of limited number of markers with known sequence information. The sequences of B. napus, B. rapa and B. oleracea are highly conserved. The annotated B. rapa [26] and available B. oleracea genome sequences (http:// www.ocri-genomics.org/) along with the single nucleotide polymorphism linkage maps of B. napus [42,43] will greatly facilitate the genomic alignment of Brassica species and allow candidate gene prediction in B. napus. In this study, a comparative genetic mapping was performed for the first time between the TN linkage map and the B. rapa and B. oleracea genomes. Most chromosomal regions in TN linkage groups can be aligned linearly with the B. rapa and B. oleracea genomes (Figure 4), which provided useful sequence information to develop molecular markers for QTL fine mapping in the future. Of the 19 P homeostasis genes in Arabidopsis located in the QTL intervals determined in this study, twelve were refracted to the corresponding segments of B. rapa or B. oleracea genomes (Figure 4 and Table 4). Among these genes, GPT1, MGD2 and SIZ1 were associated with two major loci for SY-LP, which co-located with yield-related traits between 77.1-95.0 cM on A2. GPT1 functioning as an importer of glucose-6-phosphate in Arabidopsis [44] was associated with PN detected under LP conditions by Ding et al. [19]. MGD2 encoded a type-B monogalactosyldiacylglycerol synthase involved in phospholipid metabolism and was induced under phosphate (Pi) deprivation [45]. AtSIZ1 encoded a plant small ubiquitin-like modifier (SUMO) E3 ligase and was a focal controller of Pi starvation dependent responses [46]. These three genes would be the most likely candidates for target QTLs that confer P deficiency tolerance and near-isogenic lines were developed to confirm the function of these genes. With deep sequencing already performed on the parents of TN DH population, sequence of orthologous genes in B. rapa and B. olereacea will provide useful information in the assembly and dissection of alleles between the parents, Tapidor and Ningyou 7.

Plant materials and Field trials
A B. napus DH population with 188 lines was used in this study. The DH population was derived from a cross between Tapidor and Ningyou 7 by microspore culture [47]. Ningyou 7 was characterized as a P-efficient cultivar with better growth and higher P acquisition than Tapidor under low P (LP) and optimal P (OP) conditions in pot culture [48].
Three field trials were conducted at an experimental agricultural research station in Qichun, Hubei Province in China (115u459N latitude, 30u199E longitude) which has sandy paddy soil. The first trial coded Tri. Seeds were sown in a nursery bed in the field in middle September and seedlings were transplanted 30 d after sowing. Plants were harvested in the following middle May. Standard agricultural practices were followed for field management.

Measurement of phenotypic traits
In each plot, six individuals from the middle row were used to determine height to the first primary branch (FBH) measured from ground level to the base of the lowest primary branch, plant height (PH) measured from ground level to the tip of the main inflorescence, number of primary branches (BN) measured as the number of primary branches arising from main shoot and seed number per pod (SN) measured as the average number of wellfilled seeds from 100 well-developed pods sampled from the primary branch in the middle of each plant studied. All representative individuals from each plot were harvested by hand at maturity stage to investigate seed yield per plant (SY) and seed weight of 1,000 seeds (SW). Pod number per plant (PN) was calculated using the following formula: PN = SY61000/ (SW6SN). Relative first primary branch height (RBH) was defined as the ratio of FBH to PH in each sampled plant. Tested traits were denoted with 'trait name and treatment name' such as PH-LP which indicates plant height under low P conditions. Linkage map construction and alignment with B. rapa and B. oleracea genomes About 53 GBMs designed from 46 genes regulating P homeostasis in Arabidopsis [49] were superimposed on the framework of the existing well-defined G-map [50]. In total, the new genetic linkage map had 798 molecular markers with the average distance of 2.57 cM between adjacent markers. Linkage map was constructed using JoinMap4.0 [51] and mapping was performed according to the method described by Long et al. [50]. Order of markers on the new linkage map was consistent with previously published maps [25,50].
The physical location of markers along with primer sequences in the TN DH linkage map present on B. rapa or B. oleracea genome was determined by in silico PCR. Primer sequences of markers were first aligned with B. rapa and B. oleacea genome sequences using BlastN with default parameters on the web site (http:// brassicadb.org/, http://www.ocri-genomics.org/). Physical position of the primer alignment on B. rapa or B. oleracea genome was determined by matching primer pairs that amplified a 200,2000 bp fragment. Physical position of markers was determined by aligning the product sequence of the TN DH linkage map against the B. rapa and B. oleracea genomes with evalue #1e 210 and at least 50% of the sequence length covered. Position of these markers on B. rapa and B. oleracea genomes were integrated and sorted according to their position on corresponding chromosomes and then were compared with the TN DH linkage groups for map alignment. For markers with multiple alignment positions, the best position was determined based on neighboring markers on the same linkage group.

Statistical analysis and QTL detection
Data analysis was conducted using SAS 8.1 (SAS Institute, Cary, NC, USA). Fisher's least significant difference (LSD) was used to test the significance of means at 0.05 levels. Correlation analysis was conducted to determine the relationship between the tested traits.

Supporting Information
Table S1 Seed yield and yield-related traits in 188 lines of a double haploid mapping population (TN DH) derived from a cross between a P-efficient cultivar, Ningyou 7 and a P-inefficient cultivar, Tapidor under low (LP) and optimal phosphorus (OP) conditions in three field trials.

(XLSX)
Table S2 Quantitative trait loci (QTLs) associated with seed yield and yield-related traits under low (LP) and optimal phosphorus (OP) conditions and distribution of gene-based markers (GBM) in the QTL intervals. (DOC)