Development of novel EST-SSR markers for Rhododendron longipedicellatum (Ericaceae) and cross-amplification in two congeners.

PREMISE OF THE STUDY: To investigate the genetic background and population characteristics of Rhododendron longipedicellatum (Ericaceae), a newly discovered and critically endangered species, expressed sequence tag-simple sequence repeat markers were developed, and transferability was tested in two congeners, R. molle and R. simsii. METHODS AND
RESULTS: Based on the transcriptome sequences of R. longipedicellatum, 102 primer sets were designed; 48 primer sets were successfully amplified, with 15 showing polymorphisms in 150 individuals from five extant populations of R. longipedicellatum. The number of alleles per locus ranged from four to 18, and the levels of observed and expected heterozygosity for the 15 loci varied from 0.255 to 0.913 and from 0.306 to 0.851, respectively. All 15 loci were found to amplify in R. molle and R. simsii.
CONCLUSIONS: These polymorphic SSR markers can be used in conservation genetic and phylogeographic studies to elucidate the rarity and origin of R. longipedicellatum.

Rhododendron L., renowned for its horticultural and ecological value, is the largest genus in Ericaceae and is one of the most widespread woody plants in the Northern Hemisphere. Its more than 1025 species are distributed throughout Asia, Europe, and North America, and two species extend to eastern Greenland and Queensland, Australia (Chamberlain et al., 1996; Fang et al., 2005; Cai et al., 2016). Wild Rhododendron species serve as potential genetic resources for the development of new cultivars, and more than 25,000 Rhododendron cultivars have been bred around the world. However, no evergreen rhododendron has a yellow‐flowered cultivar (Ureshino et al., 2016). At present, the breeding of flower color in rhododendrons tends to favor pure‐colored flowers internationally, especially pure yellow rhododendrons (Lan et al., 2012).

Rhododendron longipedicellatum Lei Cai & Y. P. Ma (subg. Rhododendron, sect. Vireya, subsect. Pseudovireya) is an unusual evergreen shrub, with brilliantly pure yellow flowers having no blotches or spots. Unlike all other wild Rhododendron species, whose flowering times occur between March and June in the Northern Hemisphere, the natural flowering time of R. longipedicellatum extends from the last 10‐day period of November to the first 10‐day period of February. However, R. longipedicellatum has a very limited distribution, with only five relict populations found in Malipo County, Yunnan Province, China, and with the largest population comprising about 350 mature plants (Cai et al., 2016). Furthermore, this species is at risk of extinction because of continued disturbance from anthropogenic activities. Therefore, genetic information from R. longipedicellatum is urgently needed for current and future conservation activities.

Expressed sequence tag–simple sequence repeat (EST‐SSR) markers are increasingly used in population genetic studies because they are codominant, multiallelic, and often highly polymorphic, and they are less susceptible to null alleles and homoplasy than anonymous SSRs are (Ellis and Burke, 2007; Yoichi et al., 2016). However, currently, only 77 EST‐SSR markers have been developed in Rhododendron (Yoichi et al., 2016; Xing et al., 2017), and only two markers (Rhob_1022 and Rhob_30843, percentage of polymorphic loci = 2.60%; Table 1) are amplified for R. longipedicellatum, which is insufficient for unraveling the population dynamics of most species within this genus. Therefore, we developed 15 EST‐SSR markers for R. longipedicellatum and evaluated their polymorphism and transferability to R. molle (Blume) G. Don (an important congeneric species with yellow flowers) and R. simsii Planch. (a widespread Rhododendron species that occurs with populations of R. longipedicellatum). These markers will provide an important genetic resource for rhododendron breeding programs worldwide.

Table 1

Genetic diversity of 17 polymorphic EST‐SSR markers (including 15 newly developed markers and two previously developed markers from congeneric species) in five extant populations of Rhododendron longipedicellatum.a

Locus	WBL (N = 30)			WJL (N = 30)			XCW (N = 30)			XL (N = 30)			ZWL (N = 30)			Total	Mean
	A	H o	H e b	A	H o	H e b	A	H o	H e b	A	H o	H e b	A	H o	H e b	A	H o	H e	B
RL6	4	0.233	0.216	4	0.233	0.216	8	0.367	0.377	6	0.367	0.398	5	0.069	0.308***	14	0.255	0.306	0.057
RL16	2	0.233	0.210	4	0.600	0.588	5	0.367	0.463	4	0.267	0.300	3	0.500	0.406	6	0.393	0.410	0.114
RL18	2	0.300	0.305	3	0.172	0.222	2	0.300	0.440*	5	0.733	0.570**	2	0.233	0.381*	6	0.349	0.410	0.742†
RL20	3	0.133	0.128	4	0.367	0.327	2	0.633	0.440**	4	0.700	0.557*	3	0.200	0.186	4	0.407	0.358	0.764†
RL26	4	0.633	0.584	4	0.900	0.715**	5	0.800	0.707	4	0.467	0.444	5	0.667	0.636	5	0.693	0.673	0.094
RL28	5	0.667	0.583	1	0.833	0.657**	4	0.300	0.271	3	0.467	0.501	6	0.467	0.496	8	0.380	0.402	0.115
RL37	5	0.200	0.190	5	0.233	0.303	6	0.500	0.477	6	0.552	0.500	3	0.333	0.294	10	0.362	0.360	0.174†
RL54	7	0.800	0.629**	7	0.667	0.741	6	0.833	0.752	6	0.846	0.657**	7	0.667	0.625	9	0.760	0.796	0.085
RL61	8	0.667	0.799*	10	0.833	0.858	12	0.800	0.807	12	0.933	0.794*	10	0.600	0.743*	18	0.767	0.851	0.054
RL74	6	0.633	0.776*	6	0.733	0.767	7	0.833	0.790	11	0.733	0.738	6	0.367	0.405	13	0.660	0.811	0.103
RL79	3	0.200	0.186	8	0.600	0.753**	6	0.600	0.489*	5	0.867	0.631***	5	0.433	0.428	10	0.540	0.535	0.021
RL89	3	0.167	0.158	7	0.533	0.529	8	0.567	0.612	3	0.133	0.128	4	0.138	0.134	12	0.309	0.345	0.013
RL98	5	0.333	0.402	3	0.500	0.454	4	0.586	0.598	4	0.621	0.537	4	0.367	0.322	7	0.480	0.478	0.031
RL99	5	0.333	0.353	4	0.700	0.577*	5	0.767	0.655*	5	0.308	0.645***	3	0.400	0.346	6	0.507	0.554	0.226†
RL100	7	0.900	0.776*	7	0.967	0.810**	8	0.833	0.775	6	0.889	0.798	10	0.967	0.842*	12	0.913	0.826	0.015
Rhob_1022c	2	0.467	0.370	4	0.667	0.609	5	0.267	0.308	3	0.267	0.246	3	0.667	0.508*	5	0.467	0.425	0.062
Rhob_30843c	7	0.500	0.661*	5	0.800	0.779	8	0.533	0.717**	6	0.714	0.812	5	0.533	0.526	9	0.616	0.730	0.114

A = number of alleles per locus; B = null allele frequency averaged over all populations using the Brookfeld 1 equation (Brookfield, 1996); H e = expected heterozygosity; H o = observed heterozygosity; N = number of individuals analyzed.

aLocality and voucher information are provided in Appendix 1.

bSignificant deviations from Hardy–Weinberg equilibrium after sequential Bonferroni corrections: *P < 0.05, **P < 0.01, ***P < 0.001.

cRhob_1022 and Rhob_30843 were developed by Yoichi et al. (2016).

†Loci with null alleles.

METHODS AND RESULTS

Three mature plants of R. longipedicellatum were collected from population ZWL (voucher specimen accession no. LTQ20160618; Appendix 1) and planted in a greenhouse at the Research Institute of Resources Insects, Chinese Academy of Forestry (Kunming, China). Fresh, tender leaves from the mature plant were gathered 1 y later and mixed in equal proportions for RNA extraction and transcriptome sequencing. Total RNA was extracted with Trizol (Thermo Fisher Scientific, Waltham, Massachusetts, USA) following the manufacturer's instructions. The cDNA library construction and sequencing was performed by staff at the Beijing Genome Institute (Wuhan, China) with a HiSeq 4000 (Illumina, San Diego, California, USA). Altogether, 58.30 Mbp raw reads were obtained and deposited into the National Center for Biotechnology Information (NCBI) sequence read archive (SRA) database (Bioproject ID: SRR6509877). The generated raw reads were filtered to remove reads containing adapters, ambiguous reads (N > 5%), and other low‐quality reads (base quality <15% or >20%), and a total of 44.85 Mbp clean reads were obtained and assembled de novo into 94,906 contigs using Trinity software (Grabherr et al., 2011). TGICL software (Pertea et al., 2003) was used to cluster similar contigs, which generated 74,092 nonredundant unigenes, with an average length of 938 bp. MISA software (Thiel et al., 2003) was used for SSR motif mining from all unigenes, and the minimum numbers of repeats were set as six, five, five, four, and four for di‐, tri‐, tetra‐, penta‐, and hexanucleotide motifs, respectively. Altogether, 20,304 SSR motifs were found, and 102 of them were selected with at least five tri‐ and tetranucleotide repeats, 10 dinucleotide repeats, or four penta‐ and hexanucleotide repeats for primer design in Primer3 software (Rozen and Skaletsky, 1999), with conditions as described by Li et al. (2011).

The EST‐SSR primers were initially screened for performance with two individuals from each of the five relict R. longipedicellatum wild populations (WBL, WJL, XCW, XL, and ZWL; Appendix 1). Genomic DNA was extracted from silica‐dried leaves with a modified cetyltrimethylammonium bromide (CTAB) method (Doyle and Doyle, 1987). PCR amplifications were performed in a final 10‐μL volume, containing 1 μL (10–30 ng) of template DNA, 5 μL of 0.7× Multiplex PCR Master Mix (QIAGEN, Hilden, Germany), 0.5 μL (10 pM) of each primer, and 3 μL of RNase‐free water. The PCR thermal profile consisted of an initial denaturation step at 95°C for 5 min; followed by 35 cycles of denaturation at 94°C for 30 s, annealing at 56–59°C for 30 s (Table 2), and elongation at 72°C for 1 min; with a final elongation step at 72°C for 10 min. All PCR products were resolved by electrophoresis in 1% agarose gels to determine whether amplification was successful. Of the 102 primer pairs, 48 (47.1%) target regions were successfully amplified.

Table 2

Characteristics of the 15 polymorphic EST‐SSR markers developed for Rhododendron longipedicellatum

Locus	Primer sequences (5′–3′)	Repeat motif	Allele size range (bp)	T a (°C)	Fluorescent dyea	BLAST top hit description [organism]	BLAST top hit accession no.	E‐value	GenBank accession no.
RL6	F: GAGCTTCACAAGTTAATTCCCC	(AGC)5	135–161	58	FAM1	5‐methyltetrahydrofolate‐homocysteine methyltransferase reductase mRNA [Monoraphidium neglectum]	XM_014050907.1	5.00E−06	MG585326
	R: ATCATCACCACCACCTTACCC
RL16	F: GTATGCTCTTCCTCCGTCACC	(GAG)6	94–114	59	FAM2	No hit	—	—	MG585327
	R: GAAGAATAATCTGCTGGGATCG
RL18	F: AGGAGCAAGGTATAAAAGCAGC	(AGG)8	85–95	57	FAM1	Uncharacterized LOC100249879 (LOC100249879) [Vitis vinifera]	XM_010648600.1	7.00E−11	MG585328
	R: GGGTTCTTTGTCTTTCTTAGCG
RL20	F: CATGTAGTCGGCCATCTTCC	(CCT)7	100–106	58	FAM2	CCMP1545 predicted protein [Micromonas pusilla]	XM_003062022.1	5.00E−09	MG585329
	R: TGGATCGGAGCTAAGTTCTACC
RL26	F: AGATGAACGTACCGATTAAGGG	(ATG)8	89–105	57	FAM2	No hit	—	—	MG585330
	R: CTCTCTGGTTTTACTGTTCTCTG
RL28	F: CATTGACGTTAAAAGCGATGG	(GTG)7	91–122	57	HEX3	Dof zinc finger protein DOF2.4‐like (LOC104900133) [Beta vulgaris subsp. vulgaris]	XM_010687486.1	4.00E−06	MG585331
	R: TCAGGACTCAGGTGACATATCC
RL37	F: CCAGTCAAGCTGTGTAACTGGT	(CAG)8	90–122	58	HEX2	AP2‐like ethylene‐responsive transcription factor ANT (LOC104591238) [Nelumbo nucifera]	XM_010250039.1	2.00E−08	MG585332
	R: CTCTTCTTTGTTTCCATGCCC
RL54	F: CGACTCCTACTTAAATACAGGGC	(GA)19	84–102	59	HEX3	Vitellogenin‐1 (LOC104803686) [Tarenaya hassleriana]	XM_010527700.1	1.00E−15	MG585333
	R: CTCGCGTGAGAAGACTAGACG
RL61	F: GTAGGTGGGTATATCAATCCTTGC	(CT)12	96–138	59	HEX1	Uncharacterized LOC101310701 (LOC101310701) [Fragaria vesca subsp. vesca]	XM_004299009.2	7.00E−08	MG585334
	R: CTCTATAATCCTTGACACAGGCG
RL74	F: TGTCGCATTTCTGTACACACG	(TC)10	75–129	59	HEX1	Oleracea protein CIA1 (LOC106309548) [Brassica oleracea var. oleracea]	XM_013746597.1	5.00E−06	MG585335
	R: GTATGAGATCTAGGGCACCAGG
RL79	F: CCATTTCATATTTCATCGACGG	(TC)10	110–133	56	ROX1	No hit	—	—	MG585336
	R: AGAGTCGTTGCATAGTCATTGC
RL89	F: GCTCTGTTTCTTGTTATGCTGC	(ATCT)6	102–137	57	ROX2	Cone Tongling01–10 microsatellite sequence [Lycoris radiata]	KP665168.1	8.00E−08	MG585337
	R: CTTCAGTCTAATTAGGGACGC
RL98	F: CCTCTTTACTTCCCATTACCCC	(AAACAA)4	88–103	57	ROX1	Uncharacterized LOC106446472 (LOC106446472) [Brassica napus]	XM_013888210.1	4.00E−06	MG585338
	R: AGCACATTGTTCTTGTTTTCCC
RL99	F: CCCTCTTTACTTCCCATTACCC	(AAACAA)4	85–104	57	ROX3	Uncharacterized LOC106446472 (LOC106446472) [Brassica napus]	XM_013888210.1	4.00E−06	MG585339
	R: AGCACATTGTTCTTGTTTTCCC
RL100	F: CCGTGTCGTAGAGGTTGTTACC	(GTTAGG)4	88–112	59	ROX2	Oryza sativa Japonica group protein Chromatin remodeling 24 (LOC9266414), transcript variant X2 [Oryza sativa]	XM_015781575.1	4.00E−09	MG585340
	R: CCCATAAACACTCCAAACATCC

T a = annealing temperature.

aPCR multiplex sets are indicated as 1, 2, or 3.

PCR fluorescent tagging was performed for further polymorphic screening. The 5′ end of each forward primer for the 48 markers was tagged with one of three fluorescent dyes (FAM, HEX, or ROX [Thermo Fisher Scientific]; Table 2), and multiplex PCR amplifications were performed for the 150 individuals of R. longipedicellatum, representing all extant populations (30 for each population), using the PCR conditions described above. Allele size for the tagged PCR products was obtained using an ABI 3730 sequencer with a GeneScan 500 LIZ Size Standard (Thermo Fisher Scientific) and GeneMapper 4.1 (Thermo Fisher Scientific). Population genetic parameters, including the number of alleles per locus, expected heterozygosity, observed heterozygosity, and deviation from Hardy–Weinberg equilibrium, were analyzed with GENEPOP software (version 3.4; Raymond and Rousset, 1995). Pairwise linkage disequilibrium in each population was tested with FSTAT software (version 2.9.3; Goudet, 1995). The Brookfield 1 equation (Brookfield, 1996) was used to estimate the frequency of null alleles.

Of the 48 candidate markers, 15 (31.3%) exhibited polymorphisms in R. longipedicellatum. Sequences containing those 15 markers were deposited in GenBank (Table 2). The number of alleles per locus ranged from four to 18, with an average of 10 (Table 1). Levels of observed and expected heterozygosity from the polymorphic loci ranged from 0.255 to 0.913 and 0.306 to 0.851, with averages of 0.526 and 0.545, respectively. Of the 15 polymorphic markers, four, five, four, six, and four in the WBL, WJL, XCW, XL, and ZWL populations of R. longipedicellatum, respectively, deviated significantly from Hardy–Weinberg equilibrium after sequential Bonferroni corrections, probably as a result of inbreeding in the remaining population (Li et al., 2018) or other evolutionary factors, and the presence of null alleles was confirmed at four loci (RL18, RL20, RL37, and RL99; Table 1). No significant linkage disequilibrium was observed among the markers. Furthermore, for cross‐species application, the 15 newly developed, polymorphic markers were tested in 30 R. molle and 30 R. simsii individuals from Kunming Botanical Garden (Kunming, China; Appendix 1). All 15 loci were successfully amplified, and only three loci in R. simsii exhibited monomorphism (Table 3).

Table 3

Genetic diversity of 15 polymorphic EST‐SSR markers developed for Rhododendron longipedicellatum in R. molle and R simsii.a

Locus	Rhododendron molle (N = 30)			Rhododendron simsii (N = 30)
	A	H o	H e b	A	H o	H e b
RL6	3	0.179	0.470***	5	0.733	0.766
RL16	2	0.321	0.275	3	0.740	0.577**
RL18	3	0.800	0.549***	1	0.000	0.000
RL20	5	0.600	0.663	3	0.800	0.533***
RL26	4	0.500	0.421	2	0.640	0.444**
RL28	6	0.733	0.733	2	0.519	0.391*
RL37	5	0.600	0.646	1	0.000	0.000
RL54	7	0.867	0.724*	5	0.733	0.771
RL61	7	0.633	0.675	4	0.846	0.576***
RL74	6	0.733	0.689	5	0.733	0.750
RL79	4	0.333	0.615***	5	0.600	0.733*
RL89	3	0.692	0.532**	7	0.333	0.670***
RL98	4	0.633	0.599	1	0.000	0.000
RL99	3	0.267	0.540***	1	0.333	0.414
RL100	5	0.700	0.616	6	0.767	0.646*

A = number of alleles per locus; H e = expected heterozygosity; H o = observed heterozygosity; N = number of individuals analyzed.

aLocality and voucher information are provided in Appendix 1.

bSignificant deviations from Hardy–Weinberg equilibrium after sequential Bonferroni corrections: *P < 0.05, **P < 0.01, ***P < 0.001.

CONCLUSIONS

We developed 15 highly informative EST‐SSR markers for R. longipedicellatum, which can be used in population genetic diversity, genetic structure, and phylogeographic studies to facilitate development of scientific conservation measures in R. longipedicellatum. The markers may also be valuable for population and evolutionary studies of congeneric species and closely related taxa.

Species	Population code	N	Collection localitya	Geographic coordinates	Altitude (m)	Voucher no.b
R. longipedicellatum Lei Cai & Y. P. Ma	WBL	30	Malipo, Yunnan	23°09′33.7″N, 104°56′48.7″E	1312	MH20150614
	WJL	30	Malipo, Yunnan	23°09′52.4″N, 104°56′34.8″E	1316	LTQ20161205
	XCW	30	Malipo, Yunnan	23°09′47.3″N, 104°56′45.1″E	1248	MH20141124
	XL	30	Malipo, Yunnan	23°10′1.9″N, 104°56′51.1″E	1183	DZL3637‐1
	ZWL	30	Malipo, Yunnan	23°09′59.8″N, 104°56′22.0″E	1270	LTQ20160618
R. molle (Blume) G. Don	YZC	30	Kunming, Yunnan	25°08′24.6″N, 102°44′27.9″E	1953	LXF20170322
R. simsii Planch.	YSH	30	Kunming, Yunnan	25°08′25.0″N, 102°44′31.6″E	1951	LXF20170615

N = number of individuals sampled.

aCollection localities in China.

bVoucher specimens are deposited at the Herbarium of the Kunming Institute of Botany, Chinese Academy of Sciences (KUN), Kunming, China.