Complete chloroplast genome sequence of the mangrove species Kandelia obovata and comparative analyses with related species

Sampling and DNA sequencing

Samples were collected in Xinyin National Wetland Park, Danzhou, China (19°30′N, 109°30′). Voucher specimens were deposited in the herbarium of Hainan Normal University under accession number Yang Y-201803, and replicate specimens (XY2019061201) were sent to Traditional Chinese Medicine Herbarium of Hainan Province for conservation. Five fresh leaves were collected from five healthy trees of K. obovata and then stored on ice for the return to the laboratory. Total genomic DNA was extracted from mixed fresh leaf tissues using the CTAB method (Doyle & Doyle, 1987). The extracted total genomic DNA was dissolved in 50 µL of TE buffer. After quality and concentration were analyzed by agarose electrophoresis and spectrophotometry (Beijin Puxi T6, China), a final DNA concentration of >30 ng/mL was used for Illumina sequencing. Library preparation and sequencing were performed at TGS-Shenzhen, China. The genome was sequenced on an Illumina HiSeq X Ten platform (Illumina, San Diego, CA, United States) with 150 bp paired-end reads.

Chloroplast genome assembly and annotation

To ensure accurate and reliable analyses, raw data were proofread and assembled with NOVOPlast2.7.1 (Dierckxsens, Mardulyn & Smits, 2016) and verified by SOAPdenovo2 (Luo et al., 2012). A partial rbcL gene sequence of Ceriops tagal (MH240830) in Rhizophoraceae (Chen et al., 2019) was used as seed, and the chloroplast genome sequence of C. tagal was used as the reference sequence. The chloroplast genome of K. obovata was annotated by the Dual Organellar GenoMe Annotator program (DOGMA; http://dogma.ccbb.utexas.edu/) (Wyman, Jansen & Boore, 2004). The initial annotations and putative start, stop, and intron positions were checked manually based on comparison with the chloroplast genomes of Erythroxylum novogranatense (NC030601), as Erythroxylaceae is sister to Rhizophoraceae (Setoguchi, Kosuge & Tobe, 1999). Additionally, tRNA genes were identified by the tRNAscan-SE 1.21 program (Schattner, Brooks & Lowe, 2005). Physical maps were drawn using the web tool Organellar Genome DRAW (OGDRAW) v1.2 (Lohse, Drechsel & Bock, 2007).

Simple sequence repeat (SSR) analysis

Distributed throughout the genome, SSRs are repeat sequences with a typical length of 1–6 bp that are generally considered to have a higher mutation rate than neutral DNA regions. The distributions of SSRs in the chloroplast genome were predicted by using the microsatellite search tool MISA (Kurtz et al., 2001) with the following parameters: ≥10 for mononucleotide repeats, ≥5 for dinucleotide repeats, ≥4 for trinucleotide repeats, and ≥3 for tetranucleotide repeats, pentanucleotide repeats, and hexanucleotide repeats.

Comparative Genome analysis

To investigate the sequence divergence of the chloroplast genome among the analyzed mangrove species, the whole chloroplast genome sequences of the mangrove species K. obovata and C. tagal of Rhizophoraceae, L. littorea of Combretaceae, and S. alba of Lythraceae; the semi-mangrove species B. racemosa of Lecythidaceae; and the land species E. novogranatense of Erythroxylaceae were analyzed using the mVISTA program in the Shuffle-LAGAN mode (Frazer et al., 2004). The K. obovata annotations were used as references. The differences in the chloroplast genome length, LSC length, SSC length, GC content, encoding gene types and gene numbers among these 6 species were analyzed. The LSC/IR/SSC boundaries among the species were determined by comparative analysis to explore the variation in these angiosperm chloroplast genomes.

Phylogenetic analysis

To understand the phylogenetic position of Rhizophoraceae in Malpighiales, we selected the chloroplast genome sequences of 22 other species published in the NCBI Organelle Genome Resource database and used 54 chloroplast protein-coding genes atpA, atpB, atpE, atpF, atpH, atpI, ccsA, clpP, matK, ndhA, ndhE, ndhG, ndhH, ndhI, ndhJ, ndhK, petA, petD, petG, petL, petN, psaA, psaB, psaC, psaI, psaJ, psbA, psbC, psbD, psbF, psbH, psbJ, psbL, psbM, psbN, psbT, rbcL, rpl14, rpl16, rpl33, rpl36, rpoA, rpoB, rpoC1, rps11, rps14, rps15, rps18, rps2, rps3, rps4, rps8, ycf3, ycf4) to construct a phylogenetic tree; Linum usitatissimum (NC_036356) of Linaceae were used as the outgroup. After downloading the chloroplast genome sequences from the NCBI database, the common 54 protein-coding genes among the 22 complete chloroplast genomes were aligned by MUSCLE (Edgar, 2004). Phylogenetic trees were constructed using maximum likelihood (ML) and maximum parsimony (MP) methods with MEGA7 software (Kumar, Stecher & Tamura, 2016) (https://www.megasoftware.net/), and the GTR+I+G model was used for the ML analysis. One thousand bootstrap replicates were used to calculate the local bootstrap probability of each branch.

Genome organization and gene features

In total, we obtained 5,542 Mb short sequence data with a Q20 of 98.62% using the Illumina HiSeq X Ten platform. We obtained the same chloroplast genome sequence of K. obovata by the two methods (NOVOPlasty and SOAPdenovo2). The K. obovata chloroplast genome is a typical double-stranded circular DNA molecule with a quadripartite structure. The length is 160,325 bp (Fig. 1), with a pair of IR regions 26,670 bp in length that separate an LSC region of 91,156 bp and an SSC region of 15,829 bp. The GC content of the chloroplast genome is 35.23% (Table 1), and the GC contents of the LSC, SSC and IR regions are 32.34%, 29.05% and 41.99% respectively. As the rRNA genes rrn23, rrn16, rrn5 and rrn4.5 are located in the IR region, the IR region shows a higher GC content. The above-described content of the K. obovata chloroplast genome is similar to the contents of other chloroplast genomes in Malpighiales (de Santana Lopes et al., 2018; Li et al., 2017). However K. obovata, C. tagal and E. novogranatense, all of Malpighiales, have larger genome sizes and LSC sizes and lower GC contents than those mangrove species from different orders (Table 1).

The chloroplast genome contains a total of 128 unique genes, including 80 protein-coding genes, 38 tRNA genes, 8 rRNA genes and 2 pseudogenes. Among them, 19 of these genes occur in IRs, and contain 8 protein-coding genes (rps19, rpl2, rpl23, ycf2, ndhB, rps7, rps12, ycf1), 7 tRNA genes (trnI-CAU, trnL-CAA, trnV-GAC, trnI-GAU, trnA-UGC, trnR-ACG, trnN-GUU) and 4 rRNA genes (rrn23, rrn16, rrn5, rrn4.5). In total, 15 genes with introns were found. Thirteen of these genes contain one intron, and two of these genes (clp3 and ycf3) contain two introns (Table 2). Rps12 is a trans-spliced gene with a 5′ exon located in an LSC region and two 3′ exons located in IR regions, similar to most other plant chloroplast genomes (Asaf et al., 2016; Cheng et al., 2017; Raman & Park, 2016; Jiang et al., 2017). The start of the nucleotide sequence of a protein-coding gene usually begins with ATG. However, there are some exceptions in the K. obovata chloroplast genome in which the first nucleotide has changed from A to G or the second nucleotide has changed from T to C, such as rps19 and cemA, which begin with GTG, and ndhD which begin with ACG. The K. obovata chloroplast genome with gene annotations was submitted to GenBank under the accession number MH277332.

Simple sequence repeats (SSR) analysis

With the development of next generation sequencing (NGS) technologies, SSR development has become quicker, more efficient and cheaper than before, even in species for which background genetic information is lacking (Davey et al., 2011; Zalapa et al., 2012). In our research, using the microsatellite identification tool MISA, we identified 108 SSR loci in the K. obovata complete chloroplast genome sequence, including 92 mononucleotide SSR loci (A/T), 15 dinucleotide SSR loci (AT/TA), and 1 trinucleotide SSR loci (AAT) (Table 3). A total of 103 of the 108 SSRs are located in intergenic regions, and 5 SSRs are located in gene-coding regions. Furthermore, 83, 6 and 19 SSRs, were discovered within the LSC, SSC and IR regions, respectively (Fig. 2).

Comparative analysis of the chloroplast genome sequences of six species

In this study, the chloroplast genomes of several mangrove species and one land species were analyzed by the mVISTA program. Considerable similarities in genome composition and size were identified among the species (Fig. 3). The coding regions of the two mangrove species in Rhizophoraceae were almost identical, whereas the non-coding regions are more variable. The mangrove species in Malpighiales shows a closer relationship with the land plant E. novogranatense, which belongs to the same order, than to the mangroves species from different orders.

IR contraction and expansion

In this study, we aligned the positions of the LSC, IRA, SSC and IRB borders and the adjacent genes among members of mangrove species and found that the studied locations are generally similar to those of all previously reported chloroplast genomes (Jo et al., 2016; Wei et al., 2017; Zhou et al., 2017; Yu et al., 2018). However, in the analyzed species, two copies of the K. obovata rpl22 gene are located in the boundaries of the LSC/IRA junction and IRB region; as the rpl22 gene in the IRB region has no open reading frame that encodes a functional protein, we regarded this gene as a pseudogene. However, in the other species, rpl22 is only located in the LSC region, and there is no rpl22 in the IRB region. The IR extended into the ycf1 genes, creating long ycf1 pseudogenes with variable lengths. The length of the ycf 1 pseudogene is in the IRA region is 1,353 bp, whereas the ycf1 in the SSC/IRB junction contains intact open reading frames (ORFs). There are two copies of rps19 in the IR region of the K. obovata chloroplast genome, but in the other species, there is only one copy, located in the junction of the LSC/IRA region. In all the analyzed genomes, the trnH gene is the first gene in the LSC region, although its distance from the IRB/LSC junction ranges from 1 bp to 79 bp. Comparisons among the species revealed that the K. obovata chloroplast genome has lost the ndhF gene (Fig. 4).

Phylogenetic analysis

To analyze the phylogenetic position of K. obovata in Malpighiales, we used the common 54 protein-coding genes among the 22 complete chloroplast genome sequences to infer phylogenetic relationships. The phylogenetic trees generated by the ML and MP methods have similar topologies (Fig. 5). Both show that K. obovata is most closely related to C. tagal, with 100% bootstrap supports; both of these taxa belong to Rhizophoraceae. These two species are more closely related to the land species E. novogranatense (Erythroxylaceae) than to the other species, suggesting that Erythroxylaceae is sister to Rhizophoraceae. Our study will provide valuable genetic information for genome-scale phylogenetic studies in mangrove plants.