Introduction

Mangroves are woody shrub and tropical plants that grow well in the inter-tidal zones of tropical to sub-tropical latitudes1. Globally, mangrove covers about 200,000 km2 area2,3 and belongs to family Acanthaceae that comprises of more than 400 species of polyphyletic group of trees4. These trees possess unique physio-morphological adaptations and tolerance against hypersaline environment, tidal cycles, and soil chemistry5,6,7. Importantly, mangroves serve as a hub of exponential ecological resource for habituating diverse marine life by providing a protective sanctuary to bread, shelter and grow living organisms and as a sink for continued carbon emissions. However, overexploitation of mangroves for wood and environmental pollution have drastically affected the tree population as well resulted in loss of genetic diversity8. The evolutionary history, tolerance-based mechanisms, genetic divergence and sub-speciation are the key aspects that require in-depth studies.

Among mangroves, Avicennia marina (Forssk., Vierh., gray mangroves) is one the keystone species of the genus and well-distributed species across different latitudes mostly through the dispersal of diasporas by sea and wind9,10. Gray mangroves have been divided into categories based on their habitat and importance to the community structure11. Currently, at least three sub specific or allopatric varieties of A. marina viz. var. australasica, euclayptifolia and marina have been categorized to date12. The true mangroves are further distinguished as major and minor mangroves where variations also exist in habitat exclusively in aquatic and terrestrial or both13. However, plant scientists did not reach a consensus classification14,15. DNA based molecular markers system are available to discriminate among population and species and perform phylogenetic analysis16,17,18, however, current advancement of next-generation sequencing methods can help to understand genome architecture, structure, content, divergence and evolutionary history that could solve many key questions related to population structure and taxonomy19,20,21,22,23.

In this case, the chloroplast genome offers a highly conserved sequence due to uni-parent inheritance, haploid, and non-recombinant nature24,25. To date, about 4000 chloroplast genomes have been sequenced26. The chloroplast genome size for angiosperm is around 110-165 kb containing unigenes from 90 to 11027,28. Almost all angiosperm chloroplast genomes consist of 4 regions viz. large single-copy region (LSC), followed by inverted repeats (two complementary to each other) and a single small copy region (SSC)29. Due to these reasons, chloroplast genomes sequences are most recommended and used for evolutionary studies, barcoding and phylogenetic analysis30,31. In case of A. marina, though recently Frilis et al.32 reported the first whole-genome data with genome size of 456 million base pair, however, still the sequence datasets for Acanthaceae (± 400 species) is only 15 organelle genomes. The genus suffers from lack of a comprehensive comparative genomic assessment to understand the phylogenomic and evolutionary history of A. marina and related species. Hence, in the current study we aimed to sequence A. marina growing in Oman and perform comparative chloroplast genome analysis with Andrographis paniculate, Aphelandra knappiae, Clinacanthus nutans, Echinacanthus attenuates, Echinacanthus lofouensis, Echinacanthus longipes, Echinacanthus longzhouensis, Justicia flava, Justicia leptostachya and Strobilanthes cusiaspecies from Acanthaceae family.

Results

Chloroplast genome structure of A. marina

The length of complete chloroplast (cp) genome of A. marina is 150,279 bp and it exhibits a typical quadripartite structure with a pair of inverted repeats 25,117 bp that separate a large single-copy region 82,522 bp and a small single copy region 17,523 bp (Fig. 1; Table 1). The cp genome contains overall GC contents of 38.6%. We identified 22 introns containing genes (14 protein coding genes and 7 tRNA genes), among these genes 19 genes are single intron and two genes (ycf3 and clpP) contain two introns. The largest intron was found in the gene trnK-UUU with 2,414 bp in size and the smallest intron with 486 bp in gene trnL-UAA. Additionally, the largest exon was found in the gene ndhB with 775 bp and the smallest exon was found in petB genes and was only 6 bp in size. Among the coding genes, rps12 was unequally divided, with its 5′ exon being located in the LSC region and one copy of the 3′ exon and intron being located in each of the IR regions, as in other angiosperms33,34.

Figure 1
figure 1

Gene map of the Avicennia marina chloroplast genome. Genes drawn inside the circle are transcribed clockwise, and those outside the circle are transcribed counterclockwise. The stars indicate the intron-containing genes. Genes belonging to different functional groups are colour-coded. The darker grey in the inner circle corresponds to GC content, and the lighter grey corresponds to AT content.

Table 1 Genomic details and statistics of 11 mangrove species of the family Acanthaceae.

The protein coding genes included larger and smaller subunit proteins such as rpl2, 14, 16, 20, 22, 23, 32, 33, 36 and rps2, 3, 4, 7, 8, 11, 12, 14,15 ,16, 18, 19 respectively (Table 2). Numerous groups of photosynthetic genes were also annotated including photosynthesis I, photosynthesis II, Cytochrome b6/f complex, ATP synthase and rubisco as shown in Table 2. The genome also contains a total of 132 genes including 87 protein coding genes, 8 ribosomal rRNA genes and 37 transfer tRNA genes. The total size for the protein coding region was found 79,321 bp (52% of the cp genome), for the ribosomal RNA (rRNA) the size is 9,054 bp (6.02%) and for the transfer RNA (tRNA) it is 2,839 bp (1.88). GC contents in the coding region of the chloroplast genome was found as 43.4%. Similarly, 36.8, 32.7, 43.6, 53.1, 55.3 and 38.4% of the GC contents were found in LSC, SSC, IR, tRNA, rRNA and protein coding genes, respectively. The GC content in the protein coding genes at first, the second and third position of codon were 46.06, 37.87 and 32.25% respectively as shown in Table 3. The AT distribution in chloroplast genome of A. marina was found 36.7%. The AT in LSC, SSC, IR, tRNA, rRNA and protein coding genes were 63.2, 67.3, 56.4, 46.9, 44.8 and 61.6% respectively. Within the protein coding genes, the AT/U was found 54.8, 62.1 and 67.8% on first, second and third position of codons as shown in Table 3. In A. marina the most common amino acid was found as leucine (10.7%) and the least common amino acid was cysteine (1.2%).

Table 2 Genes annotated in the chloroplast genome of A. marina.
Table 3 Base composition of the A. marina and other related species cp genome.

Comparative analysis of A. marina cp genome with related species

In order to further analyze the characteristics of the A. marina chloroplast genome, its assembled genome was compared with the chloroplast genomes of 10 other species of the same family. The results revealed that A. marina (150,279 bp) cp genome size is slightly bigger as compared to A. paniculata (150,249 bp), J. leptostachya (149,227 bp) and S. cusia (144,133 bp). However, it has a slightly smaller cp genome size when compared to the species A. knappiae, C. nutans, E. attenuates, E. lofouensis, E. longipes, E. longzhouensis and J. flava (Table 1).

The GC contents of A. marina (38.6%) were found almost similar for species (A. knappiae, E. longipes, E. lofouensis and E. longzhouensis), however GC contents were found higher when compared to the species A. paniculata (38.3%), C. nutans (38.4%), E. attenuates (38.3%), J. flava (38.2%), J. leptostachya (38.2%) and S. cusia (38.2%). The number of genes, rRNA, tRNA, SSC size and the number of protein coding regions were found almost similar in these studied cp genomes Table 1. The LSC analysis showed that it is almost similar in all the species except E. longzhouensis where the size was found slightly smaller (79,203 bp) but significantly higher in S. cusia (92,666 bp) when compared to A. marina (82,522 bp).

The synteny of A. marina cp genome with ten other species from Acanthaceae was analyzed by mVISTA. Gene divergence was determined by determining pairwise alignment of A. marina with related species. A. marina chloroplast genome was used for reference to determine variation and sequence identity in the chloroplast genomes of related species. The results showed high sequence similarities among the cp genomes of several species, especially in protein-coding and IR regions (Figure S1).

For example, the intergenic regions between psbA-matK, rps16-psbL, atpA-atpF, atpH-atpI, rpoC1-rpoB, psbE-petG, petN-psbM, psbD-rps14, ycf3—rps4, atpE-rbcL, accD-ycf4, psbL-petL, clpP-psbN, petD-rpl36, rpoA-rps11, rpl22-rps8, rpl16-rps3, ycf15-ndhB, rps19-rpl23, ndhF-ccsA, ndhD-psaC, ndhI-ndhG, ndhA-ycf1 and ndhB-ycf2 were found highly divergent. Besides these intergenic regions some divergence was observed in protein-coding genes (rps16, rpoC1, rpoC2, atpA, matK, atpF, clpP, rps12, psbN, psbB, psbT, psbH, petB, accD) in LSC region, (petD, rpl16, rpl22, rpl36, rps3, rps8ycf2) in IR region and (ycf1, ycf15,ndhA, ndhB, ndhF and ndhH) in SSC region. The gene divergence was found more prominent in LSC and SSC region as compared to the IR region. In pairwise sequence divergence analysis, A. marina exhibited highest divergence (0.084) with J. flava and showed lowest divergence with A. knappiae (0.059) (Supplementary Data 1). The most divergent genes were rpl22 (J. leptostachya (0.149), rps15 in E. attenuates (0.135), ndhF A. paniculate (0.139), E. lofouensis (0.212), and matK J. leptostachya (0.139) (Fig. 2). Similarly, lowest pairwise divergence was found in genes such as, ndhB (0.003) in E. longzhouensis, psbL (0.009) in E. lofouensis, petN (0.011) in E. longipes. psaJ (0.015) in E. attenuates.

Figure 2
figure 2

Heatmap of pairwise gene divergence in the chloroplast genome sequence of A. marina. The highly divergent genes are represented with light and green color depending on the divergence frequency.

Simple sequence repeat (SSR) analysis of A. marina cp genome

We determined SSRs in the cp genome of A. marina as well as in the cp genomes of the other ten related species (Fig. 3A). A total of 126 SSRs were found in the chloroplast genome of A. marina. Among the predominant SSRs nucleotides such as mononucleotide, dinucleotide, trinucleotide, etc., the mononucleotide SSR was found the most abundant not only in A. marina but also in the other species of the current study. In this report the presence, distribution and types of SSRs were studied (Fig. 3B). The SSRs in the IR, LSC and SSC regions of the chloroplast genomes of above-mentioned species have been analyzed and were found distributed such as mono, di and trinucleotidies. However, most of the SSRs in these species were found mononucleotides in the IR region. For example, mononucleotides SSRs in A. marina, A. paniculata and S. Cusia were 4, in the A. knappiae and C. nutans, the number of SSRs were found 6, in E. lofouensis, E. longipes, J. flava and J. leptostachya have 8 SSRs as mononucleotides, in the E. attenuates and E. longzhouensis 9 and 7 mononucleotide SSRs were found. Similarly, J. flava contain 2 dinucleotide SSRs in the inverted repeat regions, in the IR region the J. leptostachya contain 1 penta and hexa nucleotide SSR each. SSRs in LSC region were found mostly as mononucleotide and trinucleotide in all the above-mentioned species. For example, In the A. marina 90 mononucleotides SSRs were found. In all other 10 species, in LSC, mononucleotide SSRs were found in the range of 56 to 77. For example, the E. longipes contain 2nd higher number (77 mononucleotide SSRs) after A. marina in LSC region. Similarly, the E. longzhouensis has 76 mononucleotide SSRs, A. paniculata 70, C. nutans 68, E. attenuates and S. Cusia 65, A. knappiae and J. flava 59 and J. leptostachya 56 SSRs in LSC region. In some species dinucleotide SSRs were found in LSC region such as 2 in A. marina, 1 in J. flava, 5 in J. leptostachya and 4 in S. cusia. Similarly, in LSC region A. marina contain 2 trinucleotides SSRs while the highest trinucleotides SSRs were 8 in the S. Cusia.

Figure 3
figure 3

Simple sequence repeats (SSRs) analysis in the chloroplast genome of A. marina. (A) Numbers of SSR types in complete chloroplast genome, (B) Number of SSRs in LSC, IR and SSC regions. (C) frequency of SSR motifs.

SSRs were also found in SSC region most prominently as mononucleotides and trinucleotides. For example, in A. marina 24 mononucleotides SSRs were found, A. paniculata 26 mononucleotides SSRs, E. longzhouensis 25 mononucleotides SSRs, E. attenuates 24 mononucleotides SSRs, J. flava 22 mononucleotides SSRs, E. lofouensis, E. longipes J. leptostachya and S. Cusia were fund to have 21 mononucleotides SSRs each. In A. marina no trinucleotide SSRs were found, however only few species were found to have 1 or 2 trinucleotides SSRs (Fig. 3C).

Repeat sequence analysis

In the chloroplast genome of A. marina and related species forward, reverse, palindromic and tandem repeats were analyzed. Total of 25 forward repeats, 14 palindromic repeats and 15 tandem repeats were found in A. marina cp genome. Total of 16, 23, 37, 20, 17, 18, 19, 24, 21 and 16 forward repeats were found in A. paniculata, A. knappiae, C. nutans, E. attenuates, E. lofouensis, E. longipes, E. longzhouensis, J. flava, J. leptostachya and S. cusia respectively (Fig. 4). Furthermore, tandem repeats in these species were 15, 14, 62, 19, 14, 15, 14, 19, 21 and 21 respectively. Among these repeats in A. marina, 16 of the forward repeats were found in the 15–29 bp, 5 repeats in 30–44 bp, 1 repeat was found in 45–59 bp, 75–79 bp and 2 repeats were found in > 90. Similarly, in the tandem repeats, 9 repeats were found in 15–29 bp, 2 repeats were found in 30–44 bp, 1 repeat was found in 60–74 bp, and 75–79 bp and 2 repeats in > 90 bp. Among the palindromic repeats, 5 repeats were found in the 15–29 bp, 6 repeats were found in 30–44 bp, and 1 repeat was found in 45–59 bp, 75–79 bp and > 90 bp.

Figure 4
figure 4

Distribution and frequency analysis of repeated sequences in chloroplast genome of A. marina. (A) total number of forward, palindromic and tandem repeats. (B) distribution of palindromic repeats in different ranges of length. (C) distribution of forward repeats in different ranges of length. (D) distribution of tandem repeats in different ranges of length.

The highest number of forward repeats (37) were found in the C. nutans and the lowest number of forward repeats (16) were found in the A. paniculate and S. cusia. Similarly, the highest tandem repeats (62) and the lowest number of tandem repeats (14) were found in the A. knappiae, E. lofouensis and E. longzhouensis. Additionally, among all the species, the highest number of palindromic repeats (25) were found in A. paniculate and S. cusia and the lowest palindromic repeats (14) were found in the A. marina (Fig. 4A–D).

Contraction and expansion of IR region

In the current study, contraction, and expansion in the 4-junctions (JSA, and JSB, JLA, JLB) between IRa and IRb, LSC and SSC regions of A. marina and 10 related species from the family Acanthaceae were analyzed. Size of IR regions in all the species were found almost same (25,000 bp to 25,700 bp) except in the specie S. cusia in which the IR region was found significantly lower (16,328 bp) as compared to the other species (Table 1). Some of the genes present on the borders of IR, SSC and LSC region were duplicated. For example, at JLB border, the rps19 gene was found in the LSC region in the species A. marina, A. paniculata, A. knappiae and C. nutans (Fig. 5). However, this rps19 gene was found in the IRb region about 15 bp away form JLB border in E. attenuatus, E. lofouensis and E. longipes, 19 bp in E. longipes, 97 bp in E. longzhouensis, 102 bp in J. flava and 85 bp in the J. leptostachya. Similarly, the gene ndhF was found partly extended in the IRb regions in some species such as E. lofouensis 70 bp, E. longipes 93 bp, E. longzhouensis 70 bp, J. flava 117 bp, J. leptostachya 122 bp and S. cusia 44 bp extended in the IRb region at JLB junction. However, the same gene (ndhF) was found in the SSC region in A. marina about 78 bp away for JSA border. On the other hand, this ndhF gene is partly extended to IRa region in some species where it is located about 40 bp, 41 bp 43 bp and 109 bp in the IRa regions in (A. paniculate, A. knappiae, C. nutans and E. attenuatus). Additionally, the ycf1 gene was found at the IRs and SSC border in E. lofouensis, E. longipes, E. longzhouensis, J. flava, J. leptostachya and S. cusia were it extended to IRa region with (797, 820, 797, 812, 817 and 771) bp respectively (Fig. 5).

Figure 5
figure 5

Distance between the adjacent genes and junctions of Large single Copy (LSC), two inverted repeats (IR) and Small Single Copy (SSC) in A. marina and 10 other related species. boxes below and above the lines represent the genes located on the border of LSC, IR and SSC. The figure represents the changes in the length and location of certain genes across these borders.

Phylogenetic relationships among A. marina and related species

To determine the phylogenetic position of the A. marina, complete chloroplast genome was performed and compared with 24 other related species to build the phylogenetic tree. The phylogenetic position of A. marina is established using maximum likelihood (ML), maximum persimony (MP) and neighbour-jouining (NJ) methods in this study by utilizing 65 shared genes and complete cp genomes of related plant species. The phylogenetic trees constructed on three different methods shows same result and A. marina formed a clade near to A. paniculata and A. knappiae genomes from same family Acanthaceae. Similarly, Gesneriaceae was found the closest family to Acanthaceae (Fig. 6; Figure S2).

Figure 6
figure 6

Phylogenetic trees were constructed for twenty-four species from four families using three different methods maximum likelihood (ML), maximum parsimony (MP) and neighbor -joining (NJ) by using 65 shared genes. Numbers above the branches are the bootstrap values of ML, MP and NJ respectively. Red star represents the position for Avicennia marina.

Discussion

A. marina belongs to the family Acanthaceae which consist of more than 400 species distributed around the world35. The phylogenetic analysis is very important to correctly identify the taxonomic position of plants. Recently chloroplast genome sequencing has played a major role in the phylogenetic analysis of plants36,37. Though the chloroplast genomes are highly conserved when the size and genomic architecture are considered, however, the genes located on the borders of IR/SC varies tremendously in terms of size and type which make chloroplast distinctive for phylogenetic analysis38,39. The chloroplast genome sequence of all the 10 related species were also found conserved, however when compared with each other and with the chloroplast genome of A. marina, differences were observed for the genes located on the borders of IR/SC regions.

The size of complete chloroplast genome of A. marina (150,279 bp) is greater than the previously reported J. laptostachya (149,227 bp)40 and S. Cusia (144,133 bp)41, and almost same size of A. paniculata (150,249 bp)42 and J. flava (150,888)43. However, the size of A. marina chloroplast genome that we sequenced was found lower than the previously reported chloroplast genomes of the related species such as A. knappiae44, C. nutans45, E. attenuates , E. lofouensis , E. longipes46 and E. longzhouensis43. Our study also confirmed the LSC, SSC and IR region that were almost similar in size to the previously reported chloroplast genomes of Eucalyptus globulus47, Coffea arabica L.48 and Camellia japonica L.49. In our study, we have found that the GC contents (38.6%) of A. marina and on the 3rd position of codons, the GC contents (32.25%) were found lower as compared to AT/U contents (67.8%) which support previous reports of chloroplast genome sequences of C. gileadensis (37.9%) and C. foliacea (37.8%)50. The intron containing genes (14 genes) were found in which only 2 genes were found to contain two introns. The introns are very important in the gene expression regulations studies and it has been observed that when present on specific sites/positions, it can positively regulate exogeneous gene expression51. Thus, introns can be valuable tools in order to improve the transformation efficiencies. The intron sequences of chloroplast DNA also has key role in the phylogenetic analysis52.

Similarly, the size of protein coding genes (79,321 bp), rRNA (9,054 bp), tRNA (2,839 bp) and the number of genes (132), protein coding genes (87), rRNA (8) and tRNA (37) were found in our study that coincides with the previously reported studies on chloroplast genomes50. For example, in the same family (Acanthaceae), Ding, et al.42 and Huang, et al.44 previously reported the size of protein coding genes, rRNA, tRNA, number of genes, number of protein coding genes, number of rRNA and number of tRNA in A. paniculate and A. knappiae which is similar to what we found in the A. marina. Similarly, almost same results were observed in the recently reported chloroplast genomes of C. gileadensis50, Teucrium species53, Vachellia nilotica and Senegalia senegal54.

In the chloroplast genome, microsatellites or SSRs distributed in the genome with the length of almost 1 bp to 6 bp. Previously many studies reported SSRs at different positions in chloroplast genomes55,56. In the present study, SSRs distributed at different locations such as IR regions, LSC and SSC region were identified in the A. marina. Total of 126 SSRs were identified in A. marina and almost similar number of SSRs were found in the chloroplast genomes of related species such as E. Longipes (122) and E. longzhouensis (120). Other studies on chloroplast genomes have also confirmed the presence of uneven numbers of SSRs in the chloroplast genomes at different locations50,57,58,59.

Additionally, the forward, tandem and palindromic repeats in the A. marina were analyzed. The number of forward repeats were found 25, reverse repeats 11, palindromic repeats 14 and 15 tandem repeats were found in A. marina. The number of forward, reverse, tandem and palindromic repeats were found almost similar in the chloroplast genomes of the related species such A. paniculate42 and44. The pairwise alignment of gene in A. marina showed divergence of various genes in the coding regions as well as in the intergenic regions. The most divergent genes include rps16, rpoC1, rpoC2, psbH, petB and accD. Similarly the most prominent divergence in the intergenic region was found in psbA-matK, rps16-psbL, atpA-atpF, atpH-atpI, rpoC1-rpoB, petN-psbM, psbD-rps14, ycf3—rps4, atpE-rbcL, accD-ycf4, psbL-petL, rps12-psbN, petD-rpl36, rpoA-rps11, rpl22-rps8, rpl16-rps3, ycf15-ndhB, rps19-rpl23, ndhF-ccsA, ndhD-psaC, ndhI-ndhG, ndhA-ycf1 and ndhB-ycf2. The gene divergence previously studied52,60,61,62 in the chloroplast genomes strongly confirms our findings of gene divergence.

The difference in the size of chloroplast genome is considered a common evolutionary practice which can be attributed to the contraction and expansion of the inverted repeats (the most conserved regions) in the genome63. In majority of the plants, the junctions or borders of the genome particularly in the quadripartite structure are very conserved. However, inversion at the borders or junction can be found in some species as previously reported64 as well as the loss of genes and the expansion and contraction which is the common cp genome event in the angiosperms65,66. The results in our study shows that rps19 gene at JLB junction in A. marina has extended to IRb (15 to 102 bp) while the gene ndhF partly on the IRa region (40 to 178 bp) shows similar results patterns as previously studied by Cheon, et al.67.

The phylogenetic analysis is very important for the evolutionary and taxonomic studies. Many phylogenetic analysis are now based on the chloroplast genomes68,69. Previously, the phylogenetic analysis of mangrove (A. marina) were mostly based on the RAPD and other molecular markers70,71. Sahu, et al.72 reported that the role of multiple gene in the mangrove phylogeny. However, in this study we used both whole cp genomes and concatenated 65 protein coding genes to infer the phylogenetic position of A. marina. Both data sets showed same results and A. marina is closely related to A. paniculate and A. knappiae in family Acanthaceae. The present study provides a valuable analysis of the complete plastome of A. marina and related species, which may facilitate species identification and both biological and phylogenetic studies.

Conclusion

Complete chloroplast genome sequence of A. marina was found highly conserved in its structure and order of genes distribution as compared to the other mangrove species in Acanthaceae. The results showed location and distribution of SSRs as well as the sequence divergence among the chloroplast genome of A. marina and related mangrove species. Among genes, rpl22, ndhF, rps15 and ndhA were found the most divergent genes in the mangrove species. Additionally, the phylogenetic analysis shows that A. marina was closer to A. knappiae and A. paniculate species. It can be concluded from this study that complete chloroplast genome sequence may provide a better understanding of identification and phylogenetic studies of plant as compared to other strategies.

Materials and methods

Plant sample, Chloroplast DNA extraction and sequencing

A. marina seedlings were donated by the Center for Marine Conservation, Ministry of Environment, Sultanate of Oman. After shipping it to the greenhouse the leaves were collected in liquid nitrogen and ground to a fine powder. The powder samples of leaves were processed for chloroplast DNA extraction. Leaves were collected from the mangrove plants and were ground into fine powder using liquid nitrogen. Chloroplast DNA was extracted using the protocol of Khan et al.73. Manufacturer's instructions (Life Technologies, Carlsbad, CA, USA) were followed for the preparations of genomic libraries. Ion Shear Plus Reagents kit was used to shear chloroplast DNA into 400 bp fragments while Ion Xpress Plus gDNA Fragment Library kit was used for library preparation. Quantity and quality of libraries were checked using the Qubit 3.0 fluorometer and Bioanalyzer (Agilent 2100 Bioanalyzer system; Life Technologies, Carlsbad, CA, USA). The template was amplified using Ion 520 & 530 OT2 Reagents and enriched by using the OneTouch instrument and Ion OneTouch ES enrichment system. The final chloroplast DNA sample was loaded into Ion S5 Sequencing (supported by Mawarid, Oman’s Animal & Plant Genetic Resource Center, Ministry of Higher Education, Researcher & Innovation, Oman) Chip followed by sequencing through the Ion Torrent S5 protocol34,53.

Genome assembly and annotation

A total of 2,256,283, raw reads were generated for the chloroplast genome of A. marina. The cp genome of A. marina was mapped with other chloroplast genomes such as Andrographis paniculate74, using Bowtie2 (v.2.2.3)75 in Geneious Pro (v.10.2.3) software76. For the assembly, the mean coverage of A. marina was 99X. IR region was identified by using the MITObim (v.1.8) software50,77 was used for the sequence length adjustment. FASTQC78,79 was performed after sequencing to check the read quality. For biases reduction in the analysis, an in-house script was used to filter reads if less than 90% of the bases that made up the read were below Q20. Adopter sequences were removed using the Trimmomatic (v0.36)80. Bowtie2 in Geneious Pro (v.10.2.3)76 were used for mapping only high quality reads.

Dual Organellar Genome Annotator (DOGMA)81 was used for the genome annotation. For the position’s identification of transfer and ribosomal RNA and coding genes BLASTX and BLASTN were used. Furthermore, tRNAscan-SE version 1.2182 was used for the tRNA gene annotation. In order to compare the cp genome of A. marina with the related cp genomes, tRNAscan-SE and Geneious software were used. Intron boundaries, stop and start codons were adjusted by comparing it to the previously reported cp genome of the related specie from family Acanthaceae. OGDRAW83 was used for the illustrations of structural features of A. marina cp genome. MEGA-X software84 was used for the determination of deviations in synonymous codon usage and relative synonymous codon usage with amino acid composition influence. The pairwise gene divergence was determined by mVISTA85 in Shuffle-LAGAN mode.

Identification of repetitive sequences and SSRs

The REPuter software86 was used for the identification of forward, reverse and palindromic repeats, the sequence identity with 90% and minimum 15 bp was the basic criteria in the identification. Microsatellite analysis of contig sequences was carried out with the MIcroSAtellite (MISA) identification tool87. The parameters (unit_size, min_repeats) were defined as follows: 1–10, 2–8, 3–4, 4–4, 5–3, and 6–3; the minimum distance between two SSRs was set to 100. For the identification of tandem repeats, Tandem Repeat Finder version 4.08788 b was used for tandem repeat identification.

Sequence divergence and Phylogenetic analysis

The pairwise distance among the cp genes shared by A. marina and related species as well as the complete cp genome of A. marina was determined. The ambiguous and missing genes annotation were identified and removed by the comparative sequence analysis. Complete cp genomes were aligned by using the MAFFT version 7.22289 and Kimura's two-parameter (K2P) model90 was used for calculation of pairwise sequence distance. The divergence of the new A. marina cp genomes from related species of family Acanthaceae was assessed using mVISTA85 in Shuffle-LAGAN mode and by employing the new A. marina genome as reference.

For the determination of phylogenetic position of A. marina, 23 published cp genomes were downloaded from NCBI from four different families (Gesnerianceae, Lamiaceae, Bignoniaceae and Acanthaceae) respectively. Complete cp genomes and a separate partition containing only the 65 shared genes (concatenated) were used to infer the phylogenetic position of A. marina. MAFFT version 7.22289, with default parameters were used for the alignment of both complete genomes and shared genes data sets. Maximum likelihood (ML) and neighbour-joining (NJ) methods were used to infer the phylogenetic trees with MEGA-X91 and parameters were adjusted with a BIONJ tree with 1000 bootstrap replicates using the Kimura 2-parameter model with gamma-distributed rate heterogeneity and invariant sites. Maximum parsimony (MP) by using PAUP92 using previously described settings33,93.