Сhloroplast genomes (plastomes) of land plants are generally conserved in their size, gene content and order. The evolutionary origin of chloroplasts, as well as mitochondria, traces back to ancient endosymbiotic bacteria that consequently greatly reduced the genome, and contain only a small proportion of the ancestor’s genes [1, 2].

Plant plastomеs possess a quadripartite structure composed of large single-copy (LSC) and small single-copy (SSC) regions divided by two parts of inverted repeat (IR) [2, 3]. Plastome size of higher plants is usually around 150,000 bp in length and comprise approximately 120–130 genes, among which about 75 genes encode proteins of photosystem I and II, as well as for other proteins, involved in photosynthesis [see, for example, [4]], while other genes encode ribosomal RNA and proteins and transfer RNA. The highest deviations of the gene order and content in plastomes of land plants have been reported in non-photosynthetic species, in which extremely reduced plastomes were found - up to 11 Kbp [5]. However, in photosynthetic plants lineage-specific gene losses, as well as translocations that change the gene order were also observed [6–8].

IRs are usually regarded as the most stable part of the plastome. Indeed, the substitution rate of the IR sequences is the lowest compared to single copy regions [3]; the plastomes where IRs are absent (e.g. IRLC clade of Fabaceae) exhibit elevated substitution rates [9] and vice versa, genes that were translocated from SC to IR slow down their substitution rate [10]. IRs typically range in size from 15 to 30 kbp and contain a core set of genes consisting of four rRNA genes (4.5S, 5S, 16S and 23S rRNA), tRNA genes (trnAUG, trnI-GAU, trnN-GUU, trnR-ACG and trnV-GAC). The IRs of many land plants also contain a number of other genes as a result of lineage-specific expansions and contractions [3]. Few evolutionary lineages demonstrate large-scale expansions (exceeding several kbp and containing numerous genes) of the IR [3, 11]. In particular, several overlapping inversions affect size, gene content and order of the IR in leptosporangiate ferns, the clade that includes most fern species [8, 10].

Another evolutionary event that can affect plastome structure is horizontal gene transfer (HGT). HGT between nucleus, mitochondria and plastids has been shown to occur with a high rate and contributed significantly to the plant genome evolution by relocating and refashioning of the genes and consequently contributing to genetic diversity. Transfers of DNA fragments from the mitochondria or plastids to the nucleus are the most common reported ones [12–16]. The mitochondrial genomes are also often invaded by plastome-derived sequences; the presence of DNA from nuclear genomes also has been shown in a number of lineages of flowering plants [17–24] and ferns [25]. Translocations of mtDNA fragments to plastid genome are much rarer. Currently, only two cases, all from flowering plants, are known. In Daucus carota (order Apiales, Apiaceae) plastome a ~1.5 Kb region with high similarity to Vitis vinifera (order Vitales) mitochondrial sequence was found [26]. This region did not contain any typical plastid genes. Characterization of Daucus mitogenome and screening of plastomes of other Apiaceae suggest that it was inserted from the mitochondrial genome to the plastome in a common ancestor of the genus Daucus [27]. Another example is the horizontal gene transfer is a 2.4-kb segment of mitochondrial DNA into the rps2–rpoC2 intergenic spacer of the plastome of Asclepias syriaca (Apocynaceae) [28]. Thus, unlike the mitochondrial genomes, which are affected by insertions of plastid and nuclear sequences, the plastomes of flowering plants are infrequently profited by DNA transfer from the other cell compartments [22, 29, 30].

Extant ferns are non-seed vascular plants for which 45 families are currently known (with approx. 280 genera), of which more than half belong to the order Polypodiales, in line with the classification of [31]. The majority of Polypodiales species fall into two large sister clades - Eupolypods I and Eupolypods II, and the remaining to families Pteridaceae, Dennstaedtiaceae, Saccolomataceae, Lindsaeaceae (basal clade) [32]. Wolf et al. [8] for the first time examined structure of the plastome across few widely ranged representative fern taxa but no one has represented the Polypodiales. Zhu et al. [3] analyzed the evolutionary rate and shifts in IR boundaries of land plants, including seven ferns but also no Polypodiales were considered. Raman et al. [33] characterized Cyrtomium falcatum plastid genome and found some differences with congeneric species C. devexiscapulae in tRNA gene content and start codons. For 24 fern samples, five of which were Polypodiales species (other represented ten extant orders), the part of LSC region (rpoB-psbZ) was analyzed and considerable genomic changes for distant species (belonging to different orders) were found [34]. A comparison of a few taxon-wide fern plastomes (of Lycopodiophyta, Psilotopsida, Equisetopsida, Marattiopsida and Polypodiopsida) showed that some lineages have experienced multiple IR changes including expansions and inversions while others demonstrated the stasis [35]. Recently many new sequences of fern plastomes (including Polypodiales species) were released. However, the corresponding study reports only the results of phylogenetic analysis of these sequences, without detailed analysis of their gene content and structure [36].

This clearly indicates that the diversity of plastome structures in ferns is insufficiently explored. With this premise, we characterized four additional plastome sequences from Polypodiales, three from Dryopteris and one from Adianthum, and performed comparative analysis of all available fern plastomes.

We sequenced and assembled de novo new plastome sequences for three Eupolypods I species: Dryopteris filix-mas, Dryopteris blanfordii and Dryopteris villarii. The plastomes have typical quadripartite structure, are similar in their size (148,568, 152,945 and 148,727 bp, respectively). A total of 130 genes were annotated by DOGMA for D. filix-mas, including 91 protein-coding genes (5 of them are duplicated in IRs), 25 tRNA genes (5 of them are duplicated in IRs) and 4–5 rRNA genes (all of them are duplicated in IRs). Plastomes of Dryopteris species, including D. blanfordii, D. villarii, D. filix-mas and previously reported D. decipiens were identical to each other in gene content of LSC and SSC, but differed in IRs (Fig. 1).

We also sequenced one the plastome of Adianthum hispidulum, from Pteridaceae, the basal group relative to eupolypods. The complete plastome sequence of A. hispidulum was 151,327 bp in length, consisted of an LSC (83,188 bp), SSC (21,459 bp), and IRs (23,340 bp). The gene content in the plastome of A. hispidulum slightly differs from that of congeneric species A. capillus-veneris (published in [37]) – in A. hispidulum trnT-UGU gene (located in IR) is completely absent while in A. capillus-veneris it is represented by a pseudogene (Fig. 2).

Comparative dataset of IR sequences which includes all plastomes for Polypodiales available in the public sequences databases and our new data comprises 45 sequences: 31 belong to Eupolypods I, 9 – to Eupolypods II and the remaining 5 - to Pteridaceae or Dennstaedtiaceae. We have made re-annotation of published Polypodiales plastomes. The IR/SSC border of all species lies within the ndhF gene and IR/LSC - within ndhB gene. The similar IR borders were previously defined for many groups of monilophytes (Psilotales, Ophioglossales, Equisetales, Marattiales and Polypodiopsida) [35]. Commonly, chloroplast genomes of Polypodiales carry sequences (~ 246 bp) with high similarity to ycf68 ORF in the IR regions, though they were not previously annotated by the authors. The ycf68 is a putatively functional gene located in the trnI-GAU intron, which is present in many land plant chloroplast genomes but often is not annotated as its function is still unknown [38, 39]. Another putative gene, ORF42, was annotated in the trnI-GAU intron for the all species included in the analysis. ORF42 was found previously in the intron region of trnA-UGC plastid gene of some flowering plant species, for example, Veratrum patulum O. Loes. (Melanthiaceae) [40] and Pelargonium × hortorum L. H. Bailey (Geraniaceae) [11]. The sequence with high similarity to ORF42 was found in the mitochondrial genome of Phaseolus; presumably as a result of plastid-to-mitochondrion lateral gene transfer [19].

The IR structures of Polypodiales are shown in Fig. 2 and Additional file 1. The number of genes normally varies from 14 to 16. The plastomes of Polypodiales mainly accumulated the gene number variability within IR in two areas: from 0 to ~3 Kbp and from ~7 to ~11,5 Kbp regions of IR. These regions also demonstrate lower sequence similarity. The other two regions - from ~3 to ~7 Kbp and from ~11,5 Kbp till the end of IR - were largely conservative in gene number and sequence.

The 0 to ~3 Kbp region contains tRNA (pseudo)gene (Fig. 2). Though it is annotated as functional two-exon trnT-UGU gene in many species (e.g. 33, 37), the tRNA structure prediction with tRNAscan-SE does not support its functionality. In other species it is annotated as pseudogene (trnT-UGU or trnL-CAA) or completely missing. Manual check however shows that pseudogene is present in all species analyzed (Fig. 2) except for Adianthum hispidulum where it was completely lost.

In the ~7 to ~11,5 Kbp variable region three genes show partial or full deletions or duplications (trnI-GAU, ycf68, rrn16) (Fig. 2) (Additional file 1). Some species of Polypodiales have partially or completely lost trnI-GAU and ycf68 (located in intron of trnI-GAU) namely: O. sensibilis (Onocleaceae) and three of four Dryopteris species (D. blanfordii, D. villarii, D. filix-mas, Dryopteridaceae). In D. filix-mas trnI-GAU and ycf68 are completely lost. D. blanfordii has partially lost the trnI-GAU gene (the intron and one of exons were deleted), but has a duplication of the large part of rrn16 gene. D. villarii has lost intron and one exon of trnI-GAU but no rrn16 duplication. As result of deletions/duplications the IR size of Dryopteris species varied: IR of D. filix-mas and D. villarii were ~1450–1570 bp shorter but D. blanfordii, on the contrary, had IRs 647 bp longer (Cyrtomium species were used as reference). Surprisingly D. decipiens, reported in [36], has no deletions in this region. O. sensibilis (Onocleaceae) plastome has a deletion of ycf68-trnI-GAU region. Onoclea and Dryopteris species belong to different clades: Eupolypod II and Eupolypod I. Therefore, the losses of ycf68-trnI-GAU regions are independent events.

Due to the deletions and duplications IR size in Polypodiales varies from 22 Kbp (Ceratopteris richardii [KM052729]) to 26,9 Kbp (Cystopteris protrusa, [KP136830]. In particular, we found the most of the large insertions (370 bp and more) in the highly variable intergenic spacers mentioned above: intergene 14 (between rrn16 and rps12), intergene 16 (between rps7 and psbA) and intergene 19 (between ycf2 and trnN-GUU), see Fig. 3. For the insertion’s length, localization and similarity to known high plant sequences see Table 1.

Region	All Polypodiales	Eupolypods 1	Eupolypods 2	Dennstaedtiaceae + Pteridium
Intergene1	1,27	1,43	1,22	1,65
ndhB	2,46	3,50	2,78	2,83
Intergene2	11,52	12,55	10,86	14,59
trnT-UGU	16,59	26,18	3,11	31,63
Intergene3	2,68	3,32	2,20	3,45
trnL-CAA	2,20	3,24	2,44	2,53
Intergene3p	5,43	8,53	6,47	5,13
trnR-ACG	85,25	0,00	138,59	0,00
Intergene4	10,11	13,42	12,65	10,99
rrn5	0,15	0,22	0,27	0,11
Intergene5	7,25	0,46	8,02	7,12
Intergene6	2,35	2,49	3,46	0,93
rrn23	0,78	9,81	1,77	0,30
Intergene7	2,97	2,51	4,30	2,26
trnA-UGC	0,00	0,00	0,00	0,00
Intergene8	22,93	10,48	59,98	11,34
orf42	530,31	1221,25	3,38	3,77
Intergene9	0,74	1,05	0,86	0,73
trnA-UGC2	0,02	0,04	0,00	0,00
Intergene10	2,67	16,92	1,71	1,70
trnI-GAU	0,00	0,00	0,00	0,00
Intergene11	2,07	3,40	2,61	2,22
ycf68	2,87	2,55	8,39	1,07
Intergene12	68,77	41,48	245,28	15,63
trnI-GAU2	2,31	5,57	0,00	0,00
Intergene13	3,01	4,05	3,76	3,10
rrn16	0,42	0,18	0,97	0,20
Intergene14	9,37	27,47	11,23	8,26
rps12	1,88	2,13	2,77	1,44
Intergene15	1,43	1,48	1,24	1,59
rps7	2,20	2,32	2,65	1,83
Intergene16	7,31	7,69	7,79	7,57
psbA	1,42	1,89	1,78	1,30
Intergene17	2,57	4,21	2,63	2,91
trnH-GUG	0,35	0,86	0,00	0,01
Intergene18	3,14	1,92	3,36	4,20
Intergene19	5,47	5,07	4,97	6,56
trnN-GUU	0,00	0,00	0,00	0,00
Intergene20	3,02	5,42	3,27	4,03

The bolded rows marked intergenic spacers there LRT for the branch site model A was significantly better than the model B. (Model A was used for the ribosomal RNA and non-coding regions and supposes that substitutions in a DNA sequences of different tree branches are homogenous, model B – the same but supposes non-homogenous substitutions, see also materials and methods)

An unusual 1663 bp insertion was found rrn16 and rps12 genes of Woodwardia unigemmata (coordinates 93,404...95067). W. unigemmata is a fern of family Blechnaceae (Eupolypods II), whose plastome was sequenced by Lu et al., 2015 [41], and no genes were annotated in this region previously. The large part of W. unigemmata insertion has high sequence similarity to the insertion of Lepisorus clathratus (Polypodiaceae, Eupolipods I) located in the same region of IR, between rrn16 and rps12 genes. About 1160 bp of W. unigemmata and L. clathratus insertions demonstrated 65–78% similarity to each other. No similar sequences were found in IR of other Polypodiales, except for small (about 150 bp) sequence in Matteuccia struthiopteris (Onocleaceae, Eupolypods II) which has 62–66% similarity to L. clathratus and W. unigemmata insertions (further called WL-sequences).

Surprisingly the WL-like sequences were found in LSC regions of species from distant taxonomic groups of ferns, outside Polypodiales (Table 2). Firstly, small part of WL sequence has similarity to the fragment of LSC plastome of Plagiogyria species (Plagiogyria is a single genus in monotypic family Plagiogyriaceae, Cyatheales, see for example «Flora of China» [42]). To be more precise, the 772 bp part of WL sequence has 67% identity to the region of the P. glauca plastome (coordinates 29,128...29895, KP136831) and to P. japonica plastome (partial sequence, coordinates 4503–5273, HQ658099). In both species sequences with homology to WL lie within trnD-GUC-psbM intergenic spacer. Other ferns that also contain WL-like sequences in LSC are Ophioglossum californicum (KC117178) [34] and Mankyua chejuensis (KP205433). Both these species belong to Ophioglossales (basal ferns); they are distant from Polypodiales (core leptosporangiate ferns) and from Cyatheales. The fragments with similarity to WL sequences in Ophioglossales ferns have about 255–275 bp length and are located in the different regions of LSC. In O. californicum (KC117178) it is found in intergenic spacer between the trnT-GGU and trnfM-CAU genes, in M. chejuensis (KP205433, JF343520) it is also located in LSC region but in the other intergenic spacer (between trnL-UAA and rps4). It is annotated as ORF295 (Fig. 4). Altogether, this suggests that there is a translocation of DNA fragments from LSC to IR (or vice versa) during evolution of fern plastomes.

№	Clades	Species	Insertion IR coordinates		Insertion size	Location	Genbank accession number	Description
			start	end
1	1Eu*	Dryopteris blanfordii	7615	10,000	2385	Large insertion between trnI-GAU and rrn16	This work	Contains duplicated rrn16 and trnF-GAA genes
2	1Eu	Polypodium glycyrrhiza	9546	11,764	2218	intergene14 ( rrn16-rps12)	KP136832	~ 700 bp part of insertion has homology to insertion of L. clathratus
3	1Eu	Lepisorus clathratus	10,863	14,073	3210	intergene14 (rrn16-rps12)	KY419704	Part if insertion has homology to part of WL and another to insertion of P. glycyrrhiza
4	1Eu	Cyrtomium devexiscapulae	22,596	22,968	372	intergene19 (ycf2-trnN-GUU)	KT599100	No significant homology to known sequences was found
5	1Eu	Cyrtomium falcatum	22,570	22,942	372	intergene19 (ycf2-trnN-GUU)	KP189363	No significant homology to known sequences was found
6	2Eu	Cystopteris protrusa	6100	6286	186	rrn23	KP136830	Duplication inside IR (homologous to 6403–6560)
7	2Eu	Cystopteris protrusa	6897	9498	2631	rrn23	KP136830	~ 950 bp fragment of insertion has 75–76% homology to mitochondrion of different mosses
8	2Eu	Asplenium pekinense	10,447	12,688	2316	intergene14 (rrn16-rps12)	KY427331	~ 200 bp of insertion has homology to rrn16-rps12 spacer of L. clathratus
9	2Eu	Rhachidosorus consimilis	10,780	11,683	903	intergene14 (rrn16-rps12)	KY427356	No significant homology to known sequences was found
10	2Eu	Woodwardia unigemmata	11,015	12,678	1663	intergene14 (rrn16-rps12)	KT599101	WL sequence
11	2Eu	Matteuccia struthiopteris	13,885	14,255	370	intergene16 (rps7-psbA)	KY427353	Part if insertion has homology to part of WL sequence
12	2Eu	Asplenium prolongatum	22,083	23,414	1331	intergene19 (ycf2-trnN-GUU)	KY427332	The most of insertion has homology to part of rrn16-rps12 spacer of Polypodium glycyrrhiza
13	Pt	Ceratopteris richardii	14,185	14,499	314	intergene18 (trnH-ycf2)	KM052729	No significant homology to known sequences was found
14	Pt	Cheilanthes lindheimeri	22,714	25,096	2382	intergene19 (ycf2-trnN-GUU)	HM778032	No significant homology to known sequences was found

^* 1Eu Eupolypods I, 2Eu Eupolypods II, Pt Pteridaceae

Interestingly, the part of the WL sequence has high identity (67%) to the 203-bp fragment of about 5,9 Kb intergenic spacer between tRNA-CGA and tRNA-TTT genes of mitochondrial genome of Asplenium nidus (partial sequence, coordinates 3188…2986, AM600641) (Table 2). Asplenium is a genus that belongs to the same clade of Polypodiales - eupolypods II as Woodwardia but to the other family - Aspleniaceae. We performed similarity search of WL-sequence against fern mitochondrial contigs available in Utah State University Repository [45; http://digitalcommons.usu.edu/fern_genome/]. Two hits were found in Plagiogyria formosana - 244 bp (identity 72%, contig №439, coordinates 5336–5585) and 277 bp (identity 72%, contig №93, coordinates 1–277). No hits were detected in the contigs of another available fern species (D. conjugata, Pteridium aquilinum, C. richardii, Polypodium glycyrrhiza, C. protrusa).

In contrast to previous observations on the stability of the IR region, we found high variability in IR sequence and gene content in Polypodiales ferns. There are two hypervariable regions – one located at the beginning of IR, 0–3 Kb, and the second is 7–11 Kb region. These regions are the subject to the similar evolutionary changes occurred independently in the different clades. The first region in most species contains tRNA pseudogene. The members of both Eupolypods I and Eupolypods II demonstrated independent deletion of trnI-GAU-ycf68 region (i.e. Dryopteris and Onoclea). Polypodiales (together with Salviniales and Cyatheales) belong to a clade called core leptosporangiates [43]. Their plastomes acutely differ from those of eusporangiate ferns (Psilotales or whisk ferns, Ophioglossales, Marattiales, Equisetales) [8, 33]. It should be noted however that comparative analysis of fern plastomes is obfuscated by the uncertainty of the annotation of tRNA genes. This concerns, in particular, the intron-containing trnT-UGU, which was reported in IR (between ndhB and trnR-ACG) of several fern plastomes [37, 43, 44] and thought to be specific feature of core leptosporangiates. But intron-containing trnT-UGU was not found in other fern lineages [45, 46] or in any other plants outside ferns; only intronless trnT-UGU is present. This is unusual given that plastid tRNA genes, in contrast to protein-coding genes, have highly conserved exon-intron structure. Gao and coworkers supposed that tRNA genes may be lost repeatedly independently during evolution of ferns and probably the loss of trnT-UGU is the one of those events [45]. Our analysis which included manual reannotation and check using tRNA prediction program tRNAscan-SE however does not support the functionality of intron-containing IR-located trnT-UGU in any Polypodiales species where it was reported. In contrast, we found the trnT-UGU pseudogene in the IR of almost all Polypodiales. Most likely, this pseudogene is difficult to be recognized and therefore results of automatic annotation could be interpreted as gene loss. Notably, Gao and co-workers [45] compared the sequences of putative intron-containing and intronless trnT-UGU and it can be seen that the former are unusually divergent, much higher than expected for a functional tRNA gene. We conclude that intron-containing IR-located trnT-UGU is an artifact caused by the shortcomings of the automatic annotation. Moreover, two parts (“exons”) of this pseudogenes can be recognized by automation annotation programs, such as DOGMA, as different tRNAs - trnT-UGU and trnL-CAA.

In the second hypervariable region, 7–11 Kb, we found an unusual insertion (the WL-sequence) in two unrelated Polypodiales – Woodwartia unigemmata and Lepisorus clathrathus. Smaller insertion with high similarity to the WL-sequence was found in the same region in Mattheucia sthruttiopteris. In addition, the insertions with high similarity to WL-sequence were found in plastomes of Ophioglossales (basal ferns, distant from Polypodiales) but in different position – in LSC region.

The WL-sequence has high similarity with the region of mitochondrial genome of Asplenium nidus. This has two possible explanations: that it is either the sequence of mitochondrial origin, which was integrated in the plastome, or the sequence of plastid origin, which was integrated into mitochondrial genome and lost from the plastid genomes of most ferns, with exception of W. unigemmata and L. clathrathus. By now we can’t make a conclusion about which of these two hypotheses is true, due to the unavailability of fern complete mitochondrial genome sequences.

In any case, it is a result of the horizontal genome fragment transfer between mitochondria and plastids. Horizontal genome fragment transfer is a phenomenon, commonly observed in the pro- and eukaryotes. In plants, the presence of three genomes within a cell compartments (mitochondria, chloroplast and nucleus) leads to different possible types of intracellular genome fragments exchange: between organelles and nucleus and between mitochondria and chloroplasts, bidirectional [47]. The transfer of genetic material from organelles to the nucleus seems to be a continuing evolutionary process of the prokaryotic ancestors’ genome reduction [48, 49]. Many reports asserted that plant mitochondrial genomes are unusually prone to the introgression of alien sequences compared to chloroplast and nuclear genomes [47, 50]. There are only few data on mitochondria of ferns. No complete mitochondrial genome assemblies are available, only contigs. Multiple regions with strong sequence similarity to plastid DNA were detected by [51] but they didn’t relate to the plastome sequences in the total genomic contigs of six ferns species Dipteris conjugata (Gleicheniales), Plagiogyria formosana (Cyatheales), Pteridium aquilinum (Dennstaedtiaceae), Ceratopteris richardii (Pteridaceae), Polypodium glycyrrhiza (eupolypods) and Cystopteris protrusa (eupolypods). Authors speculated that the plastome-like sequences reside within the nuclear or mitochondrial genomes [42]. Assuming this is the case, it implies that the horizontal transfer of organelle genome fragments are not rare events in the evolution of ferns.

In this study we investigated the structure and evolutionary stability of IRs of plastomes in Polypodiales ferns. The two regions of IRs were found to be highly variable: (i) the sequences between ndhB and trnR-ACG genes (~3 Kbp) and (ii) the fragment including the rrn16 gene and flanking vicinity regions (~4,5 Kbp). This blinking of trnI-CAU, trnT-UGU, ndhB and rps12, trnI-GAU, ycf68, rrn16 genes related to these regions was observed in different species. The plastomes of three Dryopteris species demonstrate dynamic process of trnI-GAU elimination/rrn16 duplication.

Two Polypodiales species - W. unigemmata and L. clathratus - have an unusual sequence in the IR region. It demonstrates similarity to LSC spacers trnL-rps4 of Ophioglossales and pbsM-trnD of Cyatheales and with the part of mitochondrial genome of Asplenium (Polypodiales). We suppose these features are a consequence of intraplastomic rearrangements as well as of the transfer between the chloroplast and mitochondrial genomes during the evolution of ferns.