Next-Generation Genome Sequencing of Sedum plumbizincicola Sheds Light on the Structural Evolution of Plastid rRNA Operon and Phylogenetic Implications within Saxifragales.

The genus Sedum, with about 470 recognized species, is classified in the family Crassulaceae of the order Saxifragales. Phylogenetic relationships within the Saxifragales are still unresolved and controversial. In this study, the plastome of S. plumbizincicola was firstly presented, with a focus on the structural analysis of rrn operon and phylogenetic implications within the order Saxifragaceae. The assembled complete plastome of S. plumbizincicola is 149,397 bp in size, with a typical circular, double-stranded, and quadripartite structure of angiosperms. It contains 133 genes, including 85 protein-coding genes (PCGs), 36 tRNA genes, 8 rRNA genes, and four pseudogenes (one ycf1, one rps19, and two ycf15). The predicted secondary structure of S. plumbizincicola 16S rRNA includes three main domains organized in 74 helices. Further, our results confirm that 4.5S rRNA of higher plants is associated with fragmentation of 23S rRNA progenitor. Notably, we also found the sequence of putative rrn5 promoter has some evolutionary implications within the order Saxifragales. Moreover, our phylogenetic analyses suggested that S. plumbizincicola had a closer relationship with S. sarmentosum than S. oryzifolium, and supported the taxonomic revision of Phedimus. Our findings of the present study will be useful for further investigation of the evolution of plastid rRNA operon and phylogenetic relationships within Saxifragales.

1. Introduction

The genus Sedum comprises more than 420 recognized species, which is the most species-rich member of the family Crassulaceae [1,2]. Some species, formerly classified as Sedum, are now assigned to the segregate genera Hylotelephium and Rhodiola [3,4,5,6]. The family Crassulaceae, together with 14 other family members, has been classified in the order Saxifragales. Recently, increasing research efforts have been focused on the evolution of Saxifragales, however, phylogenetic relationships within the order are still unresolved due apparently to a rapid, ancient radiation [7,8,9,10,11,12].

Extensive genes were transferred from plastids to the nucleus during evolution. In most land plants, the plastid genome (plastome) is a circular biological macromolecule with a typical quadripartite structure [13,14,15]. In higher plants, compared with hundreds or thousands of tandem repeats in nuclear ribosomal RNA genes [16,17,18], typical plastid rRNA genes are characterized by a pair of inverted rrn operons, which show the gene order of rrn16, rrn23, rrn4.5, and rrn5 [19]. With the rapid development of next generation genome sequencing, more and more complete plastid genomes (plastomes) have been deposited in a public database. Recently, a total of over 3000 reference sequences of plastomes were available in GenBank. The plastomes have been widely accepted as a popular tool for phylogenetic studies [7,20,21,22,23,24,25,26].

Thus far, 43 plastomes have been obtained in Saxifragales, as shown in Table 1. Currently, there is no report of a plastome for S. plumbizincicola, a well-known Zn/Cd hyperaccumulator, which was newly discovered from lead and zinc mining areas in Zhejiang province, China [27]. In this paper, we sequenced the plastome of this species using next-generation genome sequencing. Together with the public sequences, we performed a comparative analysis of plastomes within Saxifragales. Consequently, the aims of this research were (1) to investigate general features of the S. plumbizincicola plastome, (2) to examine the structural evolution of the plastid ribosomal RNA operon, and (3) to clarify phylogenetic relationships within the order Saxifragales.

Table 1

Species of plastomes examined in this study.

Family	Species	Accession No.	Reference
Altingiaceae	Liquidambar formosana	NC_023092.1	[7]
Cercidiphyllaceae	Cercidiphyllum japonicum	NC_037940.1	[37]
Crassulaceae	Phedimus kamtschaticus	NC_037946.1	[38]
Crassulaceae	Phedimus takesimensis	NC_026065.1	Unpublished
Crassulaceae	Rhodiola rosea	NC_041671.1	[39]
Crassulaceae	Sedum oryzifolium	NC_027837.1	Unpublished
Crassulaceae	Sedum plumbizincicola	MN185459.1	This study
Crassulaceae	Sedum sarmentosum	NC_023085.1	[7]
Daphniphyllaceae	Daphniphyllum oldhamii	NC_037883.1	[8]
Grossulariaceae	Ribes fasciculatum	MH191388.1	[8]
Haloragaceae	Myriophyllum spicatum	NC_037885.1	[8]
Hamamelidaceae	Chunia bucklandioides	NC_041163.1	[40]
Hamamelidaceae	Corylopsis coreana	NC_040141.1	[41]
Hamamelidaceae	Fortunearia sinensis	NC_041487.1	[42]
Hamamelidaceae	Hamamelis mollis	NC_037881.1	[8]
Hamamelidaceae	Loropetalum subcordatum	NC_037694.1	[43]
Hamamelidaceae	Parrotia subaequalis	NC_037243.1	Unpublished
Hamamelidaceae	Sinowilsonia henryi	NC_036069.1	Unpublished
Paeoniaceae	Paeonia brownii	NC_037880.1	[8]
Paeoniaceae	Paeonia decomposita	NC_039425.1	[44]
Paeoniaceae	Paeonia delavayi	NC_035718.1	[45]
Paeoniaceae	Paeonia jishanensis	MG991935.1	[46]
Paeoniaceae	Paeonia lactiflora	NC_040983.1	[47]
Paeoniaceae	Paeonia ludlowii	NC_035623.1	[45]
Paeoniaceae	Paeonia obovata	NC_026076.1	Unpublished
Paeoniaceae	Paeonia ostii	NC_036834.1	Unpublished
Paeoniaceae	Paeonia rockii	NC_037772.1	[48]
Paeoniaceae	Paeonia suffruticosa	NC_037879.1	[8]
Paeoniaceae	Paeonia veitchii	NC_032401.1	Unpublished
Penthoraceae	Penthorum chinense	NC_023086.1	[7]
Iteaceae	Itea chinensis	NC_037884.1	[8]
Saxifragaceae	Bergenia scopulosa	NC_036061.1	[49]
Saxifragaceae	Chrysosplenium aureobracteatum	NC_039740.1	[50]
Saxifragaceae	Heuchera parviflora	KR478645.1	[51]
Saxifragaceae	Heuchera richardsonii	NC_042923.1	Unpublished
Saxifragaceae	Heuchera villosa	NC_042924.1	Unpublished
Saxifragaceae	Mitella diphylla	NC_042925.1	Unpublished
Saxifragaceae	Mitella formosana	NC_042926.1	Unpublished
Saxifragaceae	Mukdenia rossii	NC_037495.1	Unpublished
Saxifragaceae	Oresitrophe rupifraga	NC_037514.1	[52]
Saxifragaceae	Saxifraga stolonifera	NC_037882.1	[8]
Saxifragaceae	Tiarella cordifolia	NC_042927.1	Unpublished
Saxifragaceae	Tiarella polyphylla	NC_042928.1	Unpublished
Saxifragaceae	Tiarella trifoliata	NC_042929.1	Unpublished

2. Materials and Methods

2.1. Sample Collection and DNA Extraction

The fresh leaf samples of S. plumbizincicola (code AHNU-KPBK001) were collected from Panjiacun (29°35′16″ N, 118°35′19″ E) in Zhejiang Province, east China. Genomic DNA extraction was conducted using the Plant Genomic DNA kit (Tiangen, Beijing, China), following the manufacturer's instructions. The library was constructed using a TruSeq DNA PCR-Free Library Prep Kit (Illumina, San Diego, CA, USA) and sequenced on the Illumina Hiseq X Ten (Illumina, San Diego, CA, USA) with the strategy of 150 paired-ends and an insert size of 350 bp.

2.2. Genome Assembly, Gene Annotation, and Sequence Analyses

The paired-end reads were first checked with Fastqc [28] and then trimmed for quality using Trimmomatic 0.39 [29]. After that, obtained clean reads were filtered and assembled with GetOrganelle 1.5.2 [30] using the chloroplast genome of S. sarmentosum [7] as reference. The chloroplast genome was annotated with GeSeq [31]. The secondary cloverleaf structures of tRNAs were identified using tRNA-scan SE web server [32]. The secondary structures of rRNAs were predicted by comparison with those of other plant species [33].

2.3. Phylogenetic Analysis

To resolve the phylogenetic relationships among Saxifragales species, two phylogenetic approaches were applied: the maximum likelihood (ML) method in RAxML GUI 1.5b2 [34], as well as the Bayesian inference (BI) method in MrBayes 3.2.7a [35]. With exclusion of the termination codons, 79 protein-coding genes (PCGs) and 4 rRNAs of 37 Saxifragales species were used to construct an evolutionary tree. A phylogenomic study by Yang et al. [36] revealed a sister group relationship between Saxifragales and Rosids. We selected therefore two Vitales species within Rosids (Vitis heyneana, NC_039796; V. vinifera, NC_007957) as outgroups.

For ML analyses, we performed analyses with thorough bootstrap for ten runs and 1000 replicates under the GTRCAT model using RAxML GUI. For BI analyses, the best-fit models for 83 genes were first selected based on Bayesian information criterion (BIC) values in ModelGenerator 0.85 [53], then two simultaneous runs with eight independent Markov chains were run for 10,000,000 generations (sampling every 1000 generations).

3. Results and Discussion

3.1. General Features of S. plumbizincicola Plastome

Based on Bowtie2 mapping, in total 19,610,999 reads (21.5% of total reads) were mapped to the reference genome (S. sarmentosum, NC_023085), with a 1969× mean coverage (min, 1286×, max, 3664×, standard deviation, 71). The assembled complete plastome of S. plumbizincicola (accession number: MN185459.1) is 149,397 bp in size, with a typical circular, double-stranded, and quadripartite structure of angiosperms. The plastome has two identical inverted repeats (IRs, 25,565 bp) separated by a small single copy (SSC, 16,669 bp) and a large single copy (LSC, 81,598 bp), as shown in Figure 1. Approximately 52.0%, 4.3%, and 1.83% of the genome encodes for proteins, rRNAs, and tRNAs, respectively. Whereas, the remaining 41.87% are non-coding regions, including introns, intergenic spaces, and pseudogenes.

Figure 1

Chloroplast genome annotation map for Sedum plumbizincicola. Genes lying outside the circle are transcribed in a clockwise direction, whereas genes inside are transcribed in a counterclockwise direction. Different colors represent different functional groups. The dashed darker and lighter gray in the inner circle denote G + C and A + T contents of chloroplast genome, respectively. LSC, SSC, and IRs mean long single copy, small single copy, and inverted repeat regions, respectively.

Along with new data from this study, we comparatively investigated the structures and properties of plastomes from 44 species, representing 11 families in Saxifragales, as shown in Table 1. The size of plastomes of Saxifragales ranges from 147,048 bp (Phedimus kamtschaticus) to 160,410 bp (Liquidambar formosana), as shown in Table S1, and the total of G + C content varies from 36.40% (Myriophyllum spicatum) to 38.55% (Paeonia brownii).

The plastome of S. plumbizincicola contains 133 genes, including 85 protein-coding genes (PCGs), 36 tRNA genes, 8 rRNA genes, and four pseudogenes (one ycf1, one rps19, and two ycf15). Dong et al. [8] reported that infA and rpl32 have been lost from three species of Paeonia plastome (Paeonia brownii, P. suffruticosa, and P. obovata). In this study, comparative analysis showed that these two gene loss events occurred in all eleven plastomes of Paeoniaceae. A possible explanation is that the two functional genes have been transferred to the nucleus [8,47,54,55,56,57,58,59,60]. Furthermore, Dong et al. [8] observed that the intron of rpl2 was completely lost in Saxifraga stolonifera. There are currently about 640 species in 33 genera recognized within the family Saxifragaceae [61]. Interestingly, in the current study, the intron of rpl2 was detected in all families in Saxifragales, except for 13 species from the examined 8 genera representing the major lineages of Saxifragaceae, as shown in Table S1, which indicates an early loss of this intron within this lineage. Besides Saxifragaceae, nine other independent losses of rpl2 intron were reported in dicotyledons [62,63,64,65,66,67,68,69,70,71,72,73,74,75]. The two most probable mechanisms of loss of the rpl2 intron are homologous recombination and gene conversion [64,76,77].

3.2. Structure Analyses of Plastid Ribosomal RNA Operon

3.2.1. Structure of 16S rRNA

Similar to most other plants, the size of S. plumbizincicola rrn16 is 1490 bp. In all Saxifragales species examined, the sizes of rrn16s are the same as that of S. plumbizincicola, except for the family Paeoniaceae, with an insertion (U) between positions 576 and 577 nts. As shown in Table S2, the G + C content of the rrn16s of Saxifragales ranges from 56.5% (Rhodiola rosea) to 56.9% (Fortunearia sinensis, and Sinowilsonia henryi). The average G + C content for typical land plants is 56%, whereas this value falls from 52% to 28% for holoparasitic angiosperms, with an increasingly greater number of mutations [78].

We next examined the predicted secondary structure of 16S rRNA in S. plumbizincicola. The structure is similar to the models proposed for other plants [78,79,80], including three main domains organized in 74 helices. In total, 72 mismatched pairs have been detected, and most of them (58/72) are G-U wobble pairs, as shown in Figure 2. Furthermore, we also detected that the position 123 nt of 16S rRNA is cytosine (123-C), whereas other Saxifragales species examined are uracil. To avoid a potential sequencing error, we confirmed the mutation U123C by transcriptomic data of S. plumbizincicola (accession number: SRR5118122-SRR5118124). For further analysis, the 16S rRNAs from 3125 reference plastomes of land plants deposited in GenBank were investigated. The survey results indicated that only 13 species had the special 123-C, including two hyperaccumulator plants, Alpinia oxyphylla and Curcuma longa [81,82]. In contrast with non-canonical base pairing (G-U), we particularly observed that the mutation U123C of 16S rRNA can form stabilized base pairing (C-G) in helices H120, as shown in Figure 2. However, the underlying biological mechanisms of the mutation U123C of 16S rRNA are still unknown.

Figure 2

The predicted secondary structure model of 16S rRNA of S. plumbizincicola. Roman numbers refer to domain names.

3.2.2. Structure of 23S rRNA and 4.5S rRNA

As can be seen from Table S2, the size of rrn23 spans from 2089 bp (Sedum) to 2857 bp (Paeonia suffruticosa), and the G + C content ranges from 55.0 (Corylopsis coreana, Loropetalum subcordatum, and Chrysosplenium aureobracteatum) to 55.4% (M. spicatum), with an average value of 55.1%. In contrast to rrn23, the rrn4.5 of Saxifragales is remarkably conserved in size (103 bp), with a mean G + C content of 56.7%. The rrn4.5 and rrn23 genes are separated by 98–99 bp intergenic spacers (IGd), with G + C content between 57.1% and 60.2%, as shown in Table S2.

The predicted secondary structure of 23S rRNA in S. plumbizincicola is similar to the models of Gutell [80,83], containing 149 helices and six domains, as shown in Figure 3. Moreover, a total of 135 mismatched pairs with 101 G-U wobble pairs were found in the structure. We then comparatively analyzed 23S rRNA secondary structures of all investigated taxa in Saxifragales. Remarkably, as shown in Figure 4, the hairpin loops near helix H550 were more divergent than others, including nucleotide substitutions and indels. In particular, these divergent hairpin loops may have potential phylogenetic implications. For instance, all species of Crassulaceae are characterized by six nucleotides (5’-CACUGG-3’) in these hairpin loops. In addition, in contrast to S. plumbizincicola, P. suffruticosa had an extra 46 nts insertion between the helices H1684 and H2037 of 23S rRNA. Our study further shows that the extra insertion may form two additional helices, as shown in Figure 5. Notably, 4.5S rRNA is a unique component of plastid ribosomes from nonvascular (bryophytes) to vascular plants (pteridophytes, gymnosperms, and angiosperms), which is located on the large subunit. Several previous studies of 4.5S rRNA have failed to find known homologues in other types of ribosomes [84,85,86]. In ongoing follow-up research, 4.5S rRNA has been identified as structurally homologous to the 3’ terminus of bacterial, cyanobacterial, and green algal 23S rRNA [19,84,87,88,89,90]. Based on sequence identity analysis, 4.5S rRNA of S. plumbizincicola and 3’ terminus of Escherichia coli 23S rRNA (accession number: J01695) share 62.9% nucleotide identity. Interestingly, despite a considerable amount of nucleotide substitutions and indels between these two regions, their secondary structures exhibited similar topology, as shown in Figure 6. This finding confirms once again that 4.5S rRNA of higher plants is associated with fragmentation of 23S rRNA progenitor.

Figure 3

The predicted secondary structure model of 4.5S rRNA and 23S rRNA of S. plumbizincicola. Roman numbers refer to domain names. The bases of 4.5S rRNA are presented with blue and the bases of 23S rRNA are presented with black.

Figure 4

The predicted secondary structure models of H550 of 23S rRNA among Saxifragales species.

Figure 5

The predicted secondary structure models of the extra insertion between the helices H1684 and H2037 within Paeonia suffruticosa 23S rRNA compared with S. plumbizincicola. The color red indicates that these bases belong to the extra insertion of P. suffruticosa.

Figure 6

The secondary structure models of the 4.5S rRNA of S. plumbizincicola and 3’ end of Escherichia coli 23S rRNA. The variations between structures are presented with red. The color blue indicates the bases belonging to 23S rRNA of S. plumbizincicola.

3.2.3. Structure of 5S rRNA and Evolutionary Implications of Its Putative Promoter

Structurally, 5S rRNA is the smallest RNA component of the large ribosomal subunit in all known organisms [91]. In the S. plumbizincicola plastome, rrn5 and rrn4.5 are physically linked by the intergenic region (IGe), with the size 219 bp, as shown in Table S2. Besides, the predicted secondary structure of S. plumbizincicola 5S rRNA is similar to that of other published studies [92,93], harboring five helices, as shown in Figure 7. Furthermore, our comparative sequence analysis identified a perfectly conserved 121-bp rrn5 among Saxifragales, with medium G + C content (about 52%), as shown in Table S2. In this study, we also used the 5SRNAdb (http://combio.pl/rrna/) to survey the G + C content of plastomic rrn5. A total of 839 sequences were downloaded and analyzed. The mean G + C content is 50.73%, with the lowest in Euglena viridis (32.26%) and the highest in Staurastrum punctulatum (59.84%). The survey shows that there is a great variability in G + C content of rrn5 for photosynthetic euglenoid and green algae.

Figure 7

The predicted secondary structure model of 5S rRNA of S. plumbizincicola.

Based on similarity of nucleotide sequences, Audren et al. [94] found that a prokaryotic type promoter, which is closely related to the bacterial consensus, was located upstream of the rrn5 and downstream of the stem-loop structure from spinach. However, the putative promoter is inactive both in vivo and in vitro, likely due to the high GC content of the sextama box (TTGGGG) [94,95]. A number of studies have demonstrated that the 5S rRNA gene is transcribed with the other ribosomal genes within the same operon [19,94,96,97]. Notably, the spinach putative promoter was also detected in the similar region from all 44 Saxifragales species. As shown in Figure 8, it contains a sextama box (−35 region, T100T100G100G100G100G100) and a pribnow box (−10 region, C57A100A100T100A100T86) separated by 8–29 bp within Saxifragales, as shown in Figure 8. Interestingly, we found the sequence of putative rrn5 promoters have some evolutionary implications. For example, all spacers between -35 and -10 boxes from 44 investigated species share the 16 common nucleotides (CCTCACAATCACTAGC), except for Liquidambar formosana (CCTCTAGC). Due to nucleotide insertion, deletion, and substitution, the ancestral sequence was then further evolved to different apomorphies in diversified lineages within Saxifragales.

Figure 8

The putative promoters of rrn5 among Saxifragales species.

3.3. Phylogenetic Implications

To investigate the evolutionary relationships among the order Saxifragales, we performed phylogenetic analyses using 83 plastid genes of 44 species. Two species of Vitaceae (V. heyneana and V. vinifera) were employed as outgroups. After alignment, the concatenated sequences are 74,751 bp long. The trees derived from ML and BI analyses display the same topology, as shown in Figure 9. According to the Angiosperm Phylogeny Group (APG) system IV [98], the order Saxifragales comprises 15 families, 11 of which were chosen for the phylogenetic analyses. The order Saxifragales can be generally divided into two clades: core Saxifragales clade (maximum likelihood bootstrap [BS] = 100 and bayesian posterior probability [PP] = 1.0) and Paeoniaceae plus the woody clade ([BS] = 89 and [PP] = 1.0). The former clade is subdivided into two subclades: one containing Crassulaceae, Haloragaceae, and Penthoraceae, and the other comprising three families of Saxifragaceae alliance (Grossulariaceae, Saxifragaceae, and Iteaceae). The latter clade includes Paeoniaceae, Altingiaceae, Cercidiphyllaceae, Daphniphyllaceae, and Hamamelidaceae. In general, the framework of relationships within Saxifragales generated from this study agrees with those reported by Jian et al. [11], Moore et al. [99], and Soltis et al. [12].

Figure 9

Nucleotide-based phylogenetic tree of 44 Saxifragales species. This analysis is based on 79 protein-coding genes (PCGs) and 4 rRNAs. The maximum likelihood bootstrap (BS) and bayesian posterior probability (PP) values for each node are indicated; * indicates 100% bootstrap or 1.00 PP. The bootstrap percentages < 50% and posterior probabilities < 0.5 were omitted.

In the present study, we found that S. plumbizincicola had a closer relationship with S. sarmentosum than S. oryzifolium. Furthermore, Sedum is sister to (Phedimus + Rhodiola). Species of Phedimus, previously treated as members of Sedum, have been classified as a separate genus [100,101]. Our data support this taxonomic revision of Phedimus.

Within Saxifragaceae alliance, Iteaceae is sister to (Grossulariaceae + Saxifragaceae), with strongly supported nodes ([BS] = 100 and [PP] = 1.0). Furtherly, Saxifragaceae can be divided into two subclades: heucheroid and saxifragoid [61,101]. Within the heucheroid, two genera, Heuchera and Tiarella, have been suggested as polyphyletic by several chloroplast markers [102,103]. Our present study based on nearly whole plastome sequence data supported this view. However, both morphology and nuclear internal transcribed spacers (ITS) data have indicated that the two genera are monophyletic [102,104,105]. This incongruence between chloroplast and nuclear gene trees may be due to chloroplast capture [105,106,107,108,109,110].

Our results also accepted the monophyly of the woody clade, which is sister to the family Paeoniaceae. It is noteworthy that deep-level relationships within Hamamelidaceae are strongly supported. Nevertheless, the closest relatives of this family and relationships among these woody families are still unresolved in our analysis. This might partially be attributed to an ancient, rapid radiation [11]. Therefore, further detailed analyses need be conducted to evaluate the relationships within the woody clade.

4. Conclusions

In the present study, we first sequenced and analyzed the plastome of S. plumbizincicola. The genome structure and gene order were revealed, including 85 PCGs, 36 tRNA genes, 8 rRNA genes, and four pseudogenes. Next, we focused on the analyses of the primary and secondary structures of plastid rRNA genes. Notably, we found the sequence of putative rrn5 promoter has some evolutionary implications within the order Saxifragales. Based on the 83 plastid genes from 44 species, phylogenetic analyses demonstrated that S. plumbizincicola had a closer relationship with S. sarmentosum than S. oryzifolium. Our findings reported here shed light on the structural evolution of plastid rRNA operon and phylogenetic relationships within Saxifragales.