Sequencing results
The initial assembly strategy assumed using only PE and PacBio long reads, which should overcome issues with plastid derived regions as well as most of long repeats which form the most sequenced plant mitogenomes not exceeding 1-10 kbp range [2]. However, obtained PacBio sequencing results (2 Gbp) revealed mitochondrial reads up to 2 kbp and coverage ranged from 4-10x, which only resolved part of issues related with long repeats and plastid derived regions. Better results were obtained with the nanopore sequencing, where 24-hour runs enabled the obtaining over 8 Gbp of data. This resulted in 60x coverage of single copy mitochondrial regions and enabled proper orientation of repeat flanking contigs.
Mapping PE read onto assemble mitogenome revealed mean 531x coverage for single copy regions and relevelant multiply for repeat regions (up to 5000 x coverage).
Mitogenome structure
Obtained results do not support the master chromosome hypothesis either in circular or linear form. Hybrid assembly using ONT and Illumina reads enabled to obtain three mitochondrial contigs (named chromosomes chMt1-chMt3) flanking by overlapping repeats A - D, which were also found inside chromosomes chMt1 and chMt2 (Figure 1).
The main role in shaping the structure of Pulsatilla mitogenome are up to 87 kbp- long repeats (Figure 1). The longest of them, Repeat D is located at 5’ end of chromosome chMt1 and 3’ end of chromosome chMt2. This repeat also contains remaining long repeats including almost 39 kbp- long Repeat C, 13 kbp Repeat B and 8 kbp- long Repat A (Figure 2).
The Repeat C is flanking the 5' end of chromosome chM3t and appears also inside chromosome chMt1. The 5’ end of chromosome chMt2 is flanked by Repeat B, which is also an internal part of chromosome chMt1, while remaining ends of chromosomes are flanked by Repeat A.
In the case of Pulsatilla patens the linearity of chromosomes chMt1 and chMt2, which contain large repeats A, B and C was supported only by six and eleven nanopore reads respectively, but junctions were positively validated by long-range PCR amplification. However designed PCR primers to single copy regions closest to flanking repeats did not amplify any detectable products, which may support circular master chromosome hypothesis. The fact that single circular molecule is an oversimplified representation of the plant mitochondrial genomes and that they rather exist in vivo as a mix of circular, linear and branched forms is well documented [8, 13, 14]. However, a circular structure can usually be observed at the level of sequence assembly due to the presence of multiple repeats. This is not the case for P. patens where internal repeats are also present but their location does not allow to reconstitute the master circle chromosome. It’s also indirectly supported by unequal short and long reads coverage distribution along the main flanking repeat. The linear of P. patens mitochondrial chromosomes is also supported by features usually found in other species with linear single or multichromosomal mitochondrial genome: a terminal inverted repeat and the genes of RNA and DNA polymerases [15]. The multichromosomal, linear mitogenomes were not previously reported in early eudicots, since Aconitum, Anemone and Nelumbo support master circular chromosome hypothesis [16, 17].
Multichromosomal architecture of the mitochondrial genome was previously reported in several tracheophytes species and can form circular, linear or branched molecules. The diversification of mitogenome structures appear in the earliest tracheophyte lineages, while the evolutionarily older land plants, bryophytes exhibit single circular mitogenome molecules [18]. Starting with lycophytes, where single circular mitochondrial chromosome were confirmed only for Huperzia [19, 20, 21], non-circular and multichromosomal mitogenomes are widely distributed among different fern, gymnosperm and angiosperm lineages. Two circular mitochondrial chromosomes were found in early divergent fern Psilotum nudum, while in Ophioglossum californicum supports master circle chromosome hypothesis [22]. The sequenced gymnosperm mitogenomes of genera Cycas, Ginkgo, Welwitschia and Taxus could be also assembled into single circular molecule [22, 23, 24], but the larger, over 5 Mbp mitogenomes of Picea species are described as multichromosomal and linear [25, 26]. The mitogenomes of angiosperms are structurally variable at genus or even species level. The largest and smallest known angiosperm mitogenomes belong to the species of Silene, which exist in single or multichromosomal, but usually circular forms [10]. Single and mulichromosomal circular mitogenomes were reported for variety of Allium cepa [27], and linear, branched or circular mitochondrial molecules were found in single individual of Lactuca sativa [8]. Out of three mitochondrial chromosomes of Solanum tuberosum two could be assembled as circular, but the third exist in the linear form [12].
The presence of master circular mitogenome structure can not be completely excluded. The assembly of Chrysantheum nakingense mitogenome using similar approach confirmed the master chromosome structure relying only on 4 nanopore reads, but modules were not flanked by long inverted repeats [14], as it was in the case of P. patens.
Gene content
Annotation of Pulsatilla patens mitogenome enables the identification of a set of genes (rRNA, tRNA and CDS) that is typical for angiosperms (Table 1). According to classification of [3] sequenced mitogenomes contain all “core” genes (atp1, atp4, atp6, atp8 atp9, ccmB, ccmC, cmmFc, ccmFn, cob, cox1-3, matR, mttB, nad1-7, nad9) and most of genes defined as “variable” (rpl1, rpl5, rpl10, rpl16, rps1-4, rps7, rps11-14, rps19, sdh4). These genes were also found in the three remaining Pulsatilla species including P. alpina, P. pratensis and P. vernalis. From the “variable” gene list only rps10 and sdh3 were not found in the Pulsatilla mitogenome. The fragments of sdh3 including exons 1 and 2 were found at chromosome chMt3, but they contain two internal stop codons and cover ca. 40% of functional genes. Both genes were not detected in RNA-seq analysis, so they were unlikely to be transferred to the nucleus. The presence of this gene pair (rps10 and sdh3) in plant mitogenomes seems to be correlated, since out of 18 cases of rps10 or sdh3 absence, in 15 cases both are missing [3]. Despite using deep sequenced libraries enabling nuclear genome coverage above 40X, the presence of these genes was not confirmed even in nuclear genomes.
The repeat units in the P. patens mitochondrion contain three protein coding genes: nad3, nad5 (exons 1-3) and rps12. All three genes are located in the Repeat A, over 48 kbp long duplicated region, resulting in two identical copies in the mitogenome. The different copy of rps12 is located in a single copy region with pairwise identity 81.9%. The RNA-seq analysis confirmed the expression of both rps12 variants. These genes are also duplicated or even triplicated in few known mitogenomes including Cycas taitungensis (AP009381), Oryza minuta (NC_029816) and Vitis vinifera (NC_012119). In the case of the mitogenome of Daucus carota (NC_017855), the rps12 gene is triplicated, but opposite to Pulsatilla mitogenome, all three copies are identical. The partial duplication of nad5 gene was also described in mitogenomes of Beta vulgaris (BA000024), Cynanchum auriculatum (MH410146), Rhazya stricta (KJ485850) and Tamarindus indica (MN017227).
Pulsatilla patens mitogenome consists of two large ORFs encoded polymerases, which were recently reported in plant mitogenomes [8, 14]. The 2,850 bp long Pp_DNA_pol_B (encoding DNA polymerase) gene and 2,154 bp Pp_RNA_pol (encoding DNA-dependent RNA polymerase) are located at the chromosome chMt1.
The length of Pp_DNA-pol_B gene is similar in size to previously reported Ac_DNA_pol_B gene (2,814 bp) found in Actinidia chinensis but the homolog of Pp_RNA_pol gene of this species is remarkably shorter counting 1,743 bp.
Mining available angiosperm mitogenomes towards polymerases encoding genes revealed their presence in over 50 accessions in the case DNA polymerase and 18 accessions in the case of RNA polymerase assuming blastx identity threshold above >50%.
Both genes were reported from several mitogenomes, but the expression of them remained unconfirmed in most of the studies. Angiosperm mitochondrial genomes contain DNA derived from exons of the nuclear genome [28, 29, 30], but the examples of expression are rather rare. The expression of nuclear derived gene arf17 was confirmed in mitogenomes of Arabidopsis [31], but it's rather a rare exception. The RNA-Seq analysis confirmed the presence of mRNA of polymerase coding genes in all three analysed samples. However, since availability of nuclear genomic resources of P. patens is limited to NOR regions [32, 33], the nuclear origin of these transcripts can not be ruled out.
The presence of polymerase coding genes seems to be limited to P. patens. Both genes are absent in the three other Pulsatilla species, including P. alpina, P. pratensis and P. vernalis as well as in more distant species of Ranunculales, Anemone maxima and Aconitum kusnezoffii. Moreover, among eight analysed specimens of P. patens, the presence of polymerase like genes was not confirmed in individuals from populations STR6 and a115. The frequency of intraspecific presence polymerase-like genes was not analysed in previous studies, but the lack of these genes in samples from southern populations [34] could indicate its recent transfer from nucleus.
Table 1. Gene contents, localization and RNA editing of Pulsatilla patens mitogenome
Category
|
Gene
|
Location
|
Gene structure
|
RNA editing events
|
RNA editing frequency
|
CDS Length nc/aa
|
GC content
|
Complex 1
|
nad1
|
M1(e1);M2(e2-3);M3(e4-5)
|
e1;e2-i1-e3;e4-i2-e5
|
45
|
3,578732106
|
978
|
43,6
|
|
nad2
|
M1(E1-E2); M3(E3-5)
|
e1-i1-e2;e3-i2-e4-i3-e5
|
43
|
2,810457516
|
1530
|
41
|
nad3
|
ChM3
|
continous
|
19
|
5,322128852
|
357
|
42
|
nad4
|
ChM1
|
e1-i1-e2-i2-e3-i3-e4
|
65CDS + 3 introns
|
4,36827957
|
1488
|
42,7
|
nad4L
|
ChM1
|
continous
|
12
|
4,395604396
|
273
|
37,4
|
nad5
|
ChM1 (E1-2), CHM2 (E3-E4)
|
e1-i1-e2,e3-i2-e4
|
39
|
1,937406855
|
2013
|
41,8
|
nad6
|
ChM1
|
continous
|
17
|
2,791461412
|
609
|
41,2
|
nad7
|
ChM2
|
e1-i1-e2-i2-e3-i3-e4-i4-e5
|
40CDS + 6 introns
|
3,410059676
|
1174
|
44,6
|
nad9
|
ChM2
|
continous
|
11
|
1,919720768
|
573
|
42,9
|
Complex 2
|
sdh4
|
ChM1
|
continous
|
8
|
1,93236715
|
414
|
38,9
|
Complex 3
|
cob
|
ChM2
|
continous
|
22
|
1,842546064
|
1194
|
42,5
|
Complex 4
|
cox1
|
ChM1
|
continous
|
30
|
1,890359168
|
1587
|
44
|
cox2
|
ChM1
|
e-i-e
|
18-cds, 1 intron
|
0,5323868678
|
3381
|
49,4
|
cox3
|
ChM1
|
continous
|
20
|
2,506265664
|
798
|
45,2
|
Complex 5
|
atp1
|
ChM1
|
continous
|
3
|
0,1937984496
|
1548
|
44,3
|
atp4
|
ChM1
|
continous
|
12
|
2,116402116
|
567
|
42,5
|
atp6
|
ChM1
|
continous
|
25
|
3,267973856
|
765
|
38,7
|
atp8
|
ChM1
|
continous
|
9
|
1,851851852
|
486
|
40.1
|
atp9
|
ChM1
|
continous
|
1
|
0,4444444444
|
225
|
37,3
|
Cytochrom c biogenesis
|
ccmC
|
ChM1
|
continous
|
40
|
5,376344086
|
744
|
44,4
|
ccmB
|
ChM1
|
continous
|
42
|
6,763285024
|
621
|
42,7
|
ccmFC
|
ChM1
|
e1-i-e2
|
25
|
1,081782778
|
2311
|
49,7
|
ccmFN
|
ChM1
|
continous
|
34
|
1,960784314
|
1734
|
47,3
|
Ribosome large subunit
|
rpl2
|
ChM1
|
e1-i-e2
|
1
|
0,06600660066
|
1515
|
51
|
rpl5
|
ChM1
|
continous
|
10
|
1,792114695
|
558
|
44,1
|
rpl10
|
ChM1
|
continous
|
6
|
1,19760479
|
501
|
42,9
|
rpl16
|
ChM1
|
continous
|
10
|
2,096436059
|
477
|
45,3
|
Ribosome small subunit
|
rps1
|
ChM1
|
continous
|
5
|
0,7278020378
|
687
|
42,8
|
rps2
|
ChM1
|
continous
|
10
|
1,481481481
|
675
|
38,2
|
rps3
|
ChM1
|
e1-i-e2
|
16
|
0,9523809524
|
1680
|
43,2
|
rps4
|
ChM3
|
continous
|
31
|
2,425665102
|
1278
|
39,3
|
rps7
|
ChM1
|
continous
|
3
|
0,6711409396
|
447
|
43
|
rps11
|
ChM1
|
continous
|
4
|
0,8888888889
|
450
|
44,2
|
rps12
|
ChM1,ChM3
|
continous
|
11
|
2,91005291
|
378
|
45,2
|
rps13
|
ChM2
|
continous
|
6
|
1,709401709
|
351
|
39,9
|
rps14
|
ChM2
|
continous
|
2
|
0,7575757576
|
262
|
40,2
|
rps19
|
ChM1
|
continous
|
4
|
1,388888889
|
288
|
39,6
|
Translocation pathway
|
mttB
|
ChM1
|
continous
|
43
|
5,449936629
|
789
|
44,7
|
Maturases
|
matR
|
ChM3
|
continous
|
19
|
0,9639776763
|
1971
|
52,5
|
HGT
|
DNA pol
|
chm1
|
continous
|
0
|
0
|
2850
|
38,7
|
|
fasciclin-like
|
chm3
|
continous
|
0
|
0
|
396
|
50,8
|
|
RNA pol
|
CHM1
|
continous
|
0
|
0
|
2154
|
33,6
|
Blasting longer orfs (>300 bp) also revealed other genes that are seldom or non-reported from mitogenomes of angiosperms. The 396 bp long orf from chMt3 chromosome was identified as fasciclin-like arabinogalactan protein 21 with the greatest similarity to sequences found in genomes of Medicago truncata (99% query coverage, 96.2 percent identity) and species of genus Arachis (97% query coverage and 83-85.4% percent identity). Unannotated gene was also found in the three chromosomes of Cicer arietinum (100% QC, 98.7 PI). However, all blastx (and blastn in case of Cicer) hits refer to nuclear genomes of Fabales species, suggesting that this gene was recently acquired by intracellular HTG from nucleus or was previously unannotated in known mitogenomes. Blasting against known mitochondrial genomes revealed the presence of this gene (97-100% query coverage and 99-95% percent of identity) not only in mitogenomes of Fabales (Lotus japonicus, Tamarindus indica, Styphnolobium japonicum and Medicago truncata), but also in Paraprenanthes diversifolia and Lactuca species of Asterales. The phylogenetically scattered presence of this gene in angiosperm mitogenomes supports nuclear origin and its transfer to mitogenome independently in at least three angiosperm orders: Fabales, Asterales and reported for the first time Ranunculales.
The chromosome chMt1 of mitochondrial genome of Pulsatilla patens contains orf encoding RNA-dependent RNA polymerase (RdRP) and an expression of this gene was confirmed by RNA-seq analysis. The RdRp gene is recognized by its conserved protein domain family, pfam05919 and is presumably required for replication of mitoviruses [35] and was found in 40 out of 50 analysed plant mitogenomes [36]. Mitoviruses have small, nonsegmented RNA genomes and all are through to replicate in the host mitochondria [37, 38]. They were previously known from fungal hosts [39] and their detection in plant transcriptomes was considered as a contamination from fungal pathogens or as a result of HGT [36]. However, the recent studies suggest an existence of the contemporary plant mitoviruses [35]. The completeness of RdRP gene in plant mitochondrial genomes varies from nearly complete versions of the RdRP to remnants barely detectable in sequence searches [36]. The RdRP gene found in the Pulsatilla mitogenome is 292 aa long and belongs to the longest among those found in plants. The further studies on Pulsatilla transcriptomes are required to confirm the presence of a complete mitovirus genome, which can be up to 4.4 kbp long [35].
The 328 aa long orf encoding reverse transcriptase with retrotransposon gag protein domain pfam03732 were found on chromosome chMt3 (blastx up to 90% of query coverage and >50% identity), however blasting against known plant mitogenomes did not provide any hits. The discovery of this gene encouraged us to perform LTR search in the assembled mitogenome, which are often considered as HGT vectors among nuclear and mitochondrial genomes [40]. Several orfs longer than 500 bp remain unidentified, as no significant hits were found in public databases, despite their expression end RNA editing events. However, the genomic resources of early eudicots are poorly explored.
Primary annotation of rpl2 gene did not reveal presence of stop codon at the end of the second exon. In the most angiosperms the stop codon of mitochondrial rps2 is generated by C->U RNA editing process [41, 42, 43]. However, despite deep reads coverage, the editing required to develop stop codon in rpl2 was not detected even at low frequency. This gene in various angiosperm lineages is known by different intermediate stages of its transfer to the nuclear genome [44]. It has been proposed that part of rpl2, located between the intron and 3’ end of the gene, was transferred to the nuclear genome in the ancestor of core eudicots, most likely preceding the generation of a stop codon TAA. That left mitochondrial rpl2 shortened to form present in most angiosperm mitogenomes [45]. Following fission, the 5’ part of rpl2 was transferred to the nucleus in Fables, lettuce and probably several other eudicot lineages [46]. The mitochondrial rps2 of P. patens is longer than in core eudicots and corresponds to structure found in the early angiosperms: Magnoliales [47] and Nymphaeales [48]. The lack of variant in the rpl2 and surrounding regions suggests, that process if rpl2 transfer to nuclear genome appear after divergence of Ranunculales, but before divergence of Proteales, since mitochondrial genome of Nelumbo nucifera contains shortened rpl2 known from core eudicots [17].
Phylomitogenomics relationships of Pulsatilla patens
The phylogenomic analyses based on protein-coding mitochondrial genes resulted in trees with similar topologies (Figure 3, Figure S1) and obtained topologies are congruent with previously published papers [49].
The phylogenetics analysis of 39 protein-coding mitochondrial genes reveals Aconitum, Anemone and Pulsatilla (Ranunculales) as early divergent eudicots with Nelumbo nucifera (Proteales) as it’s closest relative (Figure 3). However, the phylogenomics analyses did not resolve Proteales-Ranunculales relationships with a maximal support. Previous studies using mitochondrial dataset [49], resolved Nelumbo as a basal for eudicots, however, it naturally did not include any member of Ranunculales. The phylogenetic position of Pulsatilla patens is congruent with plastome datasets, which resolved Ranunculales as an earlier divergent than Proteales [50]. Also, the analysis of 1,594 nuclear loci resolved Ranunculales as an older group than Proteales [51] and provided results congruent with chloroplast and mitochondrial data. The main aim of this analysis was the phylogenetic placement of mitogenomes containing polymerase genes (DNApol, RNApol) which are scattered along the whole tree.
Phylogenetic relationships of mitochondrial DNApol and RNApol genes in many cases did not reflect the phylogenetic position of species (Figures S2-S3), suggesting it’s independent acquisition via horizontal gene transfer, which was reported several time in plant mitochondrial genomes [4, 52, 53]. The sequences of polymerase genes extracted from mitogenomes of close related species usually group together, but in some clades even orders are mixed. The phylogenetic position of Pp_RNA_pol gene remained unresolved (Figure S3) and Pp_DNA_pol gene was clustered together with genes coming from evolutionary unrelated genera like Theobroma, Actinidia and Prunus (Figure S2).
Repeats
Based on the currently known genomic resources, these large repeats seem to be limited to the genus Pulsatilla. The BLAST against known mitogenomes revealed only presence of very short (up to 400 bp in case of Ligustrum sp.) parts of these repeats in other plant species. Each of three P. patens mitochondrial chromosomes is flanked with the same pair of inverted repeats, which seem to be unique among known mitogenomes. Chromosomes of Solanum tubersoum are also flanked with repeats, however in this case, the repeats were present in only two copies and their length ranged from 1,208 to 11,915 bp [12].
Beside mentioned large repeats, which shape the P. patens mitogenome structure only two repeats longer than 1,000 bp were found. The repeats S1 (1,889 bp) and S2 (1,594 bp) were located at chMt2 and present in two copies. Additional two repeats shared among chromosomes fall within 500-1,000 bp range: repeats S3 (848 bp long) were found at chMt1 and chMt2 and repeat S4 (697 bp) were found at chMt2 and chMt3.
Chloroplast derived regions
DNA transfer from plastome to mitogenome is well documented in almost all vascular plant lineages [29, 54]. This transfer is usually not limited to a single gene, but in most cases it comprises a cluster of genes, which seems to be conserved among angiosperms [55]. The Pulsatilla patens, which belongs to early divergent angiosperm order Ranunculales, reveals over 35 kbp of plastome derived regions in its mitogenome, called MTPT (plastid-derived mitochondrial DNA), which comprise 3.6% of total mitogenome length (Table S1). The pt-derived DNA is distributed among all the mt-chromosomes: 5, 3, and 1 transfer were found in chMt1, chMt2 and chMt3 respectively. Out of 9 unique cp-mtDNA regions, 5 were transferred from chloroplast inverted repeat regions, including longest, over 19 kbp. In total, almost the whole IR region (except gene ycf1) was transferred to mitochondrion in scattered form mainly over chMt1 with one part located on chMt2. The length of the longest MTPT is unusual, since the most previously identified chloroplast derived regions fall within 200-4,000 bp range [56] and the largest, 12.6 kbp long MTPT was found in Zea mays mitogenome [57]. In the other known mitogenome of Ranunculales, Hepatica maxima, the MTPT fragments are up tp 7.1 kbp long, with total length of 16 kbp [16] The presence of an almost complete IR region in P. patens mitogenome suggests that it was derived from single transfer that was subsequently split and partially duplicated by rearrangement. Similar scenario was described in Silene conica, where 35 kbp region was transferred in a single event but later, most of this MTPT was removed by a series of deletion events [57]. The mitogenome of P. patens contains 26.7% of plastome (excluding one copy of IR) which falls between Cucumis melo (22.7%) and Bambusa oldhamii (40.9%) Due to lack of mitogenomic data on early divergent eudicots it’s difficult to infer the evolutionary significance of the amount of transferred cpDNA.
Among MTPT gene clusters identified by Wang et al. [55], three were not found in P. patens mitogenome: atpE-rbcL, psbE-psbF and ycf1-trnN. Another transferred cluster specific for eudicots, psaA-psaB [55] was splitted among chromosomes chMt1 and chMt2 (Table S1). The rpoC1-rpoC2 cluster found at chromosome chMt3 was not specific for Pulsatilla, as blasting it against known plant mitogenomes confirmed its presence in six mitogenomes belonging to different orders: Phoenix dactylifera (Lililes), Sapria himalayana (Malpighiales), Solanum aethopicum (Solanales), Ziziphus jujuba (Rosales), Spondias tuberosa (Sapindales) and Ammopiptanthus mongolicus (Fabales). Recent studies hypothesize that plastome sequences were initially acquired by intracellular gene transfer and then were transferred among plant lineages via mitochondrial horizontal gene transfer [56], however in the case of Pulsatilla mitogenome this scenario in rather unlikely. The longest MTPT (chloroplast transfer 4) contains 11 genes which order is characteristic for all known Pulsatilla chloroplast genomes and intergenic spacer sequences is specific for this genus [32, 33, 58, 59]. Separation of MTPT from plastome reads with the whole genome sequencing approach could be challenging, especially while only short reads sequencing results are available, which could be mapped on both regions making them look like each other. Application of two long reads sequencing platform enabled presence of core mtDNA and MTPT in single reads, which confirmed the mitochondrial origin of MTPT which fall within the range of cpDNA variation at intraspecific level [32] suggesting it’s acquisition via intracellular gene transfer.
RNA editing
According to previous studies RNA editing is obligatory for few mitochondrial genes by creating initiation or termination codons [12]. In the Pulsatilla patens mitochondrial genome, the C->U editing is required to create start codon in nad1 and stop codons in atp6, rps11, ccmFC, and cox2 genes.
The mitochondria of flowering plants usually have 300–500 sites that are subjected to RNA editing [60]. Despite considering only editing events with frequency above 0.5, the number of affected sites (907) is bigger than in Solanum [12] and Populus [61], where 799 and 355 editing events were detected respectively. High number of editing sites in CDS of P. patens mitogenome can not be explained by tissue used for RNA extraction, since leaf (as in case of our study) and root tissues are considered as less frequently edited than flower tissues [62]. The lower number of predicted RNA-editing sites was given for Hepatica maxima [16] due to presence of two genes missing in P. patens (rpl10 and sdh4). However, predicting-based estimation of RNA-editing sites usually overstates the amount of RNA-seq verified editing sites [63].
The number of observed edited sites varies from single in rpl2 up to 65 in nad4 (Table 1, Fig. S4). Considering gene category, the highest number of edited base pairs was found in the genes of Complex 1 (nad) and Cytochrome C biogenesis (ccm), which is congruent with data obtained for Populus [61] as well as for early land plants [63]. Despite deep transcriptome sequencing, the reverse U->C editing, which was documented in many early divergent tracheophytes lineages [64], has not been observed at assumed frequency in P. patens mitogenome.
In our study, the high number of editing sites can not be explained by possible false-positive results often observed in many surveys, where DNA and RNA libraries were prepared from different individuals.
Interspecific variation in protein- coding genes
The comparative analyses of P. alpina, P. patens and P. pratensis protein coding genes (CDS) revealed 50 SNPs (43 dN and 7 dS) and 23 (18 dN and 5 dS) differentiating P. patens from P. alpina and P. pratensis respectively (Table XS2). Eight out of 38 compared genes did not reveal any interspecific variation (nad1, nad3, nad6, nad9, cox3, ccmC, rpl2, mttB) and in the next seven genes only single SNPs were detected (nad4, nad4L, nad5, nad7, cox1, ccmB, rps12, matR). The interspecific comparison revealed that the genes of complexes I-IV are least variable, while the most of interspecific variation is accumulated in the genes of Complex V and genes of small ribosomal subunit. Among the genes of Complex V the lowest variation was found in atp9 (only two synonymous mutations), while in the remaining number of nucleotide substitution among atp genes varied from 7 (atp8) to 12 (atp1).
The genes belonging to other groups revealed low variation, except from ccmFC where 16 substitutions were found.
Most mutations found in protein- coding genes were substitutions, however the indels were also identified in five genes including 2aa indel in atp4, 5aa in ccmFC, 3aa in rps1, 4aa indel in rps3 and the largest, 30aa indel in rps4. With single exception of atp4 (3 bp shorter CDS in P. pratensis than P. patens), the rest indels were found in P. alpina.
The mitochondrial genes are usually considered as slower evolving than chloroplast genes [65, 66], but the comparative analyses at genus level are limited to only a few genera.
The number of mutations of Pulsatilla mitochondrial CDS did not deviate from chloroplast CDS, especially when comparing atp and rps genes. Chloroplasts CDS of atpA, atpB, atpE and atpI revealed 3, 2, 1 and 1 substitutions respectively despite analysing six individuals from three Pulsatilla species [32]. Also the genes of small ribosomal subunit of chloroplast genome were less variable, with only two substitutions in rps11 and rps15. With exception of the most variable ycf1 (39 substitutions), the numbers of substitutions in the remaining, relatively variable chloroplast CDS, ranged from five to seven [32].
The raw numbers of indels and substitutions could be biased by differences in genes length, therefore the pi nucleotide diversity is often used to estimate relative differences in the variation among genes. The pi values in the case of analysed Pulsatilla species (excluding mentioned above non-variable genes) varied from 0.0003 (matR, nad5) to 0.0516 for rps4 (Table S3). Nine of the ten most variable genes comprised genes from Complex V and small ribosomal subunit (Table S4), which corresponded to the raw numbers of detected SNPs.
The mitochondrial protein- coding genes of Pulsatilla genus revealed relatively high variation in comparison to other vascular plants. Only eight SNPs were found among mitochondrial genes in the genus Larix [66]. On the other hand, the number of SNPs at species level in Silene vulgaris was much higher (144 SNPs), but they mostly concentrated into three genes: atp1, atp6 and cox1 [10]. However, the genus Silene is known for its extreme mitochondrial diversity, not only at substitution level but also in structure and gene content [10].
Mitochondrial genomes as a resource for phylogenetic and phylogeography studies
The plastomes of Ranunculaceae are structurally variable [32, 68], but at the genus level the gene order seems to be conserved [32, 33, 58, 59]. Currently the data on sequenced mitogenomes of Ranunculaceae are restricted to three species, with different mitogenome structure, gene content and order. However, based on preliminary mappings of sequenced short-read libraries from Pulsatilla alpina, P. vernalis and P. pratensis, the selection of regions for phylogenetic and phylogeographic studies on generic level is quite challenging due to common presence of large repeats and MTPTs. Moreover, the linearity of preferred regions should be conserved among species. Eight single copy linear regions of minimum length 10 kbp with shared, identical gene order among analysed species were identified at chromosome chMt1 at 89K-122K, 149K-160K, 168K-179K, 123K-252K, 283K-299K, 303K-333K, 377K-399K, 409K-424K. Together 165 kbp of the total 452 kbp (36 %) of chromosome chMt1 were suitable for low-level phylogenetic and phylogeographic studies. In the case of shorter chromosome chMt2, for regions meeting above criteria was found: 83K-94K, 101K- 126K, 214K-276K, 295K-323K, which consist 37% (126K bp) of total chMt2 length. At chromosome chMt3 only two phylogeny suitable regions were found between 49K-67K and 80K-100K, which consisted 19% (38k bp) of total chMt3 length. The phylogenetic usage of mitogenomes in plant studies is mainly limited to CDS sequences, due to high evolutionary dynamics of molecule structure. The beyond-CDS mitogenomic phylogenetics is limited mainly to early land plants with a stable structure like liverworts [65] and mosses [6, 69, 70] and few studies on vascular plant taxa where both CDS an intron sequences were used [71]. The amount of shared mitogenomic clusters among Pulsatilla species seems to be high in comparison to genera like Picea, where only two shared clusters up to 9 kbp were found [26]. Mitochondrial conserved regions seem to have a great potential in phylogenetic studies, despite their taxa specific restrictions. Since their application in low- level phylogenetics is not straightforward like plastid genomes, they can provide better phylogenetic resolution than plastomes [72].