QTL and Candidate Gene Identification for Silique Length Based on High-Dense Genetic Map in Brassica napus L.

Silique length (SL) is an important yield trait and positively correlates with seeds per silique and seed weight. In the present study, two double haploid (DH) populations, established from crosses Zhongshuang11 × R11 (ZR) and R1 × R2 (RR), containing 280 and 95 DH lines, respectively, were used to map quantitative trait loci (QTL) for SL. A high-dense genetic map from ZR population was constructed comprising 14,658 bins on 19 linkage groups, with map length of 2,198.85 cM and an average marker distance of 0.15 cM. Genetic linkage map from RR population was constructed by using 2,046 mapped markers anchored to 19 chromosomes with 2,217-cM map length and an average marker distance of 1.08 cM. Major QTL qSL_ZR_A09 and qSL_RR_A09b on A09 were identified from ZR and RR populations, respectively. Both QTL could be stably detected in four environments. QTL qSL_RR_A09b and qSL_ZR_A09 were located on 68.5–70.8 cM and 91.33–91.94 cM interval with R2 values of 14.99–39.07% and 15.00–20.36% in RR and ZR populations, respectively. Based on the physical positions of single nucleotide polymorphism (SNP) markers flanking qSL_ZR_A09 and gene annotation in Arabidopsis, 26 genes were identified with SNP/Indel variation between parents and two genes (BnaA09g41180D and BnaA09g41380D) were selected as the candidate genes. Expression analysis further revealed BnaA09g41180D, encoding homologs of Arabidopsis fasciclin-like arabinogalactan proteins (FLA3), as the most promising candidate gene for qSL_ZR_A09. The QTL identification and candidate gene analysis will provide new insight into the genomic regions controlling SL in Brassica napus as well as candidate genes underlying the QTL.


INTRODUCTION
The genetic control of domesticated traits has been studied in various crop plants (Nasyrov, 2003;Aitken et al., 2008), including rapeseed, by quantitative trait loci (QTL) mapping using different pairs of crosses between distinct varieties, for specific traits especially to enhance yield (Udall et al., 2006;Bouchet et al., 2014), oil content Rout et al., 2018), and disease resistance (Wu et al., 2019). Brassica napus L. (AACC,2n = 38), as one of the essential oil crops worldwide, is a polyploid crop derived by spontaneous hybridization between Brassica rapa (AA,2n = 20) and Brassica oleracea (CC,2n = 18) (Nagaharu, 1935) <7,500 years ago followed by the process of duplication (Boulos et al., 2014;Bayer et al., 2017). Comparative genomic studies showed that the A genome of B. rapa and B. napus exhibited collinearity with fewer genomic changes and translocations (Wang et al., 2011). Silique length (SL) is one of the most effective components for yield selection in rapeseed, by which not only seed yield can be increased but also the total oil yield (Samizadeh et al., 2007;Bennet et al., 2011). Information regarding inheritance pattern and gene action for SL is useful in rapeseed breeding program. However, the understanding of silique development in B. napus is limited and only few functional genes for SL have been explored due to high ploidy level.
In general, increase of SL may enhance the source of matter, which results in larger seeds (Chay and Thurling, 1989). To increase the yield of B. napus, identification of QTL and genes for SL is of importance for the development of cultivars with long siliques and incorporating the trait to improve yield (Youssefy et al., 2012). Number of siliques per plant varied independently on SL (Chay and Thurling, 1989). However, SL is positively correlated with seed weight and seed number per silique  and could serve as one trait for indirect selection of seeds per silique in the breeding process. As one of the important yield-related traits, SL has been extensively studied in terms of QTL analysis in rapeseed (Shi et al., 2019;Zhang et al., 2011;Li et al., 2014;Fu et al., 2015;Wang et al., 2016;Yang et al., 2017).
For QTL identification, second-generation molecular markers including simple sequence repeats (SSR) and amplified fragment length polymorphism (AFLP) are of low density across the genome, resulting in low-quality mapping (Yang et al., 2017). Advanced next-generation sequencing (resequencing) technologies provide new insight to rapidly identify QTL of interest. High-density genetic maps based on high-throughput SNP markers have the advantages of high throughput, low cost, high resolution, and high positioning accuracy, which can significantly improve the efficiency and precision of genotyping and genetic map construction. QTL for complex agronomic traits such as plant height, flowering time, oil content, and seed weight have been mapped based on high-throughput SNP markers (Yang et al., 2012b;Li et al., 2014;Wang et al., 2015).
Many QTL for SL or seed size have been identified in rapeseed (Zhang et al., 2011;Li et al., 2014;Fu et al., 2015;Liu et al., 2015;Wang et al., 2016;Yang et al., 2017). Eight QTL for SL distributed on A01, A06, A07, A09, and C06 were detected from immortalized F2 (IF2) and double-haploid (DH) populations in different environments (Wang et al., 2018). One major QTL was identified for SL and thousand seed weight on A09 that shared a common linked marker suggesting tightly linked or pleiotropic QTL (Qi et al., 2014). Though many QTL have been detected, only two genes affecting SL and seed size have been cloned and characterized. An auxin response factor gene (ARF18) involved in SL and seed weight has been cloned by fine mapping and association analysis . One 55-aa deletion prevents ARF18 forming homodimer and inhibited the activity of downstream auxin genes . Another research group identified one major QTL on A09 and then fine-mapped BnaA9.CYP78A9 as underlying gene for qSLWA9 locus (Yang et al., 2012b;Shi et al., 2019). A 3.7-kb insertion of a CACTA-like transposable element (TE) in the regulatory region functioned as an enhancer to increase the expression of BnaA9.CYP78A9 in varieties with long siliques and large seeds (Shi et al., 2019), whereas the insertion was replaced by a 12.3-kb deletion in varieties with short siliques. Interestingly, the physical distance of these two genes is less than 20 kb, suggesting the complex of regulating genes for SL in the A09 chromosome region.
In the present study, we identified QTL for SL by using two independently generated double haploid populations derived from crosses Zhongshuang11 × R11 and R1 × R2. Major QTL on A09 chromosome were stably detected in all environments from both populations, and the contribution rate is between 15.00% and 20.36%. Markers with physical position revealed that the qSL_ZR_A09 of ZR population may be one novel QTL on A09. Gene structure variation and expression analysis showed that BnaA09g41180D, which is homologous to Arabidopsis FLA3, may be the best promising genes for SL control. This QTL at the specific locus qSL.A09 and the information of candidate gene will be valuable for marker-assisted selection (MAS) breeding for higher yield and further identification of functional gene related to yield improvement in rapeseed.

Plant Materials
Four semi-winter-type inbred lines developed by OCRI-CAAS including Zhongshuang11 (ZS11), R11, R1, and R2 were selected for dissecting the genetic bases of SL in B. napus. ZS11 and R1 have long siliques while R11 and R2 have regular SL. Two DH populations were developed from the crosses between ZS11 × R11 and R1 × R2. The DH populations derived from the crosses of ZS11 × R11 and R1 × R2 parents were named as "ZR population" and "RR population", respectively. In total, 280 and 95 DH lines were generated from ZR and RR populations, respectively.

experimental Design and Trait Measurement
The DH populations including DH lines and parents were grown in field using a randomized complete block design with three replications during September to May in Wuhan (114.33°E, 30.50°N), Huanggang (114.87°E, 30.44°N), and Zunyi (106.55°E, 27.42°N) each year (Supplementary Figure 1). Location-year combination was coded as environments. The four environments for ZR population were WH15 (Wuhan, 2014-15), WH16 (Wuhan, 2015-16), WH17 (2016, and HG17 (Huanggang, 2016-17). And the four environments for RR DH population were WH13 (Wuhan, 2012-13), WH14 (Wuhan, 2013-14), ZY13 (Zunyi, 2012-13), and ZY14 (Zunyi, 2013-14). Each plot contained three rows and 20 plants in each row. The distance between rows was 33 cm and between plants was 10 cm. At maturation stage when most siliques turned yellow, 10 plants were randomly sampled from the middle row of each plot and considered as representative for trait measurement. Thirty welldeveloped siliques of the main inflorescence from each plant were measured for SL.

Statistical analysis
Analysis of variance (ANOVA) was performed using the procedure MEANS of SAS 9.1 (SAS, 2004) for descriptive phenotypic statistics. The phenotypic t test of parents was calculated using the procedure TTEST of SAS 9.1. The broad-sense heritability was calculated as h n nr g g ge e 2 2 2 2 2 = + + σ σ σ σ /( / / ) , where σ g 2 is genetic variance, σ ge 2 is the interaction variance of genotype with environment, σ e 2 is error variance, n is the number of environments, and r is the number of replications. Components of variance ( σ g 2 , σ ge 2 , and σ e 2 ) were estimated using the GLM procedure of SAS 9.1.

SNP Genotyping and Linkage Map Construction of RR Population
The Brassica 60 K SNP BeadChip Array was used to genotype 95 RR DH lines and two parental lines. Genetic linkage map of RR population has been constructed in our previous study (Liu et al., 2016).

Sequencing Library Construction and High-Throughput Sequencing
The genomic DNA of 280 ZR DH lines and parental lines (ZS11 and R11) was extracted from the leaves by CTAB method. The genomic DNA was sheared into ~500-bp fragments using the S2/ E210 Ultrasonicator (Covaris, USA). About 2 μg genomic DNA from each line was prepared for sequencing library construction. Sequencing library was constructed by "NEB DNA Library Prep Kit" according to the manufacturer's instructions. The sequencing libraries were constructed by terminal repair, then added 3′A and a sequencing linker. Sequence data were generated by Illumina HiSeq 4000 (San Diego, California, USA) with paired-end (PE150). Sequence quality and adaptor trimming was conducted by SOAPnuke 1.4 (BGI, Shenzhen, China).

SNP Identification and Genotyping
The raw reads which came from the original data of the Illumina HiSeq™ sequencing platform were processed using a series of in-house C-scripts for quality control (QC). During QC procedures, low-quality reads were filtered out and clean reads were obtained. Burrows-Wheeler Aligner (BWA) (Li and Durbin, 2009) software was used to align the clean reads of each sample against the B. napus genome (https://www.ncbi.nlm.nih. gov/assembly/GCF_000686985.1/). SAMtools (McKenna et al., 2010) software was used to convert the alignment files to bam files and remove duplicate reads. The bam files were filtered by GATK toolkit and potential SNPs between all lines and the genome were detected. The SNPs identified between the parents were considered as polymorphic for a subsequent bin calling. To reduce false-positive SNPs caused by sequencing errors, the base support for each SNP in the parents and DH lines were required to be 4 or more (and no more than 1,000). All SNP markers between the parents were classified into four segregation patterns (aa × bb, hk × hk, lm × ll, and nn × np). For the DH populations, the segregation patterns aa × bb were selected for construction of genetic map.

Bin Map and Linkage Map Construction
The recombinant bin map based on SNP markers was constructed by a modified slide window method which was developed by Huang's lab (Huang et al., 2009). The genotype of each window was described as a window of 15 SNPs and a step size of 1 SNP. Windows containing more than 13 "aa" or "bb" types were genotyped as aa or bb, respectively. Fifteen adjacent SNP intervals with the same genotype across the entire ZR population were combined into a recombination bin. The linkage map was constructed from the recombination bins serving as genetic markers using the HighMap program . The Kosambi mapping function was used to calculate the genetic distance between the markers (Kosambi, 2016).

QTL Mapping and SNP Bin Map for ZR Population
Data which were collected in different years were considered as an independent environment. WinQTLCart2.5 (Wang, 2006) was used to conduct QTL mapping of ZR and RR populations with composite interval mapping method and Model 6 Zmapqtl procedure. A 1,000 permutation test of shuffling the phenotypic means with genotypic means was performed to calculate the significance threshold of the test statistics for one QTL (Doerge and Rebaï, 1996). The parameters of control markers, windows size, and walk speed were fixed to 5, 10, and 1 cM, respectively. The marker distance was set to 5 cM to determine a significant QTL. The percentage of phenotypic variance explained by QTL was predicted by the highest probability peaks. For SL trait, QTL repeatedly detected at the same location of a chromosome with the identical direction of genetic effects were considered as the same QTL.

Candidate Genes analysis
To identify candidate genes for SL, genetic and physical map integration of QTL confidence intervals was performed. We used the 1-LOD confidence intervals as the QTL confidence intervals, and the flanking sequences of the SNP markers in QTL confidence intervals were mapped to the B. napus reference genome (http://www.genoscope.cns.fr/ brassicanapus) using BLAST tool to project QTL on the physical map. Candidate genes were identified based on the annotation of the Arabidopsis genome and physical positions of SNP markers in the reference B. napus genome.

expression analysis
Total RNA was extracted with Trizol Reagent (Invitrogen, America). Reverse transcription was performed according to the instruction of FastQuant RT Kit (Tiangen, China). RT-PCR was performed as described previously using the primers listed in

Phenotypic Variation of SL in Parents and DH Lines
The 280 ZR DH, 95 RR DH lines, and their corresponding parents were grown in winter-spring seasons. The field performance of DH populations for SL across the years was analyzed. The parents bearing long siliques showed maximum SL of 8.33 ± 0.37 cm (ZR crosses) in WH16, and minimum SL of 6.65 ± 0.21 cm (RR crosses) was exhibited in WH14 ( Table 1 and Figure 1). Oneway ANOVA across environments indicated that genotypes (representing established DH populations and their parents), growing environments, and genotype-environment (G × E) interaction had a highly significant effect on SL (P < 0.0001) ( Table 2). The broad-sense heritability (h 2 ) for SL was 94.53% in RR and 97.30% in ZR populations. The normal distribution of SL variation in both DH populations suggested that SL trait is controlled by multiple genes (Table 1).

High-Density SNP Genetic Map Construction
Re-Sequencing and Genotyping The RR population map contained 2,046 marker loci with an average marker density of 1.08 cM along a 2,217.2-cM genetic map. This linkage map was described in more detail elsewhere (Liu et al., 2016). For ZR population, over 609 Gb of clean data (~2 billion reads) was generated for the parental lines and 280 DH lines. In parental lines of ZS11 and R11, 50.78 Gb clean data (~170 million reads) and 62.51 Gb clean data (~209 million reads) were produced, respectively. While 495.75 Gb clean data (~16.53 billion reads) was generated for the 280 DH lines with

QTL Mapping
In total, nine repeatedly detected QTL were identified for SL in RR and ZR populations ( Table 4). Five QTL were frequently detected in different environments in RR population and were located on A01, A07, and A09. The explained phenotypic variance ranged from 3.28% to 39.07%, and the LOD values ranged from 2.53 to 21.99. A major QTL, qSL_RR_A09b (R 2 > 14%) was detected in all environments (Table 4 and Figure 4A). Four QTL were repeatedly identified in different environments in ZR population and were located on A08, A09, C03, and C06. The explained phenotypic variance ranged from 3.01% to 20.36%, and the LOD values ranged from 2.57 to 17.50. A major QTL, qSL_ZR_A09, explaining more than 15% of the total SL variation was repeatly detected in all environments (Table 4 and Figure  4B). Another QTL, qSL_ZR_C06, was also defined in all four environments and explained about 7% phenotypic variance in all environments. Different from RR population, two QTL for SL, including qSL_ZR_A08 and qSL_ZR_C03, were repeatedly detected in two or three environments in ZR population.

Candidate Gene Identification
According to the results of QTL mapping, two major QTL, qSL_RR_A09b from RR population and qSL_ZR_A09 from ZR population, were detected in all four environments. We then mapped the flanking sequences of the SNP markers in QTL confidence intervals to the B. napus reference genome (http:// www.genoscope.cns.fr/brassicanapus/) using BLAST tool. The physical location of these two QTL was found on A09 between 27,577,590-28,622,867 bp and 28,619,958-28,994,184 bp, respectively ( Table 5). This physical location overlaps with the SL QTL reported before (Dong et al., 2018). There are two genes,   BnaA9.ARF18 and BnaA9.CYP78A9, controlling SL cloned in B. napus; both are located on A09 Shi et al., 2019). The CACTA-like TE which inserted into the upstream region of BnaA9.CYP78A9 acting as an enhancer to increase SL Shi et al., 2019) was detected in R1 (long silique), but not in R2 (short silique) (Figure 5A), suggesting BnaA9.CYP78A9 might be the responsible gene for qSL_RR_A09b in RR population. In contrast, qSL_ZR_A09 in ZR population located in the physical interval of 28,619,958-28,994,184 was diverged from the QTL reported before, suggesting it is a novel QTL for SL control.
To determine whether BnaA9.ARF18 and BnaA9.CYP78A9 were the SL genes for qSL_ZR_A09 we further compared sequence variation between the two parental lines. No SNP was detected for BnaA9.ARF18 sequence between two parents ZS11 and R11, and both ZS11 and R11 contained the 3.7-kb CACTA-like TE insertion (Figure 5A), excluding BnaA9.ARF18 and BnaA9. CYP78A9 as the responsible genes underlying qSL_ZR_A09. We then conducted candidate gene identification for this locus.
From the physical interval of qSL_ZR_A09, 78 candidate genes were found against the B. napus reference genome and 26 genes harbor SNPs or indels in gene or promoter region (Supplemental File 1, Table S3). According to the annotation from Arabidopsis, two candidate genes (BnaA09g41180D and BnaA09g41380D), respectively encoding homologs of Arabidopsis FLA3 and AGL61, which have been demonstrated to be involved in silique developmental regulation were identified. RT-PCR was performed to investigate the expression pattern of these two genes from five different tissue samples. The mRNA transcript of BnaA09g41180D was detected in flower buds, stem, and silique, but not in leaf sample. The BnaA09g41380D transcript was detected in all five  tissue samples with various intensities (Figure 5B). The mRNA transcript of BnaA09g41180D was enriched more in ZS11 than R11 in young silique sample ( Figure 5B). Fifty-one SNPs or indels were identified in the promoter region of BnaA09g41180D between ZS11 and R11. Among them, 42 including 31 SNPs and 11 Indels were located in the regulating cis-element of the promoter region. Thus, we speculated that BnaA09g41180D is the most promising candidate gene controlling SL at qSL_ZR_A09.

DISCUSSION
The B. napus genome was sequenced (Chalhoub et al., 2014) and various genetic maps had been published (Zhang et al., 2014), which enabled genotyping and provided a more effective method for exploring target regions with plenty of SNPs (Yang et al., 2012a;He et al., 2014). Traditional QTL prediction was mostly conducted by genotyping with a large number of individuals in a primary population, which is a time-consuming and painstaking work (Salvi and Tuberosa, 2005). In rapeseed, different techniques have been employed for QTL mapping (Uzunova et al., 1995;Zhao et al., 2006;Fu et al., 2007;Long et al., 2007;Radoev et al., 2008;Zhao et al., 2008;Shi et al., 2009;Yang et al., 2012b;Zhang et al., 2012;Raman et al., 2013;Wang et al., 2013;Raman et al., 2016;Yang et al., 2018). Identification of QTL by high-density SNP array analysis (Clarke et al., 2016) and re-sequencing (Fletcher et al., 2016) was proven a more efficient and quick method to ascertain the genomic regions of QTL. In this study, the ZR DH lines were re-sequenced and exhibited considerable variation in SL traits (Figure 1) and were ideal for genetic map construction (Figure 2) and QTL detection (Figure 4). The high-density genetic map developed in this study covered a total 2,198.85 cM with 14,658 recombination bin markers, indicating high-density markers for possible application in gene identification studies. Four QTL for SL were identified located on A08, A09, C03, and C06 chromosomes in ZR DH population, explaining 3.01-20.36% phenotypic variation. Five QTLs located on A01 (2), A07, and A09 (2) were also detected in RR population, explaining 3.28 to 39.07% phenotypic variance. Major QTL was referred as QTL occurring at least once with R 2 ≥ 20% or at least twice with R 2 ≥ 10% (MacCaferri et al., 2008;Shi et al., 2009). Four QTL, qSL_RR_A09b occurring at four environments with R 2 from 14.99% to 39.07%, qSL_ZR_A08 occurring at two environments with R 2 from 14.82% to 15.89%, qSL_ZR_A09 occurring at four environments with R 2 from 15.00% to 20.36%, and qSL_ZR_C06 occurring at two environments with R 2 from 10.99% to 11.67%, are all major QTL for SL.
Many QTL detected by this study indicated the polygenic nature of SL trait. Thus, different environment conditions such as rainfall, temperature, and soil may have some influence on the SL phenotype. As a matter of fact, our result on qSL_RR_A09b which explained 14.99-39.7% phenotypic variance also approved SL as a polygenic trait. This renders the difficulty in breeding for long siliques and the importance of marker assistant selection. Thus, the development of allele-specific markers for different loci is of great value for improving SL in breeding. For example, the favorable alleles of qSL_ZR_A09 and qSL_ZR_C06, the two major QTL on A09 and C06 detected in ZR population, are from different parents. These favorable alleles can be pyramided using markers specific to the alleles and thus, longer SL can be readily acquired after pyramiding the different alleles from different loci.
We identified several QTL for SL including qSL_ZR_A09 and qSL_RR_A09b in two genetic populations, revealing the possible existence of stable genetic loci on the A09 chromosome.
In previous SL QTL mapping studies, different QTL were identified on different LGs (Uzunova et al., 1995;Chen et al., 2007;Radoev et al., 2008;Xu et al., 2010;Yang et al., 2012b;Fu et al., 2015;Yang et al., 2017;Dong et al., 2018). One major effect QTL on A09 explained more than 50% variation as compared to a minor QTL revealing <10% variation (Yang et al., 2012b). Two genes, BnaA.ARF18.a and BnaA9.CYP78A9, have been demonstrated to be involved in the regulation of both SL and seed weight in rapeseed and are located on A09 Shi et al., 2019). The QTL region of qSL_RR_A09b from RR population overlapped with BnaA.ARF18.a and BnaA9. CYP78A9. The CACTA-like transposable element inserted into the upstream region of BnaA9.CYP78A9, which contributes to the increase of SL, was found in R1 with long silique but not in R2 with short silique, indicating that QTL qSL_RR_A09b is identical to the QTL for seed weight (SW) and SL identified previously (Yang et al., 2012b;Li et al., 2014;Fu et al., 2015;Liu et al., 2015;Shi et al., 2019), and BnaA9.CYP78A9 is the responsible gene for this locus.
The qSL_ZR_A09 detected from ZR population showed a moderate effect (with R 2 value from 15% to 20.3%) in all environments, which could also be considered as major effect QTL. The physical region of qSL_ZR_A09 was not overlapped with the QTL interval on A09 reported before. Both the parents (Z11 and R11) of ZR population contain the CACTA-like transposable element insertion ( Figure 5A) and had no sequence variation in BnaA.ARF18.a. Thus, qSL_ZR_A09 is very likely one novel QTL controlling SL. In the novel QTL region, two candidate genes (BnaA09g41180D and BnaA09g41380D) which are homologs of Arabidopsis FLA3 and AGL61, respectively. In Arabidopsis, modulation expression of FLA3 influenced elongation of the stamen filament and female fertility, which led to the change of silique length . Siliques of the heterozygous mutant of agl61 are smaller than the wild-type siliques and seed set was also reduced (Bemer et al., 2008). The mRNA transcript of BnaA09g41180D was increased in ZS11 compared to R11 in young silique sample. Many SNPs or indels were identified in the promoter region of BnaA09g41180D between ZS11 and R11, with most of the sequence variation in the promoter region. The variation regulating cis-element in the promoter region may lead to differential expression level of candidate genes. Thus, we considered that BnaA09g41180D is the most promising candidate gene for qSL_ZR_A09.
In summary, we have developed two DH populations and conducted high-density SNP maps to identify novel and major QTL for SL, which is one of the most important yield components. Two major QTL on A09 were identified and candidate genes were predicted. One of the major QTL identified is probably a previously identified SL gene and the other one is a newly found SL gene in B. napus. Sequence variation of these genes ensures the development of gene-specific molecular markers and thus will facilitate pyramiding of favorable alleles of these genes via MAS in the breeding program for yield increase. As both identified QTL are located on A09, this result may be helpful for breeding FIGURe 5 | Detection of specific insertion of CACTA-like transposable element (TE) and expression comparison of two candidate genes. (a) The PAV PCR marker was used to detect the presence or absence of the 12.3kb fragment and the presence or absence of the CACTA-like TE. TEa-F and TE-R specifically detect the presence of the 12.3-kb fragment and the absence of the TE, and the combination of TEp-F and TE-R specifically detects the presence of the TE and the absence of the 12.3-kb fragment. (B) Expression analysis between ZS11 and R11 in five tissue samples. 1 leaf, 2 stem, 3 flower bud with 2-2.5 mm, 4 flower bud with 3-3.5 mm, 5 young silique.
Frontiers in Plant Science | www.frontiersin.org November 2019 | Volume 10 | Article 1579 of other Brassica species including B. rapa and Brassica juncea. Specific markers developed for the main effect QTL in this study will be useful for marker-assisted selection in breeding. Main effect QTL detected in the population will also provide target for map-based cloning of genes controlling silique length. Sequence variation of candidate genes will ensure the development of genespecific markers and thus facilitate pyramiding breeding.

DaTa aVaILaBILITY STaTeMeNT
The datasets generated for this study can be found in NCBI https://www.ncbi.nlm.nih.gov/bioproject/PRJNA588593/.

aUTHOR CONTRIBUTIONS
HW and QH designed the experiment. HW, WH, LF, and DM contributed to phenotypic measurements. WW, HC, and JL contributed to data analysis. BD and MH contributed to the expression and gene analysis. HW and QZ wrote the manuscript. HC and QH revised the manuscript. All authors reviewed and approved this submission.