A Molecular Genetic Linkage Map of Eucommia ulmoides and Quantitative Trait Loci (QTL) Analysis for Growth Traits

Eucommia ulmoides is an economically important tree species for both herbal medicine and organic chemical industry. Effort to breed varieties with improved yield and quality is limited by the lack of knowledge on the genetic basis of the traits. A genetic linkage map of E. ulmoides was constructed from a full-sib family using sequence-related amplified polymorphism, amplified fragment length polymorphism, inter-simple sequence repeat and simple sequence repeat markers. In total, 706 markers were mapped in 25 linkage groups covering 2133 cM. The genetic linkage map covered approximately 89% of the estimated E. ulmoides genome with an average of 3.1 cM between adjacent markers. The present genetic linkage map was used to identify quantitative trait loci (QTL) affecting growth-related traits. Eighteen QTLs were found to explain 12.4%–33.3% of the phenotypic variance. This genetic linkage map provides a tool for marker-assisted selection and for studies of genome in E. ulmoides.


Introduction
Eucommia ulmoides Oliver (2n = 34), the single extant species of the genus Eucommia (Eucommiaceae), is strictly a dioecious perennial tree [1]. It is an economically important plant for both herbal medicine and organic chemical industry. Chemical constituents (e.g., phenylpropanoids and flavonoids) in the bark and leaves have high pharmacological activities and health care functions of lowering blood pressure and blood sugar, resisting oxidation and mutation, improving the health, strengthening the body, promoting metabolism and relieving tiredness [2][3][4][5]. The whole plant except xylem contains Eucommia-rubber which is an important raw material in the chemical industry. Eucommia-rubber is a hard rubber with thermoplasticity, and it has properties that are similar to those of plastic [6]. Historically, only the bark was officially recognized as a traditional Chinese herbal drug. In recent years, the bark of E. ulmoides also was used to produce Eucommia-rubber in China, Russia and Japan. To improve the quality and yield of the bark, height and diameter growth were the main parameters for selection [7].
Conventional breeding of E. ulmoides has mainly focused on the selection of promising plants from existing natural populations. These selected plants were propagated vegetatively and released as clones. Recently, these cultivars were used as parents in crossbreeding. However, classical breeding often takes decades to fully evaluate and release new cultivars. The ability of E. ulmoides breeders to select promising parents for crossing, and to identify progenies with favorable combinations of characters, is hampered by the limited knowledge of the genetic basis of economically important traits. The speed and precision of breeding can be improved by the development of genetic linkage maps. Such genetic linkage maps can facilitate the development of diagnostic markers for polygenic traits and the identification of genes controlling complex phenotypes. The linked molecular markers identified in quantitative trait loci (QTL) analysis could then potentially be used in breeding practice via marker-assisted selection, where the selection is based on DNA sequence rather than the phenotype.
For forest trees, given the high genetic load and long generation time, segregating populations derived from crosses between inbred lines are not available. To circumvent this limitation, a pseudo-testcross approach is generally used to construct linkage maps from full-sib populations. Combined with the pseudo-testcross strategy, molecular markers such as random amplified polymorphic DNA (RAPD), sequence related amplified polymorphism (SRAP), amplified fragment length polymorphism (AFLP), inter-simple sequence repeat (ISSR) and simple sequence repeat (SSR) have been used extensively for the preparation of linkage maps of a number of tree species [8][9][10][11][12]. In a pseudo-testcross, only dominant markers that segregate in a 1:1 ratio are used to build separate molecular maps for each parent. Considering modern marker technologies are available for full-sib populations, markers that segregate in 3:1 (dominant), 1:2:1 (co-dominant) and 1:1:1:1 (co-dominant) ratios, in addition to 1:1, can be used to integrate individual linkage maps [13][14][15]. Using genetic linkage maps, QTL analysis have been conducted for traits of leaf, growth, vegetative propagation, wood quality, resistance, yield, flowering and fruiting in tree species [16][17][18][19][20][21].
In order to construct a genetic linkage map of E. ulmoides, we produced a F 1 mapping population from the cross between a wild genotype Xiaoye and a cultivar Qinzhong No.1. The female parent Xiaoye originated from the forest at Yantuo, Lingbao, Henan. The male parent Qinzhong No.1 was one of the four earliest cultivars [22], and it was planted in the museum garden of Northwest A&F University, Yangling, Shaanxi. Xiaoye and Qinzhong No.1 were chosen as parents because they differ in important quantitative traits. For instance, Xiaoye has late budding and flowering times, low content of secondary metabolite, small leaves, and smooth bark, whereas Qinzhong No.1 has early budding and flowering times, high content of secondary metabolite, large leaves, and rough bark. Besides, Qinzhong No.1 is an excellent cultivar, fast growing and with high resistance to drought and cold. In this study, we present a genetic linkage map of E. ulmoides based on SRAP, AFLP, ISSR and SSR markers. Results from our QTL analysis for height and basal diameter, measured over four consecutive years, are reported.

Genetic Linkage Map
The genetic linkage map (DZ0901) consisted of 706 markers distributed over 25 linkage groups (LG) covering 2133 cM (Table 5 and Figure 1). The number of mapped makers per LG ranged from 5-106 with a mean of 28.2. The map size of the LGs ranged from 19.9-194.0 cM with a mean of 85.3 cM. The average map distance between adjacent markers was 3.1 cM. In addition, 165 markers distributed over 25 triplets and 45 doublets. There were 628 unlinked markers and 20 markers that successfully linked with a group but could not be ordered. Since our estimate of E. ulmoides genome length was 2403 cM, the genetic linkage map constructed in our study covered approximately 89% of the genome. Table 1. Primer sequences used in the sequence-related amplified polymorphism analysis.

Growth Traits and QTL Analysis
A high degree of genetic variation was found for height and basal diameter ( Table 6). Figure 2 showed the frequency distributions of these traits. Pearson correlation analyses showed significant correlations between height and basal diameter and moderate weak correlations over years (Table 7).     The genetic linkage map of DZ0901 was used to search for putative QTLs (Table 8 and Figure 1). Eleven height QTLs were detected. In 2010, one height QTL was located on LG18 and explained 17.1% of the phenotypic variation. In 2011, three additional height QTLs were located on LG10, LG10 and LG12, and explained 29.7%, 27.7% and 22.8% of the phenotypic variation, respectively. In 2012, three height QTLs were located on LG9, LG13 and LG22, and explained 12.6%, 12.4% and 33.3% of the phenotypic variation, respectively. In 2013, two height QTLs identified at similar genomic regions as the height QTLs in 2012 were located on LG9 and LG22, and explained 13.5% and 25.3% of the phenotypic variation, respectively. Other two height QTLs in 2013 were located on LG21 and LG24, and explained 26.6% and 27.1% of the phenotypic variation, respectively. Four basal diameter QTLs were identified at similar genomic regions as the height QTLs. In 2010, one basal diameter QTLs was located on LG18 and explained 13.4% of the phenotypic variation. In 2011, one basal diameter QTLs was located on LG12 and explained 20.2% of the phenotypic variation. In 2012, two basal diameter QTLs were located on LG21 and LG22, and explained 25.1% and 21.4% of the phenotypic variation, respectively. Three additional basal diameter QTLs were detected on LG18, LG1 and LG1, and explained 29.8%, 17.7% and 16.8% of the phenotypic variation, respectively. Four of the 18 QTLs were significant, and they were Dht0-1, Dht1-1, Dht1-2 and Dbd0-2. Other QTLs were not significant, but they had a LOD score greater than 3.0. Flanking markers and QTLs supported by Kruskal-Wallis nonparametric test were indicated in Table 8.

Marker Amplification
SRAP has been recognized as an efficient and useful marker system [12,13,23]. It has several advantages such as simplicity, high throughput, numerous co-dominant markers and easy isolation of DNA fragments for sequencing, and it targets open reading frame regions. In this mapping population, SRAP analysis was an efficient method for generating polymorphic markers. Every primer combination gave at least six polymorphic markers with an average of 12.2 per primer combination. This is comparable to the polymorphism in other tree mapping projects using SRAP analysis [12,13].
It is known that AFLP marker produces a larger number of polymorphic fragments than other techniques. In our study, the average number of polymorphic DNA fragments per primer combination was 16.4. This is comparable to the average obtained in other tree mapping projects [9,17], but lower than the average reported for mapping using interspecific crosses [24,25]. As reported on many plant species [24,25], AFLP markers were dominant in this mapping population. ISSR analysis has been used successfully to construct genetic linkage map of many tree species [9,10,12]. In these studies, ISSR markers were highly polymorphic and tended to be evenly distributed throughout genomes. Besides, the ISSR analysis was faster and easier than the AFLP analysis. However, it was less efficient with an average of 6.9 polymorphic markers per primer and had a limited number of primers. Like AFLP marker, ISSR markers were dominant in this mapping population.
SSR markers are typically co-dominant, highly polymorphic and highly reproducible across laboratories. They are also useful for comparing and combining linkage maps from different mapping populations. Furthermore, many SSR markers are transferable across related species [26,27]. Unfortunately, E. ulmoides is the single extant species of the genus Eucommia, and there are fewer available SSR primer combinations. In this study, we used 19 SSR primer combinations. Additional SSR markers are currently being added to better bridge this map with future E. ulmoides maps.

Segregation Distortion
In this study, 29% of the markers showed segregation distortion. We excluded these markers to obtain a more accurate genetic linkage map because distorted markers can affect the mapping accuracy by overestimating the map distances and causing marker clustering [9,28]. Also, the order of markers on linkage groups may be affected by segregation distortion [29]. We may have lost some information by excluding the distorted markers. However, we obtained a genetic linkage map covering approximately 89% of the estimated E. ulmoides genome with an average of 3.1 cM between adjacent markers. In a follow-up study, we intend to map these distorted markers using a larger mapping population and co-dominant markers.
Segregation distortion has been reported frequently in woody species. The percentage of markers showing segregation distortion was highly variable: 47% in spruce [28], 38% in pear [25], 29% in citrus [12], 18% in Salix [8], 9% in grape [9], 8.5% in Populus [13] and 1.8% in peach [30]. Compared to these data, distorted frequency in this study appeared to be intermediate (29%). Many biological mechanisms have been implicated in causing segregation distortion including divergence of the parental genotypes [13,25,31], chromosome loss [32], genome size differences [31], genetic load and recessive lethal alleles [33], meiotic drive locus [34], and gametic and zygotic selection [35,36]. In the present study, the female parent Xiaoye was a wild genotype from the forest in Henan province. The male parent Qinzhong No.1 was a cultivar produced by controlled breeding, and it was planted in the museum garden of Northwest A&F University. They differ in traits of growth, phenology, morphology and content of secondary metabolite. Thus, the divergence of the parental genotypes may contribute to the observed segregation distortion.

Genetic Linkage Map
We constructed a genetic linkage map of E. ulmoides based on the segregation of SRAP, AFLP, ISSR and SSR markers as a first step towards understanding the E. ulmoides genome. The total map distance was 2133 cM, and the average map distance between adjacent markers was 3.1 cM. The present map covers a significant portion of the E. ulmoides genome, which should provide adequate coverage of the genome to begin QTL analysis. E. ulmoides is a diploid species with 2n = 34. The number of linkage groups was more than the number of haploid chromosomes of E. ulmoides. The presence of more than 17 linkage groups may be due to some gaps preventing connection between groups belonging to the same chromosome. However, gaps in the genetic linkage maps, resulting in two or more linkage groups per chromosome, are common in tree species even with large numbers of markers [8,11,13,15]. In future work, more co-dominant and functional markers are needed to be added to this genetic linkage map in order to fill the gaps, integrate some linkage groups and cover the entire genome.

QTL Analysis
It is often assumed that a quantitative trait exhibits continuous variation because of the interaction of environmental effects and multiple genes of small and cumulative effects. In the present study, we were able to detect QTLs with moderate to large effect for growth-related traits. The estimated magnitude of the individual QTL effect ranged from 12.4%-33.3% of the phenotypic variance. Our results agree with other QTL studies in tree species indicating that growth-related traits may in part be controlled by a few genes with large effect. In an F 2 population based on an interspecific cross of Populus, Bradshaw and Stettler [37] reported that effects of single QTL for growth-related traits explained 24%-33% of the phenotypic variation. In an interspecific backcross family of white poplar, Zhang et al. [38] found that four QTLs for stem volume explained 35.8% of the total phenotypic variance. In a tetraploid hybrid F 2 population of Salix, most of the QTL for the different growth-related traits each explained around 12% of the phenotypic variation, with a few exceptions explaining more than 20% of the variation [39]. Furthermore, in an intraspecific cross of Salix, 11 QTLs were identified for growth-related traits with each QTL explaining 14%-22% of the phenotypic variance [8].
Four basal diameter QTLs were identified at similar genomic regions as the height QTLs. This was not surprising because of the high correlation coefficient between height and basal diameter. The Pearson correlation coefficients were 0.80 in 2010, 0.83 in 2011, 0.66 in 2012 and 0.72 in 2013, respectively (Table 7). This suggested that height and basal diameter growth in E. ulmoides had common genetic components. The clustering of QTLs controlling highly correlative growth-related traits have been reported in other tree species of Populus [37,38], Salix [8,39], Eucalyptus [18,19,40] and apple [41].
No QTL was consistently expressed over the four years. However, QTL Dht2-1, Dht2-3 and Dbd2-1 in 2012 were identified at similar genomic regions as the QTL Dht3-1, Dht3-3, and Dbd3-1 in 2013, respectively. A similar result of QTL analysis for height and basal diameter in radiata pine was reported by Emebiri et al. [42], who observed that none of the putative QTL positions detected at any one age was strongly expressed at all of the four stages of measurement and that 45% of putative QTLs significant at one age were also detected at a subsequent age. For growth-related traits, QTL instability has been reported frequently in tree species [8,15,37,40,41,43]. Verhaegen et al. [40] did not find the same QTLs over three consecutive years for growth-related traits in hybrid Eucalyptus. In rubber tree, QTLs detected during the summer were different from the QTLs detected during the winter for height and girth growth [15]. To explain this phenomenon, Verhaegen et al. [40] assumed that a set of regulatory genes may differentially control the temporal expression of the genes controlling a trait or that different sets of regulatory factors may be involved during different periods of time. Kenis and Keulemans [41] proposed that genetic control of these traits is largely influenced by environmental factors and probably changes as the tree matures.
To be able to utilize marker-assisted selection successfully in a breeding program, the molecular markers should be consistently found in various environments and show a large effect on the trait. In this study, we have considered only the first four years of a tree's life, and the phenotypic assessment was undertaken in a single environment. Therefore, further QTL analysis under different environmental conditions over the years is necessary for providing additional insights on the pattern and stability of the growth QTLs.

Plant Material
The population consisted of 152 F 1 individuals that resulted from the cross between a wild genotype Xiaoye and a cultivar Qinzhong No.1. Controlled pollination was carried out in the spring of 2009 at Yantuo, Lingbao, Henan, and seeds were collected in October and stored at 4 °C. In March 2010, seeds were sown in a substrate with humus, sand and soil (1:1:1 mix) in plastic cups. Subsequently, seedlings were transplanted to the flat in a greenhouse when they had grown to a height of approximately 20 cm. The progenies were planted in the field in March 2011 at the nursery of Northwest A&F University, Yangling, Shaanxi. The F 1 population was designated as "DZ0901".

DNA Extraction
DNA was extracted from young leaves of the 152 F 1 individuals and the two parental trees according to a modified CTAB procedure [44]. DNA quality was visually assessed on a 1% agarose gel by electrophoresis, and the concentration was determined using a NanoDrop ND-1000 spectrophotometer (NanoDrop Technologies Inc., Wilmington, DE, USA).

SRAP Analysis
SRAP analysis was performed according to Li and Quiros [23] with some modifications. Approximately 50 ng DNA was added to a mixture containing 2.5 mM MgCl 2 , 0.2 mM dNTPs, 0.4 mM of each primer, 1× PCR buffer and 1.5 U Taq DNA polymerase for a total volume of 25 μL. PCR parameters were as follows: 5 min at 94 °C, 5 cycles of 94 °C for 1 min, 35 °C for 1 min and 72 °C for 1.5 min, 30 cycles of 94 °C for 1 min, 50 °C for 1 min and 72 °C for 1.5 min, and a final extension of 10 min at 72 °C. DNA fragments were separated by electrophoresis on 8% non-denaturing polyacrylamide gel and visualized by silver staining. The SRAP primers used in this study are listed in Table 1.

AFLP Analysis
AFLP analysis consisting of genomic DNA digestion with EcoRI and MseI restriction enzymes, adapter ligation, pre-amplification, and selective amplification using EcoRI plus three and MseI plus three selective nucleotide primers were similar to those from Vos et al. [45] with modifications described by Wang et al. [46]. The following cycling parameters were used for pre-amplification: 94 °C for 2 min, 30 Table 2.

ISSR Analysis
The protocols of Zietkiewicz et al. [47] for ISSR were adapted. Reaction mixture was as described above for SRAP except that a single primer was used. Thermal cycling conditions were as follows: 94 °C for 4 min, followed by 38 cycles of 94 °C for 30 s, 45 s at the locus-specific annealing temperature and 72 °C for 1.5 min, and then a final extension step of 72 °C for 5 min. PCR products were detected as described above for SRAP. The 100 primers were from the #9 ISSR primer kit (801-900) of the Biotechnology Laboratory, University of British Columbia (UBC, Vancouver, BC, Canada).

SSR Analysis
The SSR reaction mixture was as described above for SRAP. Thermal cycling conditions were described by Deng et al. [48]: 4 min at 94 °C, locus-specific amplification cycles of 50 s at 94 °C, 50 s at the locus-specific annealing temperature and 90 s at 72 °C, and a final extension step for 10 min at 72 °C. PCR products were detected as described above for AFLP. Nineteen SSR primer combinations developed for E. ulmoides by Deng et al. [48] were used in this study.

Segregation Analysis and Map Construction
Data of segregating markers was analyzed as a "cross-pollinated" population using JoinMap 4.0 [49]. Deviation from expected Mendelian ratio was determined using a chi-square test. The marker placement was determined using a minimum LOD threshold  [50]. To estimate observed genome coverage, the expected genome length of each linkage group was calculated by multiplying the observed length by (m + 1)/(m − 1), where m is the number of markers in that linkage group, and the estimated genome length was the sum of revised length of all linkage groups [51]. Observed genome coverage was assessed by dividing the observed genome length by the estimated genome length.

Growth Traits Assessment and QTL Analysis
Height and basal diameter were measured to evaluate the growth of progenies in October from 2010-2013. The descriptive statistics, the skewness of the distributions and Pearson correlations of traits were calculated using SPSS 13.0 (SPSS Inc., Chicago, IL, USA, 2004) for Windows. QTL analysis was done using MapQTL 5.0 (Plant Research International B.V. and Kyazma B.V., Wageningen, Gelderland, The Netherlands, 2004) [52]. Kruskal-Wallis nonparametric test, interval mapping (IM) and multiple QTL mapping (MQM) were performed for each trait. In MQM, the markers closest to the QTL peaks detected by IM were used as cofactors. The limit of detection (LOD) thresholds was estimated with a 1000-permutation test. The QTLs with LOD values higher than the genome-wide threshold at p < 0.05 were considered significant. However, those QTLs with a LOD score greater than 3 and smaller than the threshold were also reported. The genetic linkage map and QTL positions were drawn using MapChart 2.2 [50].

Conclusions
In this study, we report a genetic linkage map of E. ulmoides constructed by SRAP, AFLP, ISSR and SSR markers. This genetic linkage map provided an adequate coverage of the E. ulmoides genome for QTL analysis. A saturated genetic linkage map will be constructed by adding more co-dominant and functional markers. The QTL analysis provided a better genetic understanding for growth-related traits of E. ulmoides seedlings. Projects have been initiated to use the genetic linkage map to identify QTLs controlling other biological and economically important traits, and this will allow the potential of marker-assisted selection in the improvement of E. ulmoides cultivars.