A Comparative Genomics Approach for Analysis of Complete Mitogenomes of Five Actinidiaceae Plants

Actinidiaceae, an economically important plant family, includes the Actinidia, Clematoclethra and Saurauia genus. Kiwifruit, with remarkably high vitamin C content, is an endemic species widely distributed in China with high economic value. Although many Actinidiaceae chloroplast genomes have been reported, few complete mitogenomes of Actinidiaceae have been studied. Here, complete circular mitogenomes of the four kiwifruit species and Saurauia tristyla were assembled. Codon usage, sequence repeats, RNA editing, gene transfers, selective pressure, and phylogenetic relationships in the four kiwifruit species and S. tristyla were comparatively analyzed. This research will contribute to the study of phylogenetic relationships within Actiniaceae and molecular barcoding in kiwifruit.


Introduction
According to the endosymbiosis theory, the mitochondrion is an endosymbiotic alphaproteobacterium engulfed by the archaeal-derived host cell and eventually evolves into a semiautonomous organelle [1][2][3]. Mitochondria, known as energy factories, play a crucial role in numerous metabolic processes related to energy generation, synthesis, and degradation in living cells [4]. Mitochondrial DNA is maternally inherited in most seed plants [5]. With genome sequencing technology's rapid development, various complete organelle genomes in plants have been extensively studied [6]. Nearly 7576 chloroplasts and mitogenomes of land plants have been published in the National Center for Biotechnology Information (https://www.ncbi.nlm.nih.gov/genome/browse#!/organelles/) (accessed on 10 January 2022). The number of mitogenomes in land plants published was less than 20 before 2015 (Supplementary Figure S1). In recent years, it is no doubt that the number of land plant mitogenomes has increased significantly from 2018 to 2021 (Supplementary Figure S1). However, compared to the completed chloroplast genomes (7246), only 324 completed mitogenomes were assembled (Supplementary Table S1), suggesting that the interpretation and functional annotation of the mitochondrial genome is complex in comparison to other organelles.
The intergenomic DNA transfers and highly dynamic, multipartite structures of plant mitogenomes may make it challenging to build plant mitogenomes [7]. Several articles have reported that most plant mitogenomes range from 200 to 2000 kb in size [8]. Differences in mitogenome size can be attributed to repetitive sequences and foreign DNA derived from other organisms during evolution [9]. Many intramolecular recombination events and subgenomic conformations have been found in some land plants, such as Scutellaria tsinyunensis [10], Cucumis sativus [11], Ipomoea batatas [12], and Brassica napus [13]. In extreme environments, gene loss and RNA edits may occur during plant mitogenome rearrangement [14]. Several separate chromosomes can be found in some higher plant mitochondrial genomes. The cucumber mitogenome, for instance, has three separate chromosomes [11]. The mitogenomes of Globodera ellingtonae and Camellia sinensis have two separate chromosomes [15,16]. The mitogenomes of the plant's seeds contain many repeating sequences, including simple sequence repeats (SSRs), tandem repeats, and scattered repeats. In addition, there are also many insertions/deletions (indels) and single nucleotide polymorphisms (SNPs) within mitogenomes [17,18]. SSRs and SNPs have been widely applied to identify species rapidly and for phylogenetic plant analyses, especially in Chinese herbal medicine classification [19,20]. Moreover, it is an essential feature for mitogenome evolution via intracellular transfer between the mitochondria and the chloroplast genomes [21]. Most of the transferred sequences are transferred from the nucleus to the mitochondria, but several chloroplast-derived tRNA genes are transferred to the mitochondria and perform essential functions [22]. The horizontal gene transfer (HGT) phenomenon also plays a significant role in the evolution of plant mitogenomes [9]. These findings suggested the existence of instability in higher plants' mitogenome structures. Finally, a long-reads strategy in combination with short-reads technologies (Pacbio SMRT, Oxford Nanopore, or Illumina mate-pair) were applied to solve the problem caused by this structural instability in mitogenome assembly.
Actinidiaceae, an economically important plant family, includes the Actinidia, Clematoclethra, and Saurauia genus [23,24]. Among the Actinidiaceae family of the Asterids, kiwifruit with remarkably high vitamin C content, commonly known as 'the king of fruits', is an economically important horticultural fruit tree. Kiwifruit is widely cultivated in Asia, Europe and Oceania (https://www.fao.org/faostat/zh/#data/QCL, Supplementary Figure S2) (accessed on 3 December 2021). Worldwide annual kiwifruit production increased rapidly from 2012 to 2020 and reached approximately 2 million tons in 2020 (https: //www.fao.org/faostat/zh/#data/QCL, Supplementary Figure S3) (accessed on 3 December 2021). Total kiwifruit production in Asia was the highest, accounting for 52.5%, followed by Europe (25.3%) from 2012 to 2020 (https://www.fao.org/faostat/zh/#data/QCL, Supplementary Figure S4) (accessed on 3 December 2021). It is noted that the annual kiwifruit production in China was the highest in Asia and reached up to 1.49 million tons in 2020 (https://www.fao.org/faostat/zh/#data/QCL, Supplementary Figure S5) (accessed on 3 December 2021). This may be due to the abundant kiwifruit germplasm resources in China. So far, diploid A. chinensis and hexaploid A. chinensis var deliciosa are the most commercial kiwifruit varieties. Abiotic and biotic stresses, including drought, salinity, low or high temperatures, and Pseudomonas syringae pv. actinidiae (Psa) seriously affect the yield and quality of kiwifruit [25]. After incidence of Psa, there is no remedy available to control it, except for destroying the tree to prevent the spread of the disease. Thus, Psa seriously threatens the production and development of the kiwifruit industry worldwide. It has been reported that A. eriantha var 'huate' and A. chinensis var deliciosa 'jinkui' strongly resist Psa [26,27]. Mitochondria genetic engineering would be beneficial in developing a method of resilience to abiotic and biotic stresses [25]. Hence, the complete mitogenomes of kiwifruit are sequenced, providing great promise for breeding kiwifruit cultivars with resilience to abiotic and biotic stresses.
For the last three decades, 4 kiwifruit nuclear genomes and over 29 complete chloroplast genomes from the Actinidiaceae family have been sequenced [28][29][30][31][32], while no complete mitogenome of this family has been reported previously. To elucidate the evolutionary mechanisms and structural features that underlie the Actinidiaceae family's mitogenomic diversity, the complete mitochondria genome of the diploid A. chinensis, and the tetraploid A. chinensis, hexaploid A. chinensis var deliciosa, A. eriantha and S. tristyla were sequenced and assembled in this study. The mitochondria genome with a two-chromosomal conformation was found in diploid A. chinensis, tetraploid A.chinensis and S. tristyla. The genome size (939 kb) of A. chinensis var deliciosa was significantly more extensive than other Actiniaceae species. Therefore, we hypothesized that the mitogenome may experience expansion during A. chinensis ploidy doubling. We analyzed the mitogenome structures of four kiwifruit species and S. tristyla to elucidate/unveil the genomic repeats, RNA editing sites, relative synonymous codon usage, gene transfer, and the evolutionary relationships among the Actinidiaceae family. To sum up, our study will be instrumental for genetic engineering and breeding programs. Table S3 shows details of the tested materials. Fresh leaves were wrapped in aluminum foil, flash frozen in liquid nitrogen, and stored at −80 • C for subsequent use. High-quality total genomic DNA was extracted using a DNAsecure Plant Kit (Tiangen Biotech, Co. Ltd., Beijing, China). The DNA library construction and sequencing were performed as previously reported by Emerman et al. [33].

RNA Editing Predicting and Codon Usage
We used the online PREP-Mt suite of servers (http://prep.unl.edu/, accessed on 5 March 2022) [48], with a cutoff value of 0.2, to predict the RNA editing sites of the 39 protein-coding genes of the 5 mitogenomes. The relative synonymous codon usage (RSCU) was calculated by MEGA X53 [49].

Substitution Rate Calculation and Phylogenetic Inference
Pairwise 19 protein-coding gene sequences of the mitogenomes of Actinidia were used to estimate the pairwise nucleotide substitution rates, including the non-synonymous substitution rate (Ka) and synonymous substitution rate (Ks), and the ratio of Ka to Ks. The Ka/Ks ratios were calculated by PAML (v4.9) [50] using the yn00 module with default parameters. The Ka/Ks values' heatmap was plotted using Tbtools [47]. In order to further analyze the phylogenetic position of the Actinidiaceae species, 23 plant mitogenomes from GenBank were downloaded for phylogenetic tree construction. A total of 20 orthologous mitochondrial genes were identified and extracted using PhyloSuite (v1.2.1) [51]. The corresponding nucleotide sequences were aligned using MAFFT (v7.450) [52] implemented in PhyloSuite. The phylogenetic tree was constructed using the maximum likelihood (ML) method via RAxML v8.1.5 [53] with 1000 bootstrap replicates. Furthermore, the web iTOL (https://itol.embl.de, accessed on 8 April 2022) [54] was used to visualize the phylogenetic trees.

Mitogenome Assembly, Annotation and Gene Features
The de novo assembly assembled five complete mitogenomes of Actiniaceae species. The de novo genome assembly yielded a single circular molecule for A. chinensis var deliciosa (939 kb) and A. eriantha (768 kb). In contrast, two distinct circular chromosomal genomes for A. chinensis (2×), A. chinensis (4×) and S. tristyla ( Figure 1A-D and Supplementary Table S2) were reported. A. chinensis (2×) and A. chinensis (4×) mitogenomes exhibited similar genome size (916 kb vs. 907 kb). The genome size (939 kb) of A. chinensis var deliciosa was significantly more extensive than the other Actiniaceae species. S. tristyla (482 kb) has the smallest genome size among them. Interestingly, their mitogenomes, containing similar GC content, are about 46% (Supplementary Table S2). A comparison of the annotated genes in Actiniaceae revealed that the pseudogene rps2 is absent in tetraploid and hexaploid A. chinensis (Supplementary Figure S2). The loss of the sdh4 gene in the A. arguta mitochondrion genome was also notable (Supplementary Figure S2).

Repeat Sequences and Chloroplast-Derived Region Analysis
As shown in Figure 2A, a total of 49-124 tandem repeats were found in the Actiniaceae mitogenome. The number of tandem repeats, comprising between 10 and 20 bp in the A. chinensis (2×) mitogenome, was significantly lower for the other species. Compared to the tetraploid and hexaploid A. chinensis, the number ranged from 40 bp to 105 bp in diploid A. chinensis. About half of the tandem repeats ranged from 10 bp to 20 bp in kiwifruit, whereas most of the tandem repeats ranged from 41 bp to 105 bp in S. tristyla. A. eriantha contained the fewest tandem repeats. A wealth of SSRs was identified ( Figure 2). SSRs in hexamers were discovered in all species except S. tristyla ( Figure 2C). Nearly 78% of the SSRs belonged to monomers and dimers ( Figure 2B). Tetramers and pentamer-nucleotide repeats were less frequent in the Actiniaceae mitogenome. A. chinensis mitogenomes contained a higher number of SSRs than those of S. tristyla and A. eriantha. S. tristyla had the lowest number of SSRs ( Figure 2B).

RNA Editing Sites and Codon Usage Analysis of PCGs
As shown in Figure 4, the RNA editing sites of 39 PCGs of the mitogenomes of 5 Actinidiaceae plants were predicted in this study. Three cytochrome c biogenesis genes, including ccmFn, ccmB, and ccmC, displayed the most RNA editing sites in five Actinaceae sp. plants. Interestingly, we found that the number of NAD1 gene RNA editing sites in A. chinensis var deliciosa was significantly higher than in the other species (Figure 4). Only the rpl2 gene in the A. chinensis (4×) had no RNA editing sites (Figure 4). The rpl16, rps1, and rps2 genes contained the same number of editing sites in the A. chinensis (2×) and A. eriantha, but not in the other species (Figure 4).  The codon distribution and relative synonymous codon usage (RSCU) of five Actinidiaceae species' mitogenomes were analyzed. The RSCU analysis showed that Leu, Ser, and Arg appeared the most frequently, whereas those that encoded Met and Trp were relatively less abundant in five Actinidiaceae species' mitogenomes. Five species in the Actinidiaceae family share a similar RSCU style ( Figure 5A-E).

The Synonymous and Nonsynonymous Substitution Rate (Ka/Ks) and Phylogenetic Analysis
Nineteen protein-coding genes of six Actinidiaceae mitogenomes were used to calculate the Ka/Ks ratios. As shown in Figure 6, we observed that the sdh3 gene had an abnormally high Ka/Ks ratio > 1 compared to the other genes between S. tristyla and the kiwifruit, indicating possible positive selection. The Ka/Ks values of most PCGs were less than 1, such as atp9, ccmB, ccmC, cox3, nad6 and rps12, indicating that most PCGs were under purification selection ( Figure 6). These results suggested that most PCGs may be highly conservative in the evolutionary process of Actinidiaceae.  In order to further analyze the phylogenetic position of Actinidiaceae, 23 plant mitogenomes from GenBank were downloaded for phylogenetic tree construction based on 20 PCGs. Phylogenetic analysis showed that 23 plant mitogenomes were divided into 6 categories (Figure 7). We selected V. vinifera and N. nucifera as outgroups. The phylogenetic tree strongly demonstrated that five kiwifruit species (A. chinensis (2×), A. chinensis (4×), A. chinensis var deliciosa (6×), A. eriantha and A. arguta) clustered into one clade with a 100% bootstrap value (Figure 7). It also revealed that S. tristyla was closely related to five kiwifruit species (Figure 7).

Discussion
So far, 25 Actinidia genus chloroplast genomes have been reported [32], and our group has comprehensively analyzed them. Unlike conserved genome structures and small Actinidia chloroplast genomes, Actinidia mitogenomes generally have multiple different sizes and structural variations [55]. This makes Actinidia mitogenome research relatively challenging. Here, complete mitogenomes of five Actinidiaceae plants were sequenced and assembled. Two subgenomic circles were found in A. chinensis (2×-4×) and S. tristyla. Similar results have been reported in C. sinensis var. assamica and C. assamica [16,56]. The GC content of five Actinidia mitogenomes is about 46%, similar to most other mitogenomes [57]. Among the observed size variations, the genome size in A. chinensis var deliciosa (6×) was about twice that of the S. tristyla (482 kb) and the closely related species Vaccinium macrocarpon (459 kb) [58]. The genome size of the hexaploid A. chinensis var deliciosa was nearly 20 kb larger than that of diploid and tetraploid A. chinensis, which is probably the result of a gradual increase in sequence duplication and intracellular transfer of the plastid or nuclear genome or horizontal transfer of mitochondrial DNA during evolution [59,60].
Repeat sequences widely exist in plant mitogenomes, including tandem repeats and SSRs [61]. Positive correlations between genome size and repeat sequences were identified in 38 Rosaceae mitogenomes [62]. As shown in Figure 2A, tandem repeats from 10 to 30 bp were the most abundant for the Actinidiaceae mitogenomes, with similar results in Diospyros oleifera [63]. Guo et al. [64] have reported that SSRs played a pivotal role in intermolecular recombination during evolution. The number of SSRs in the A. chinensis mitogenome (18.24%) was higher than that of S. tristyla and A. eriantha ( Figure 2B), which may cause the A. chinensis mitogenome size to be larger than S. tristyla and A. eriantha. It is consistent with the findings of previous studies [65,66]. Dimer repeats were the most abundant SSR type (about 48%) in the Actinidia mitogenome ( Figure 2C), which is commonly found in Suaeda glauca [67]. Gene transfer from the chloroplast to the mitogenome frequently occurs during long-term plant evolution [68]. A total of 9-55 kb of plastid-derived sequences was observed, which occupied 3-6% of the Actinidiaceae mitogenomes ( Figure 3). Similar results have been reported by Adams et al. [69]. Some plastid-derived protein-coding genes (cp-derived PCGs), such as rpoC1, ndhB, rps7, rps19, and rpl23, were identified in the A. chinensis (2×-4×) mitogenomes, which is commonly found in angiosperms [70]. In addition, we also found that psbJ, petL and petG cp-derived genes only exist in hexaploid A. chinensis var deliciosa (Supplementary Table S2), suggesting that special evolutionary events may have occurred during genome evolution.
Mitochondrial gene expression may be affected by RNA editing [71]. The number of RNA editing sites varies in different species [72]. Previous studies identified approximately 491 RNA editing sites within 34 genes in rice [68] and 486 RNA editing sites within 31 genes in Primula vulgaris [73]. In Actinidiaceae, the number of RNA editing sites in most PCGs was extremely conserved in Actinidiaceae, similar to other plant studies [74]. Interestingly, we also observed that the number of NAD1 gene RNA editing sites increased with ploidy in A. chinensis (2×-4×) (Figure 4). Whether the number of RNA editing sites is positively correlated with the ploidy of the kiwifruit requires further research. Relative synonymous codon usage (RSCU) refers to the relative probability of a specific codon between the synonymous codons that encode the corresponding amino acid [75]. The RSCU value showed that the codon usage pattern in the Actinidiaceae plants' mitogenomes shared a similar RSCU style ( Figure 5A-E), which was commonly found in higher plant mitogenomes [76].
Calculating the Ka/Ks ratio plays a vital role in understanding the dynamics of molecular evolution [77]. As shown in Figure 6, the Ka/Ks ratios of most PCGs were less than 1 (Figure 6), suggesting that these genes were highly conserved and had undergone neutral and negative selections in Actinidiaceae. These results were also supported in the report on Lamiales [78]. However, among the five Actinidiaceae plants, sdh3 with dN/dS values greater than 1.0 was found between the S. tristyla and kiwifruit mitogenomes, indicating that this gene may have suffered from positive selection during the evolution in Actinidiaceae. The phylogenetic trees in this study showed a close relationship between the Actinidiaceae and other Ericaceae plants (Figure 7), as Wang et al. [56] proposed. Notably, our analyses also demonstrated that the Actinidiaceae is monophyly, with the sampled five Actinidia taxa clustering in a clade as a sister to S. tristyla (Figure 7), in agreement with the result of Wang et al. [32]. Moreover, A. chinensis (2×-4×) was closely related to A. chinensis var deliciosa (6×) (Figure 7), which was consistent with the results of a previous study [79].

Conclusions
The large size variation in Actinidiaceae mitogenomes appeared due to increasing sequence duplication and intracellular transfer of the plastid. The number of RNA editing sites and codon usage in most PCGs of five Actinidiaceae plants' mitogenomes were highly conserved. Most of the coding genes had undergone negative selection, indicating the conservation of mt genes during evolution. We found that sdh3 may have suffered from positive selection during the evolution in Actinidiaceae. Kiwifruit species showed high similarities and were highly similar to S. tristyla and A. chinensis (2×-4×) was closely related to A. chinensis var deliciosa (6×). This study provides important mitochondrial genome resources for the Actinidiaceae species and has deepened our understanding of organelle genome evolution in flowering plants.

Supplementary Materials:
The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/genes13101827/s1. Table S1. Statistics information of sequenced land plant complete mitogenomes; Table S2. Features of five mitogenomes that belong to the Actiniaceae family; Table S3. Taxonomic information, GenBank accession numbers and collection place of mitogenomes used in the study; Figure S1. Statistics information of sequenced land plant complete mitogenomes from 1993 to 2022; Figure S2. Production quantities of kiwifruit by country; Figure S3. Yield quantities of kiwifruit in the world from 1999 to 2020; Figure S4. Production share of kiwifruits by region; Figure S5. Production of kiwifruit, including the top 10 producers.