Ancient Horizontal Gene Transfers from Plastome to Mitogenome of a Nonphotosynthetic Orchid, Gastrodia pubilabiata (Epidendroideae, Orchidaceae)

1. Introduction

The advancement in sequencing technology has led to a surge in research on nucleotide sequences of flowering plants and their use for evolutionary and phylogenetic studies. However, these studies have mainly focused on the plastome due to cost and time constraints in sequencing. As of February 2023, 7,588 species of angiosperm plastome sequences have been deposited in NCBI whereas only 347 mitogenome sequences have been deposited in NCBI. Research using plastome sequences has been conducted on various flowering plants, such as Amborella [1], magnolids [2], monocots [3], and rosid [4]. Despite challenges in sequencing, nuclear genomes of some flowering plants have been decoded, including Amborella [5], Oryza [6], Phalaenopsis [7], and Triticum [8].

Both mitochondria and chloroplast are essential organelles in plant function. The plastome is the most widely used genome in plant evolutionary study because it has relatively stable genome structure, size, and rich evolutionary information. In contrast to plastome, the mitogenome vary widely in size and structure across taxa even among different populations of same species [9,10]. The smallest mitogenome is 66 kb for Viscum [11], while the longest mitogenome is 11.3 mb for Silene [12]. It has been suggested that these differences result from genomic rearrangements, repetitive sequences, and/or foreign DNA insertions [13,14,15]. Due to many repeats, long length, and rearrangement of the mitogenome, it is difficult to decipher the overall structure of the mitogenome in majority of plant species. Frequent gene transfers have been reported in the mitogenome of Valeriana of photosynthetic Dipsacales [16]. Intracellular gene transfer (IGT) and horizontal gene transfer (HGT) have been reported in the parasitic plant Aeginetia [17]. Gene transfer has been reported during the symbiotic relationship between fungi and orchids that are closely related to saprophytic life styles of the Orchidaceae [18]. The trace of gene transfer has been studied in a few orchid mitochondria. These imported sequences might be due to past biological contact [10,18,19,20].

The Orchidaceae are one of the largest families of flowering plants, containing approximately 25,000 species [21]. Recent studies have divided this family into five subfamilies: Apostasioideae, Vanilloideae, Orchidoideae, Cypripedioideae, and Epidendroideae [22]. Following the first complete plastome sequencing report of Phalaenopsis aphrodite [23], plastome sequencing reports have been published for all major lineages within Orchidaceae, such as Corallorhiza [24,25], Cymbidium [26], Dendrobium [27], Holcoglossum [28], Neottia [29], and Vanilla [30]. To date, at least 321 orchid plastome sequences have been deposited in NCBI (as of February 2023). The growing number of reported orchid plastome sequences has led to various studies of the Orchidaceae family, focusing on phylogenetics, evolution, and plastome structure [3,30,31,32,33,34,35,36]. These studies often concentrate on specific lineages such as epiphytic or non-photosynthetic orchids. There are a total of 43 genera that contain non-photosynthetic mycoheterotrophs [37]. Several plastome sequences of these types of orchids have been reported, such as Aphyllorchis [29], Corallorhiza [38], Cyrtosia [33], Epipogium [39], Gastrodia [40], Hexalectris [35], Lecanorchis [41] and Rhizantella [31]. Typical photosynthetic orchid plastomes contain 83 protein coding genes [23]. However, non-photosynthetic species such as Epipogium roseum have significantly reduced plastomes, containing only 17 protein coding genes [39]. This contraction of plastomes results in low levels of phylogenetically informative sites. Genes retained in these plastomes often exhibit elevated evolutionary rates [42], making it difficult to construct accurate phylogenetic trees and interpret their results [3,41].

For the Orchidaceae family, the mitogenome sequence has only been reported for Gastrodia elata [40] and nuclear genome sequences have only been reported for Cymbidium [43], Dendrobium [44], Gastrodia [40], and Phalaenopsis [7]. Despite active research into phylogenetic and evolutionary aspects of Orchidaceae, phylogenetic trees from different organelles are often inconsistent [41,45,46]. This might be due to factors such as incomplete lineage sorting, ancient plastome capture, and/or sampling issues [41]. To address these incongruences, more mitochondrial and nuclear genome sequence data are needed.

The species Gastrodia pubilabiata belongs to the tribe Gastrodieae of the subfamily Epidendroideae of the family Orchidaceae [22]. It is widely distributed in Korea, Japan, and Taiwan [47,48], but it is a rare and endangered species. For its closely related species, G. elata, plastome sequences along with nuclear and mitogenome sequences have been reported [40,41]. The plastome of G. elata is greatly reduced compared to that of a typical photosynthetic orchid. Flowers and stems of G. pubilabiata are pale greenish in color. In a previous study on Cymbidium macrorhizon [49] which also has pale greenish stems and flowers, the plastome was almost intact condition, although it is known as a non-photosynthetic species. Only a few ndh genes were missing from its plastome. Therefore, the extreme reduction of plastome sequence in G. elata raises the question of what type of plastome sequence will be present in G. pubilabiata which has pale greenish stems and flowers.

The current state of research on the plastome and mitogenome of non-photosynthetic orchids presents several difficulties, including a lack of parsimonious informative sites in the plastome and incongruence in phylogenetic trees based on different genomes. To address these issues, this study aimed to decode both of the plastome and mitogenome of the non-photosynthetic orchid G. pubilabiata. Through this research, organelle genomes of non-photosynthetic orchids will be better understood and interactions between different organelle genomes will be confirmed. In addition, any potential foreign nucleotide sequences in organelle genomes will be identified. Through this study, if any evidence of externally imported nucleotide sequences is found, the source of sequences will be investigated and the relationship with the taxon will be determined. Results of this study will contribute to the further understanding of the evolutionary modes of mitogenome in non-photosynthetic orchid species.

2. Results

2.1. Organelle characteristics and comparative analysis

The plastome sequence of Gastrodia pubilabiata was a 30,698 bp in length (Table 1 and Figure 1). The GC content of the plastome was 24.9% with a depth coverage of 236.0x. Unlike other reported orchid plastome sequences such as Cypripedium or Habenaria, the plastome of G. pubilabiata was much shorter in length. In addition, it lacked the typical quadripartite structure of angiosperm plastomes, which typically consist of a large single copy (LSC) and a small single copy (SSC) separated by two inverted repeats (IRs). Instead, the plastome of G. pubilabiata, like that of G. elata, lacked IRs. As a non-photosynthetic, mycoheterotrophic orchid without any leaves, G. pubilabiata lost all photosynthetic related genes, including ndh, pet, psa, psb, and rpo as well as several housekeeping genes, atp genes, rrn5, ccsA, cemA, infA, matK, rbcL, ycf3, and ycf4. Finally, the plastome of G. pubilabiata was annotated. It had only 19 protein coding genes, three tRNA genes, and three rRNA genes.

As a result of de novo assembly, 44 mitochondrial contigs of G. pubilabiata were obtained (Table 1). The goal was to obtain a complete, master circle of the mitogenome of G. pubilabiata. However, it was proved to be challenging due to the presence of AT-rich regions, rearrangements, and the recombination. These 44 mitochondrial contigs were considered as the complete mitogenome with a total length of, 867,349 bp (Table 1 and Figure 2). The coverage depth was 150.0x, with the highest coverage being 627.2x (14,604 bp) and the lowest coverage being 40.1x (4,510 bp). The GC content of the mitogenome was 42.1%. With these 44 mitochondrial contigs, a total of, 38 protein coding genes, nine tRNA genes, and three rRNA genes were annotated. Compared to other monocot mitogenomes, three genes (rpl10, sdh3, and sdh4), nad1 intron1, and nad1 intron2 were missing in the mitogenome of G. pubilabiata (Figure 2).

Repeat analysis was conducted to compare the mitogenome and plastome of G. pubilabiata. Three additional organelles were also included (the mitogenome of G. elata, the plastome of G. elata, and the plastome of G. longistyla) in this analysis. The analysis revealed a total of 437 repeats (Supplementary Table S1). Of these, 279 repeats were found in the mitogenome of G. pubilabiata, 79 repeats were shared between the mitogenomes of G. pubilabiata and G. elata, 44 repeats were shared among the three Gastrodia plastomes, 34 repeats were identified within the mitogenome of G. elata, and one repeat was shared between the G. elata plastome and mitogenome (Figure 3 and Supplementary Table S1). Additionally, 89 simple sequences repeats (SSRs) were found within the plastome of G. pubilabiata, including 16 mononucleotide SSRs, five dinucleotide SSRs, four trinucleotide SSRs, 13 tetranucleotide SSRs, and 51 pentanucleotide SSRs (Supplementary Table S2).

2.2. Foreign sequences in the organelle

BLASTN and BLASTP searches were carried out for the plastome of G. pubilabiata to identify exotic regions. Results of the BLASTN showed that most regions had the highest matches with plastome sequences of green plants. However, there were also matches with algae, bacteria, and animal sequences in regions corresponding to rRNAs (rrn4.5, rrn5, rrn18, and rrn23). Results of the BLASTP for open reading frames (ORFs) in the G. pubilabiata plastome confirmed that most ORFs were of green plant origin. In some cases, ORFs located near or within rpl16, rrn16, rrn23, and rps12 showed a high proportion of matches with bacterial sequences in BLASTP search results.

During annotation of protein coding genes in the plastome of G. pubilabiata, no debris from other mitogenome or nuclear genomes was found. However, during annotation of mitogenome, sequences like plastome derived psbF genes were detected, which was 119 bp in length. Regarding this regions, similar regions were found in other previously reported mitogenomes. These sequences were aligned and used to construct a maximum likelihood tree (Supplementary Figure S1 and Supplementary Figure S2). In this ML tree, all mitochondrial psbF-like regions formed a monophyletic group, indicating that mitochondrial and plastome psbF genes evolved independently.

The plastome and mitogenome of G. pubilabiata were searched using local BLASTN against two databases of published angiosperm plastomes and mitogenomes. As a result, 23 blocks of the G. pubilabiata plastome were found to match with the mitogenome database. Of these 23 blocks, 17 were related to rrn16 and rrn23, while the remaining were related to clpP, rpl2, rpl16, and ycf1. In the G. pubilabiata mitogenome, 118 blocks were identified as having plastome-related sequences, 82 of which had a higher degree of similarity to mitogenomes than to the plastome database. Thirty-six blocks were found to be plastome-related blocks, largely located in intergenic spacers or rRNAs. Six regions were found to be coding regions in the plastome (ndhC, ndhJ, ndhK, psaA, rpoB, and rps4) (Supplementary Table S3). A phylogenetic tree was constructed using these regions after performing MAFFT alignment (Figure 4).

2.3. Phylogenetic position of Gastrodia pubilabiata

The phylogenetic position of G. pubilabiata was determined using 79 protein coding regions and four rRNAs in its plastome. Results showed that G. pubilabiata, G. elata, and G. longystyla formed a monophyletic group. All nodes in the Gastrodia clade were supported by 100% bootstrap values (Figure 5). General phylogenetic position of the Orchidaceae was found to be Apostasioideae[Vanilloideae[Cypripedioideae[Orchidoideae, Epidendroideae]]].

As previously mentioned, studies decoding mitogenomes of orchids are limited. Thus, a phylogenetic tree was constructed using previously reported angiosperm mitogenomes. Results showed that G. pubilabiata and G. elata formed a monophyletic group and the genus Gastrodia formed a sister relationship with Allium cepa and Asparagus officinalis. Most species were supported by high bootstrap values (above 90%). However, species such as Beta, Daucus, Fagus, Geranium, Luffa, Nelumbo, and Prunus showed relatively lower bootstrap values. The phylogenetic tree based on the mitogenomes was in line with previously reported phylogenetic trees based on plastomes or nuclear genomes (Figure 6).

3. Discussion

3.1. Organelle genome evolution

The plastome of Gastrodia pubilabiata had a distinct feature in that it has extremely reduced genome size and consist of only a single copy region. This same plastome structure has previously been also reported in the same genus, G. elata [40]. Other taxa without IR include Orobanchaceae [50], Cactaceae [51], Geraniaceae [52,53], Fabaceae [54], and Gymnosperms [55,56]. In Fabaceae, a representative group of taxa lacking IR, deletions of genes or introns, rearrangements of the plastome structure, gene transfer, and hypermutations have been reported [57]. The lack of IR in the legume plastome might be a result of plastome rearrangement mediated by repetitive DNA [58]. These repeats are thought to have facilitated rearrangements by enabling recombination [53].

The plastome of G. pubilabiata contains simple sequence repeats (SSRs) located around rps2-rps14, rps2-rps3, and ycf2 regions (Supplementary Table S2). The gene trnQ-UUG is located between rps2 and rps14 and the trnC-GCA gene is situated near rps2 and rps3. Recombination of repeats and tRNAs might be the primary mechanism behind genomic rearrangements in Trachelium [59]. In a typical orchid plastome, rps2 is found near the LSC-IRa junction and rps3 is located near the LSC-IRb junction. Rps2 and rps14 are 20 kb apart. The change in position of rps2, rps3, and rps14 in G. pubilabiata is believed to be due to the presence of repeats or tRNAs in the intergenic spacer (IGS) near the relevant region. However, it is uncertain whether gene reduction, in which protein coding regions between corresponding genes are lost, precedes genome rearrangement and reduction [60]. Some studies have suggested that repetitive sequences contribute to genome rearrangement in the Geraniaceae [53]. Accumulation of SSRs and repeats in the plastome of G. pubilabiata might have also contributed to gene loss or plastome size reduction.

Several changes such as ndh gene contents, IR boundary shifts, and IR absence have been reported in plastomes of various orchid taxa [3,30,40,41,60]. For non-photosynthetic taxa, changes have also been reported in gene groups other than just ndh genes [41]. This level of change is consistent with previous reports for non-photosynthetic flowering plants [61]. In the case of G. pubilabiata, most genes have disappeared. The plastome of G. pubilabiata belongs to the fourth stage of plastome degradation pattern (Supplementary Figure S3). Other orchid taxa that belongs to stages 4 of plastome degradation, including Epipogium [39], Gastrodia [40], Rhizanthella [31], and Lecanorchis [41], are also non-photosynthetic mycoheterotrophs with greatly altered plastomes. These taxa have only retained some housekeeping genes.

Flowers and stems of G. pubilabiata have a green color, which sets them apart from the other fourth stages of orchid with red, white, or yellow flowers and stems. In contrast to the other green stemmed mycoheterotrophs C. macrorhizon, which only lacks ndh genes in its plastome [49], all photosynthesis related genes have disappeared from the confirmed plastome of G. pubilabiata (Figure 1 and Supplementary Figure S3). Despite this, low level photosynthesis has been observed in flowers and stems of C. macrorhizon due to residual chlorophyll [62,63,64]. Given this and the presence of residual plastids in other non-photosynthetic taxa, it is hypothesized that there might be an alternative pathway for photosynthesis in G. pubilabiata. Further research such as RNA-seq is needed to confirm this.

Mitogenomes of flowering plants are kwown to be variable in structure due to the presence of repeats, rearrangements, and submolecules [65]. Mitogenomes of monocots can vary greatly in length and genetic composition (Table 2 and Supplementary Figure S4). In the genus Gastrodia, multiple mitogenomes have been reported [40]. The mitogenome of G. pubilabiata was found to have 44 mitogenomes (Table 1 and Figure 2). Unlike the plastome, which displays a significant difference in genetic composition between photosynthetic and non-photosynthetic taxa, it is believed that gene composition of the mitogenome is not significantly impacted by photosynthesis or the lack thereof.

3.2. Putative gene transfer

BLAST search was carried out to detect any foreign sequences in the plastome and mitogenome of G. pubilabiata. The search uncovered regions with high plastome frequency having mitochondrial features and regions with high mitochondrial frequency having plastome characteristics (Figure 7). Regions with high frequency of mitogenome even in the plastome were mostly made up of rRNA regions (Supplementary Table S3). Other protein coding regions in the plastome (accD, clpP, rpl16, and rpl2) did show BLASTN results for the mitogenome. However, it was confirmed that there were many more search results for the plastome (Supplementary Table S3). The plastome might have infiltrated to the mitogenome, or there was a coincidental similarity of nucleotide sequences, rather than the presence of corresponding regions derived from the mitogenome. Reports of the presence of mitochondrial plastome sequences (MTPTs) have been documented in flowering plants [13]. A recent study has also reported the presence of plastome originating sequences in the mitogenome of Mangifera species [66].

BLASTN results for the mitogenome were like those for the plastome. Many results showing the high frequency for the plastome contained regions having rRNA, such as 16s rRNA, and 18s rRNA, suggesting that there might have been an exchange between rRNAs of the mitogenome and the plastome. However, for some regions, BLASTN results were found for plastome regions (ndhC, ndhJ, ndhK, psaA, psbF, rpoB, and rps4), but not for mitogenome references. This suggests that gene transfer from the plastome to the mitogenome might have occurred. The mitogenome contains many photosynthesis-related genes from the plastome (ndh, psb, psa, rpo), even though G. pubilabiata is a non-phtoosynthetic plant and all photosynthesis related genes have been deleted from the current plastome (Figure 1 and Supplementary Figure S3). This suggests that transfer of photosynthesis related genes to the mitogenome might have occurred before the loss of photosynthetic function in plastome. Therefore, this gene transfer might have occurred a longer time ago than previously expected. Except for a few obligate mycoheterotrophic taxa, a previous study have shown that the divergence of mycoheterotrophs might have occurred earlier than 30 mya [41]. This suggests that photosynthesis-related genes might have disappeared after that time. Gene transfer might have occurred a longer than 30 mya. However, as pointed out in that previous study, research results on the plastome and mitogenome of related obligated mycoheterotrophic taxa will be needed to verify that hypothesis.

A phylogenetic tree was constructed to determine the origin of gene fragments (ndh, psa, psb, and rpo) in the mitogenome of G. pubilabiata with other major flowering plants. Results are shown in Figure 4. In the case of psbF, G. pubilabiata formed a monophyletic relationship with other mitochondrial psbF. It tended to be distinguished between the plastome and the mitogenome. This suggests that these gene fragments might have been transferred from the plastome to the mitogenome in the distant past. In the case of psaA, the mitogenome derived psaA displayed a long branch, which was believed to result from the absence of other mitogenome derived sequences with the same cause as the psbF phylogenetic tree. The phylogenetic tree using rpoB region was distinct from trees using psaA and psbF gene fragments. It did not show a long branch and rpoB sequences of G. pubilabiata from a monophyletic branch with Euonymus of Celastrales, which was supported by a bootstrap value of 100. This suggests that G. pubilabiata has a close relationship with Celastraceae plants in its vicinity. Similar relationships have been reported for parasitic plants such as Aeginetia with nearby Poaceae plants [17]. Although it is difficult to directly compare mycoheterotrophic orchid to parasitic plants, this orchid might be related to surrounding plants. It was cornfirmed that Celastrus and Euonymus lived some distance away, not nearby.

Unlike photosynthesis-related genes found in the mitogenome, the rps4 gene is classified as a housekeeping gene. The presence of this gene in both the plastome and mitogenome of G. pubilabiata was confirmed and a phylogenetic tree was constructed to examine its origin. The phylogenetic tree showed that both plastome and mitogenome derived sequences of the rps4 gene formed a monophyly with other orchids. However, the plastome derived sequence showed a long branch and formed a monophyly, while the mitogenome derived sequence was at the base of the clade. This difference in evolutionary rates between the two sequences was likely due to relaxed selection in degraded plastomes of mycoheterotrophic plants, which tended to drive the acceleration of gene evolution [67]. Moreover, mitochondrial datasets shows much slower evolutionary rates than the plastome datasets in the previous study [68]. Thus, it could be concluded that the rps4 gene in the plastome evolved at a faster evolutionary rate than the gene in the mitogenome.

It is known that ndh genes are typically the first to be lost in the plastome degradation process of non-photosynthetic species, followed by other photosynthesis related genes [50]. This study found that ndhC, J, and K genes in the mitogenome of G. pubilabiata formed a monophyly with other orchids, unlike photosynthesis-related gene rpo, psa, or psb, which showed relationships with other taxa of flowering plants (Figure 4). This suggested that the degradation of G. pubilabiata’s plastome was more recent or that alternative pathways were involved. Reference guided assembly analysis only found a few base sequences for ndhC, J, K in the mitogenome. No base sequence was found to support the possibility of transfer to the nuclear genome. Hence, it was assumed that ndhC, J, K genes disappeared more recently from the plastome or that the degradation of the plastome was more recent. While it was believed that non-photosynthesis in Orchidaceae developed within the past 30 mya [41], further experiments such as transcriptome analysis are needed to explore the possibility of alternative pathways.

3.3. Phylogenetic implications

Phylogenetic relationships within the Orchidaceae are presented as Apostasioideae[Vanilloideae[Cypripedioideae[Orchidoideae, Epidendroideae]]] based on the analysis of protein coding regions and rRNAs of orchid plastomes. This finding aligns with a previous study [3]. Three species of the genus Gastrodia, G. elata, G. longyphyla, and G. pubilabiata, formed a monophyletic group and the tribe Nervilieae was the sister group to Gastrodieae (Figure 5). In the case of Gastrodieae, Nervilieae, Diurideae, and Vanilleae, containing non-photosynthetic orchids, a long branch was formed, which might be due to reduced number of genes used in the analysis or elevated mutation rates in mycoheterotrophic species. A phylogenetic tree based on mitochondrial genes indicated that G. pubilabiata was part of a monophyletic group along with other Asparagales species (Figure 6). Unlike the plastome phylogenetic tree, non-photosynthetic species exhibited limited levels of long branches.

4. Materials and Methods

4.1. Plant material and DNA extraction

Plant samples were collected from Jeju Island, South Korea. To protect plant habitat, only a few flowers were taken from individuals of Gastrodia pubilabiata. Fresh samples were ground into powder using liquid nitrogen in a mortar. Ground samples were then used to extract genomic DNA using a G-spin II Genomic DNA extraction Kit (Intron, Seoul, Korea) following the manufacturer’s manual. The quality of the extracted DNA was checked with a UV/VIS spectrophotometer. DNA was then stored in Plant DNA bank in Korea (PDBK2020-0086). Due to limitations of sample collection, the voucher specimen was replaced with photographs of G. pubilabiata’s in its habitat. These photographs were deposited in Korea University Herbarium (KUS) under the voucher-specimen number of KUS2020-0086 (Supplementary Table S4).

4.2. Sequencing, Assembly and Annotations

Extracted genomic DNA of G. pubilabiata was used for sequencing using two different NGS platforms: Illumina NovaSeq (Illumina, CA, US) and Pacbio Sequel (Pacbio, CA, US). A total of 224,926,562 NovaSeq reads, and 746,492 Pacbio Sequel reads were produced. The quality of sequencing reads was examined using FastQC v.0.11.9 [69]. Raw reads from Illumina NovaSeq were about 50 gb. They were trimmed with BBduk 37.64, as implemented in Geneious 11.1.5 (length: 27 kmer) [70]. BBNorm 37.64 was used to normalize trimmed reads (target coverage level: 30; minimum depth: 12). Raw reads from Pacbio Sequel were about 5 gb. They were corrected and polished using trimmed NOVAseq reads with CANU v.1.9 [71].

Error corrected Pacbio reads were used for de-novo assembly using FALCON v.0.3.0 [72,73] and CANU v.1.9. Mitogenome of Allium cepa (NC030100) and 18 mitochondrial sequences of G. elata (MF070084 ~ MF070102) were used as references to verify putative mitochondrial contigs among de novo assembly results obtained using FALCON and CANU. Trimmed NovaSeq reads and selected contigs from FALCON and CANU results were used for hybrid de novo assembly using SPAdes v.3.14.1 [74]. Mitochondrial contigs were selected from SPAdes results by performing local BLASTN [75] and using BANDAGE v.0.8.0. [76]. The final selection of mitogenome was based on their coverage, which was obtained using hybridSPAdes [77]. The depth-coverage of mitogenomes of G. pubilabiata is represented in Supplementary Table S5.

Plastome sequences of G. elata (NC037409) and Habenaria radiata (NC035834) were used as references to identify plastome contigs from results of FALCON and CANU. Selected contigs were assembled using Geneious assembler in Geneious 11.1.5 to obtain complete plastome sequence. Trimmed NovaSeq reads were then mapped to the complete plastome sequence of G. pubilabiata’s using Geneious assembler, and the mapped NovaSeq reads were extracted. Extracted NovaSeq reads were used to perform de novo assembly using Geneious assembler to validate G. pubilabiata’s complete plastome sequence.

The complete plastome sequence of G. pubilabiata was annotated using BLASTN, tRNAscan-SE 2.0 [78], ORF finder, and the find annotation function in Geneious 11.1.5. Plastome sequences of G. elata and H. radiata were used as references. Alternative start codons (ACG and TTG) were also considered in ORF finder. Pseudogene judgement was performed based on the criteria described previously [41]. ORFs that were larger than 150 bp were translated into protein sequences. These sequences were used to perform a psi-blast search [79] using several reported Orchidaceae plastome sequences as references. A circular plastome map was generated using the OGdraw web server [80].

4.3. Phylogenetic analysis

Ninety-two plastome sequences were downloaded from NCBI for phylogenetic analysis (Supplementary Table S6), including 87 sequences for the Orchidaceae, four for the Asparagales, and one for Liliales. Seventy-nine protein coding genes and four rRNA genes were extracted from these plastome sequences. Each extracted gene was aligned using MAFFT [81]. Alignments were manually checked. All alignments were concatenated to a length of 85,749 bp and subjected to jModeltest in the CIPRES Sience Gateway to determine the best fit model [82,83]. The best fit model was determined to GTR + I + G. Missing genes were treated as missing data. A maximum likelihood (ML) tree was constructed using RaxML-HPC2 on XSEDE in CIPRES Science Gateway with 100 bootstrap replicates [84]. The result ML tree was visualized using Treegraph2 [85].

Thirty-eight mitogenome sequences were downloaded from NCBI (Supplementary Table S6). Thirty-five coding genes and three rRNA genes were extracted from these mitogenomes. Each gene was aligned using MAFFT. All alignments were manually checked. Alignments were concatenated to a total length of 51.862 bp and then subjected to jModeltest in CIPRES Science Gateway to determine the best model. An ML tree was constructed using RaxML-HPC2 on XSEDE in CIPRES Science Gateway with a GTR + I + G model and 100 bootstrap replicates. The obtained ML tree was graphically represented using Treegraph2.

Eighteen mitochondrial sequences of G. elata (MF070084 ~ MF070102) and two plastome sequences of G. elata (NC037409) and G. longistyla (MW879162) were used to compare organelles among the genus Gastrodia. Simple sequence repeats (SSRs) were identified using Phobos 3.3.12 [86]. REPuter was used to distinguish palindromic repeats among organelles with default options [87]. The repeat finder in Geneious was used to find repeats larger than 50 bp with a perfect match option. Repeat information were summarized and presented graphically using Circos v.0.69-9 [88].

A total of 5,199 reported angiosperm plastome sequences and 219 mitogenome sequences were downloaded from NCBI and used as a local BLAST database. Plastome and mitogenome sequences of G. pubilabiata were divided into 150 bp fragments, which were then used to perform local BLASTN to identify potential instances of IGT or HGT. Results were summarized in terms of chloroplast counts and mitochondria counts based on the number of results generated. Summarized data were visualized using ggplot2 in R [89]. Based on these data, mitochondrial plastome (MTPT) sequences were extracted to construct a phylogenetic tree. Fifty-three plastome sequences of the main angiosperm lineage were used as references for this phylogenetic tree. For psaA and rpoB genes, a total of 301 Malpighiales and Celastrales plastome sequences were used as references to construct the ML tree. All obtained ML trees were visualized using Treegraph2.