劉夢露，戴亮芳，程馬龍，張 金，李海燕，謝建坤，羅向東

(江西師范大學 生命科學學院，南昌 330022)

RhododendronmolleG. Don (Ericaceae) is one ofRhododendronspecies native to China and Japan. With flourishing branches and gorgeous flowers, it is one of elite parents for new varieties improvement inRhododendron[1]. Moreover, the plant ofR.molleis rich in ericol-in, rhodofoxin and rhodojaponin-Ⅲ. Rhodojaponin-Ⅲ is one of grayanane diterpenoids, possesses significant blood pressure lowering and heart rate slowing effects[2]. Both of ancient and modern medicine studies suggested that different organs and tissues ofR.molle(including flower, fruit, root and leaf ect.) could be used as analgesics and insecticides[2-4]. In recent years, the new research on the pharmacological effects ofR.mollehas extended to anti-schistosome[5]. Therefore,R.molleis a high valuable medicinal and attractive ornamental plants.

R.mollewas widely distributed across central and south China until the 1980s. As a result of overexploitation, habitat loss, weak fertility and lower ability of natural regeneration,R.mollehas dramatically decreased in recent years. It has suffered severe habitat fragmentation and population loss. At present,R.molleexists in a restricted number of populations with small population sizes. The remnant populations are exposed to higher risk of critically endangered or extinction[6]. Understanding the processes that determine population genetic structure and diversity is very important for the conservation and management of threatened populations, which is necessary to employ more powerful DNA markers[7]. Simple sequence repeats (SSR) or microsatellites marker is a co-dominant marker, which is regarded as one of the most effective markers to investigate the genetic diversity and population genetic structure[8].

Here, we report a set of novel polymorphic microsatellites markers for endangeredR.molleby using Restriction site Associated DNA approach and Illumina paired-end sequencing (RAD-seq), which will provide a useful tool for studying population genetics and evolutionary history as well as developing a conservation strategy for this species.

1 Materials and methods

1.1 Plant material and DNA extraction

The young leaves ofR.mollewere collected from three locations of China, Yongxiu County (29°07′N, 115°43′E, Jiangxi Province), Panan County (28°52′N, 120°26′E, Zhejiang Province) and Jinzhai County (31°23′N, 115°46′E，Anhui Province). The young leaf tissues were sampled and flash frozen in liquid nitrogen and stored at -80 ℃ until DNA extraction. DNA was extracted from the young leaf tissues following modified CTAB method[9]. The quality and concentration of the DNA were examined by electrophoresis on 1% agarose gels with D2000 DNA marker.

1.2 Methods

1.2.1RADlibraryconstructionandIlluminasequencingThe genomic restriction site associated DNA sequences (RAD-seq) libraries were created by Novogene company (Novogene，China) as described by Bairdetal.[10]. Genomic DNA ofR.mollewas digested with the restriction endonuclease (EcoR I), sheared randomly. Illumina P1and P2adaptors were ligated after end repair and addition of single ‘A’ base. After ligation, 1.0×SPRI (<100 bp) cleanup was performed and the library was amplified for enrichment of adapter ligated fragments. Then a single library with an insert size range of 300-700 bp was prepared and sequenced from both ends with 100 bp read lengths in one lane of an Illumina sequencer on the Illumina HiSeqTM2000/MiseqTM. The raw data from the images acquired after sequencing was transformed by base calling into raw reads and stored in FASTQ format. Raw reads were subjected to quality check by usingSeQCv2.1 (Genotypic Proprietary Tool)[11]. Clean reads were obtained after sequencing error rate distribution checks, GC content distribution checks and sequencing data filtering.

1.2.2AssemblyofRAD-seqandSSRmarkersdiscoveryA total of 7.513G bp clean reads was partially assembled using VelvetOptimiser to obtain a final assembly sequence[12]. The filtered assembled contigs were used as a source of sequence information for SSRs screening by SR search software with 2 and 6 as the minimum and maximum length of the repeating unit, respectively. The minimum distance between two SSR sequence is 12 bp and SSR upstream and downstream sequence length of 100 bp. SSR primers were designed using Primer 3.0[13]. And the major parameters for designing the primers were: (1) primer sequence length from 20 to 28 nucleotides, with 24 as the optimum, (2) annealing temperature from 60 ℃ to 65 ℃, an optimum annealing temperature of 63 ℃, and (3) a pair of primers annealing temperature difference of the maximum of 1 ℃.

1.2.3IdentificationofpolymorphismSSRanddataanalysis128 pairs of primers selected randomly from 11 687 pairs of primers were synthesized by Sangon Company (Shanghai). Polymerase chain reaction (PCR) amplifications in six different geographicalR.molleindividuals were performed to identify polymorphism microsatellite loci. Each PCR reaction consisted of 4 μL of ultrapure water, 7.5 μL of 2×Taq PCR Green Mix, 2 μL each of forward and reverse primers (0.5 μmol·L-15 pmol), and 1.5 μL of template genomic DNA in a final reaction volume of 15 μL. Annealing temperature gradient from 55 ℃ to 60 ℃ was set for finding optimum annealing temperature. The thermal cycles for the SSR PCR runs were programmed for an initial denaturation at 94 ℃ for 3 min; followed by 34 cycles of 94 ℃ for 45 s, 56 ℃ for 45 s, 72 ℃ for 1 min and finally, 72 ℃ for a final extension of 7 min. The PCR product was then separated on an 8% non-denaturing polyacrylamide gel and visualized by silver staining. 20 bp DNA ladder was used as standard for scoring. The characteristics and availability of polymorphic primers were further evaluated in 63 individuals from three different populations. The evaluation parameters include the number of alleles per locus (A), observed heterozygosity (Ho), expected heterozygosity (He) as well as deviations from Hardy-Weinberg equilibrium (HWE), which were calculated by the GenAlEx 6.1[14]. Null allele frequency was checked by Micro-Checker 2.2[15].

2 Results and analysis

2.1 RAD-seq for R. molle

A FASTQ format file size of 7.653G bp raw reads was generated after Illumina sequencing (Table 1). After trimming adapter sequences, containing ‘N’ and low quality bases, the initial clean data (7.513G bp) were obtained. The error rate of clean reads was about 0.04%. The Q20and Q30were to 93.73% and 88.33%, respectively (Table 1). The percentage of GC content observed inR.mollewas about 39.69% and RAD-Tag digestion ratio is to 98.48% (Table 1). Clean reads were assembleddenovoto serve as framework reference genome. About 171.534 Mb ofR.mollegenome distributed over 498 252 contigs were obtained upon assembly of sequencing data. The contig lengths ranged from 100 to 952 bp with an N50length of 470 bp (Table 1, Fig. 1). These results suggested that the RAD-seq sequencing was successful, which was attributed to the sufficient samples data, normal GC distribution and high sequencing quality in this study.

2.2 SSR discovery and evaluation for R. molle

After several filtering criteria, 349 912 of 498 252 contigs were selected to detect SSR loci, since the contigs of less than 212 bp were filtered. Approximately 3.27% of all assembled contigs (11 457 out of 349 912) contained at least one detectable SSR locus and the total number of microsatellite loci was 11 961. Finally, a set of 11 687 simple sequence repeats with primers were obtained from 11 961 microsatellite loci by Primer 3.0 software. We found that the proportion of five different SSR units was significant different. Different repeat units occurred at numbers of di-, tri-, tetra-, penta- and hexa-repeats were 6 049 (51.76%), 4408 (37.72%), 695 (5.94%), 369 (3.16%), 166 (1.42%), respectively (Table 2).

Table 1 Details of the raw reads’ quality for RAD-seq in R. molle

Fig.1 Distribution of contig length (>212 bp)

A total of 481 types of repeat motifs were identified among the 11 687 (SSRs length≥12 bp) (Table 2). Of the different perfect SSR motif units, the most type number was penta-nucleotide (32.43%), followed by tetra-nucleotide (26.82%), hexa-nucleotide (25.79%), tri-nucleotide (12.47%) and di-nucleotide (2.49%). Further analysis indicated that the number of different repeat motifs in different repeat units was distributed unevenly. The (GA)nrepeat motif was the most abundant type in di-nucleotide repeat motifs. The four other main unit type were the (AAG)nin tri-nucleotide, (AAAT)nin tetra-nucleotide, and (AAAAG)nin penta-nucleotide, (TCTGAG)n/(AAAAAG)nin hexa-nucleotide repeat motifs, and occurred at frequencies of 19.42%, 5.20%, 12.09%, 4.07% and 3.01%, respectively.

2.3 Characteristics of polymorphic microsatellites markers

128 microsatellite markers were selected randomly to detect the polymorphism among sixR.molleindividuals from three different provinces, of which 20 pairs of primer present high polymorphism. To test and characterize the genetic variability ofR.molleat the 20 loci, 63 individuals were collected from three locations of China (Yongxiu County, Panan County and Jinzhai County). All the SSR bands were clear and distinct (Fig. 2).

The characteristics of the polymorphic SSR markers ofR.mollewere shown in Table 3. The number ofAvaried from 4 to 16 with an average of 9.2 alleles per locus,HoandHevalues ranged from 0.050 to 1.000 and 0.489 to 0.908 (Table 3). Seven loci (HD26, HD34, HD35, HD55, HD62, HD88, and HD100) showed a significant departure from HWE (P< 0.05) probably due to heterozygote deficiency in all populations. The rest of thirteen loci were not deviated from HWE in at least one population. All loci were successfully amplified for three populations ofR.molle, which suggested that these genomic-SSR markers will be valuable for genetic research inR.molle, such as population genetic structure analysis, diversity analysis and conservation etc.

Table 2 Number and frequency of five perfect SSR repeat motifs in R. molle

3 Discussion

Due to the advantage of high abundance, random distribution within the genome, high polymorphism information and co-dominant inheritance, SSR markers have been widely used for construction of genetic maps, diversity analysis and genetic structure[8]. There are a number of traditional ways to develop SSR, including microsatellite enrichment, selectively amplified microsatellite polymorphic loci (SAMPL), and EST SSRs or genomic microsatellites according to public DNA databases[16]. All of these ways, however, should obtain the detailed sequence information about both sides of the repeat sequence and then design primers, which are costly and time-consuming. Recently, RAD-seq was developed and regarded as a promising and powerful tool to acquire the genomic and transcriptomic data for non-model organisms and non-sequenced genomes, which has been used for genome assembly and SSR marker development in a variety of species[17]. In the present study, 11 687 microsatellites were developed using RAD-seq. Then 20 polymorphic microsatellite loci from 128 microsatellites were identified. Compared with previous studies, the RAD-seq for SSR development forR.mollewithout an available reference genome shows more time-saving, high efficient and high throughput[18].

In this study, the RAD-seq data was reliably called for 93.70% (sequence coverage at least 1X) and average depth of coverage of 26.16X. It allowed the discovery of a set of 11 961 perfect SSR loci of ≥12 bp length inR.molle. About 3.3% contigs inR.mollepossessed at least one SSR, which was lower than previously reported SSR prevalence in the ESTs of pigeonpea (7.6%), wheat (7.4%) and flax (3.5%), higher than barley (2.8%)[19-22]. The SSR frequency also depended on the parameters used in detecting SSR markers, such as the repeat length and number of repeat unit thresholds. The abundance of SSR (SSR/kbp) inR.mollewas 14.3 compared to 20.0 in cotton[17], 14.0 inArabidopsis, 8.4 in pigeonpea and 5.4 in wheat[19-20]. Our result showed the abundance of SSR existed inR.mollewas moderate.

SSRs are microsatellites with 2-6 nucleotides tandem repeats which are distributed randomly throughout the genome of all the eukaryotes[11]. In the present study, we found that the proportion of different SSR units was significant different. The di-nucleotides were the most abundant (51.76%) repeat sequences, which was agree with other species, for example, di-nucleotides account for 48 to 67% of repeat sequences in plant libraries[23], e.g. pigeonpea (60.41%)[19],P.alata(62.5%)[24]. In ESTs of grape and soybean, however, the tri-nucleotides SSR locus was the most abundant class[25]. Yuetal.[26]reported that the di-nucleotide repeats were mainly in the noncoding regions, whereas a portion of tri-nucleotide repeats were in the coding region. In this study, RAD-seq is characterized high throughput and whole-genome scanning, which resulted in RAD-seq tags include large amounts of information representing noncoding regions and the di-nucleotides were the most abundant (51.76%) repeat sequences inR.molle.

Out of all the repeat motifs from 11 687 primer pairs, (GA)n(10.05%) was the most abundant di-nucleotide repeat in the samples, which was in agreement with the previous results obtained in jute[8], barley, maize, rice, sorghum, and wheat[27]. But Tanetal.[28]reported that TG/CA (9.2%) was the most common repeat motif for bermudagrasses. InR.molle, CTT/AAG, TCT/AGA and AAT/ATT accounted for 4.31%, 4.15% and 4.04% of tri-nucleotide repeats. Similarly, in jute the most abundant trimeric motifs were AAG/CTT, AAC/GTT, TCT/AGA[8], and inP.alata, AAG/CTT, AGA/TCT and GAA/TTC were the most abundant motifs[24]. Analysis of the most abundant di-nucleotide and tri-nucleotide repeat, we finnd that most of these repeats are rich in A/T which is prone to mutation during the SSR evolution[29]. SSR plays an important role in species formation due to the effect to regulate of gene expression[29]. An evolutionary view, natural selection of species and organism adapt to the environment result in the difference of repetitive sequences. Higher the evolution degree of biological was, more the repetitive sequences accounted for the proportion of the genome[30]. Therefore, some kind of SSR, such as GA, TG or CTT/AAG, AGA/TCT，their abundance in a species perhaps due to natural selection and the evolution of the species. In addition, it is likely that these types of repeats have already been used to design probes to capture DNA fragments containing microsatellites.

To our knowledge, this is the first report of RAD-seq to develop more than ten thousand of SSR markers ofR.mollewhich will be highly valuable in the diversity and genetic structure analysis. The 20 primer pairs will be used in future investigation of spatial genetic structure, comparative levels of genetic variation within and between populations. We expect that these markers could be utilized to facilitate efficient conservation and management strategies for this medicinal and ornamental plant. On the other hand, these developed SSR primers can be used as universal primers of the related species[31].

In this article, a dataset of 498 252 contigs were derived from 7.513G bp clean reads ofR.molleusing next-generation DNA sequencing technique. A comprehensive set of 11 687 genomic-SSR markers was developed. 128 microsatellite markers were selected randomly to detect the polymorphism. 20 out of 128 pairs of polymorphic primers were identified in 63 individuals from three different populations. The number of alleles per locus ranged from 4 to 16, expected heterozygosity (He) values ranged from 0.489 to 0.908, respectively. Our results further confirm that the RAD-seq method is an efficient and cost-effective means for SSR discovery. And these polymorphic microsatellite makers developed in this study will be useful for further investigation of population structure and diversity ofR.molleand other congener species.

:

[1] ZHONG G H, LIU J X, WENG Q F,etal. Laboratory and field evaluations of rhodojaponin-III against the imported cabbage wormPierisrapaeL. (Lepidoptera: Pieridae) [J].PestManagementScience, 2006,62: 976-981.

[2] CHENG D M, HU M Y. Research advances in the research onRhododendronmolleG. Don [J].NaturalProductResearchandDevelopment, 1999,11: 109-113.

[3] BAO G H, WANG L Q, CHENG K F,etal. Diterpenoid and phenolic glycosides from the roots ofRhododendronmolle[J].PlantaMedica, 2003,69: 434-439.

[4] ZOU H Y, LUO J, XU D R,etal. Tandem solid-phase extraction followed by HPLC-ESI/QTOF/MS/MS for rapid screening and structural identification of trace diterpenoids in flowers ofRhododendronmolle[J].PhytochemicalAnalysis, 2014,25: 255-265.

[5] 谷旭田, 趙紅梅, 陳洪亮. 羊躑躅等6種中藥對釘螺與血吸蟲毛蚴殺滅效果觀察[J]. 中獸醫(yī)醫(yī)藥雜志. 2008,3: 20-22.

GU X T, ZHAO H M, CHEN H L. Study on the effects of six kinds of Chinese herbs on killingOncomelaniaand blood fluke miracidium [J].JournalofTraditionalChineseVeterinaryMedicine, 2008,3: 20-22.

[6] 程馬龍, 戴亮芳, 羅向東, 等. 瀕危羊躑躅自然居群的遺傳多樣性和遺傳結構分析[J]. 西北植物學報, 2016,36(3): 674-680.

CHENG M L, DAI L F, LUO X D,etal. Genetic diversity and genetic structure analysis on natural populations of endangeredRhododendronmolleG. Don[J].ActaBotanicaBoreail-OccidentaliaSinica, 2016,36(3): 674-680

[7] TAN X X, LI Y, GE X J. Development and characterization of eight polymorphic microsatellites forRhododendronsimsiiPlanch (Ericaceae) [J].ConservationGenetics, 2009,10: 1 553-1 555.

[8] DAS M, BANERJEE S, DHARIWAL R,etal. Development of SSR markers and construction of a linkage map in jute [J].JournalofGenetics, 2012,91: 21-31.

[9] DOYLE J J. DNA protocols for plants-CTAB total DNA isolation. In: G. M. HEWITT and A. JOHNSTON, ed.MolecularTechniquesinTaxonomy[M]. Berlin: Springer-Verlag. 1991, 283-293.

[10] BAIRD N A, ETTER P D, ATWOOD T S,etal. Rapid SNP discovery and genetic mapping using sequenced RAD markers [J].PLoSOne, 2008,3: 3 376-3 376.

[11] MUDALKAR S, GOLLA R, GHATTY S,etal. De novo transcriptome analysis of an imminent biofuel crop,CamelinasativaL. using Illumina GAIIX sequencing platformand identification of SSR markers [J].PlantMolecularBiology, 2014,84: 159-171.

[12] LI H, DURBIN R. Fast and accurate short read alignment with Burrows-Wheeler transform [J].Bioinformatics, 2009,25: 1 754-1 760.

[13] ZERBINO D R, BIRNEY E. Velvet: algorithms for de novo short read assembly using de Bruijn graphs [J].Genomeresearch, 2008,18: 821-829.

[14] PEAKALL R, SMOUSE P E. GENALEX 6: Genetic analysis in Excel. Population genetic software for teaching and research [J].MolecularEcologyNotes, 2006,6: 288-295.

[15] OOSTERHOUT O C, HUTCHINSON W F, WILLS D P M,etal. MICRO-CHECKER: Software for identifying and correcting genotyping errors in microsatellite data [J].MolecularEcologyNotes, 2004,4: 535-538.

[16] VARSHNEY R K, GRANER A, SORRELLS M E. Genic microsatellite markers in plants: features and applications [J].TrendinBiotechnology, 2005,23: 48-55.

[17] WU K, LIU H Y, YANG M M,etal. High-density genetic map construction and QTLs analysis of grain yield-related traits in Sesame (SesamumindicumL.) based on RAD-Seq techonology [J].BMCPlantBiology, 2014,14: 1-14.

[18] ABDELKRIM J, ROBERTSON B, STANTON J A,etal. Fast, cost-effective development of species-specific microsatellite markers by genomic sequencing [J].BioTechniques, 2009,46: 185-192.

[19] DUTTA S, KUMAWAT G, SINGH B P,etal. Development of genic-SSR markers by deep transcriptome sequencing in pigeonpea [Cajanuscajan(L.) Millspaugh] [J].BMCPlantBiology, 2011,11: 1-13.

[20] PENG J H, LAPITAN N L. Characterization of EST-derived microsatellites in the wheat genome and development of eSSR markers [J].Functional&IntegrativeGenomics, 2005,5: 80-96.

[21] CLOUTIER S, NIU Z, DATLA R,etal. Development and analysis of EST-SSRs for flax (LinumusitatissimumL.) [J].Theoretical&AppliedGenetics, 2009,119: 53-63.

[22] VARSHNEY R K, GROSSE I, HAHNEL U,etal. Genetic mapping and BAC assignment of EST-derived SSR markers shows nonuniform distribution of genes in the barley genome [J].TheoreticalandAppliedGenetics, 2006,113: 239-250.

[23] LI Y C, KOROL A B, FAHIMA T,etal. Microsatellites: genomic distribution, putative functions and mutational mechanisms: a review [J].MolecularEcology, 2002,11: 2 453-2 465.

[24] PENHA H A, PEREIRA G D S, ZUCCHI M I,etal. Development of microsatellite markers in sweet passion fruit, and identification of length and conformation polymorphisms within repeat sequences [J].PlantBreeding, 2013,132: 731-735.

[25] HISANO H, SATO S, ISOBE S,etal. Characterization of the soybean genome using EST-derived microsatellite markers [J].DNAResearch, 2007,14: 271-281.

[26] YU J K, DAKE T M, SINGH S,etal. Development and mapping of EST-derived simple sequence repeat (SSR) markers for hexaploid wheat [J].Genome, 2004,47: 805-818.

[27] KANTETY R V, ROTA M L, MATHEWS D E,etal. Data mining for simple sequence repeats in expressed sequence tags from barley, maize, rice, sorghum and wheat [J].PlantMolecularBiology, 2002,48: 501-510.

[28] TAN C C, WU Y Q, CHARLES M,etal. Development of simple sequence repeat markers for bermudagrass from its expressed sequence tag sequences and preexisting sorghum SSR markers [J].MolecularBreeding, 2012,29: 23-30.

[29] KATTI M V, RANJEKAR P K, GUPTA V S. Differential distribution of simple sequence repeats in eukaryotic genome sequences [J].MolecularBiologyandEvolution, 2001,18: 1 161-1 167.

[30] 郭瑞星, 劉小紅, 榮廷昭, 等. 植物SSR標記的發(fā)展及其在遺傳育種中的應用[J]. 玉米科學, 2010,45: 457-459.

GUO R X, LIU X H, RONG T Z,etal. Development of SSR maker and its application in plant genetics and breeding [J].JournalofMaizeSciences, 2005,13: 8-11.

[31] QIN Z C, YANG J B. Development and characterization of 15 microsatellite loci forRhododendrondelavayiFranch. (Ericaceae) [J].Hortscience, 2010,45: 457-459.