Additional file 1: Fig. S1. Representative cropped agarose gels of CR1/CR2 amplicons obtained for Struthio camelus (a), Rhea pennata (b), Rhea americana (c) and Crypturellus tataupa (d). The numbering of amplicons separated in the agarose gel corresponds to the reaction numbered in Table S1 in Additional file 2. The numbers of PCR fragments, which were finally sequenced, are underlined. Lane M - GeneRuler 1 kb DNA Ladder (Thermo Scientific). Original uncropped and unprocessed gels are presented in Fig. S13a-d. Fig. S2. Sequence comparison of two copies of tRNA-Pro, ND6from, tRNA-Glu genes and control region, found in this study as duplicated, with appropriate sequences of mitogenomes previously deposited in GenBank and showing a typical avian gene order. The alignments are shown for Struthio camelus (a), Rhea americana (b), Rhea pennata (c) and Crypturellus tataupa (d). Dots indicate residues identical in the obtained copies with those in the single sequence previously published. Genes for tRNA-Pro are marked in yellow coloring, NADH dehydrogenase subunit 6 in cyan coloring, tRNA-Glu in grey shading, and control regions in red fonts. In the case of Crypturellus tataupa species only two obtained copies are compared to each other because the previously published mitogenome of these taxa is incomplete and it does not contain any of the analyzed genes or the control region. The sequences were aligned with MUSCLE [125] in MEGA [124]. Fig. S3. Sequences of mitogenomic fragments amplified and sequenced in this study for Struthio camelus (MH264503), Rhea americana (MK696563), Rhea pennata (MK696564) and Crypturellus tataupa (MK696562). Genes for tRNA-Leu are marked in bold, NADH dehydrogenase subunit 5 in magenta coloring, cytochrome b in red coloring, tRNA-Thr in green coloring, tRNA-Pro in yellow coloring, NADH dehydrogenase subunit 6 in cyan coloring, tRNA-Glu in grey shading, 12S rRNA in blue fonts, tRNA-Val in brown fonts, 16S rRNA in green fonts, tRNA-Phe are in bold and underlined, and control regions in red fonts. Fig. S4. Cropped agarose gels of potential NUMTs amplified for Dromaius novaehollandiae (a) and Casuarius casuarius (b) with the use of primers designed for detection of CR1/CR2 fragments. The numbering of lanes corresponds to the reaction numbered in Table S1 in Additional file 2. Lane M - GeneRuler 1 kb DNA Ladder (Thermo Scientific). Original uncropped and unprocessed gels are presented in Fig. S13e and f. Fig. S5. Comparison of selected sequence associated with reannotation of Palaeognathae mitogenomes. Eudromia elegans gene for tRNA-Pro with the last 200 nucleotides of the control region (a); Emeus crassus, Dinornis giganteus and Anomalopteryx didiformis control region 5′ spacers (b); Probosciger aterrinus (MH133970.1), Eolophus roseicapilla (MH133971.1) and Cacatua moluccensis (MH133972.1) ND6 pseudogenes (c); Anomalopteryx didiformis ND6 gene with control region 5′ spacer (d); Emeus crassus ND6 gene with control region 5′ spacer (e); Dinornis giganteus ND6 gene with control region 5′ spacer (f). Sequences rich in G and A residues are marked in a red box. The sequences in b and c were aligned with MUSCLE [125] in MEGA [124], whereas those in a, d, e and f using the optimal global:local algorithm (glsearch) form FASTA package version 36.3.8 g [127]. Fig. S6. Sequence comparison of two Eudromia elegans control regions. Dots indicate residues identical in the compared copies. Dots indicate residues identical in the compared sequences. Sequence repeats present in both control regions are marked in red and blue boxes. The sequences were aligned with MUSCLE [125] in MEGA [124]. Fig. S7. Influence of removing partitions with the highest substitution rate on alignment and tree parameters: the number of parsimony informative sites and mean distance (a) as well as mean support values (b). The mean phylogenetic distance was obtained from MrBayes tree. The mean support values were calculated from posterior probabilities received in MrBayes and PhyloBayes as well as SH-aLRT and non-parametric bootstrap percentages obtained in IQ-TREE. The posterior probabilities were scaled to 100%. Fig. S8. Maximum parsimony (a, c) and maximum likelihood (b, d) reconstruction of ancestral states and mapping of mitogenomic duplications onto the Palaeognathae tree topology t1 based on sequences of mitochondrial genes (a, b) and the Palaeognathae tree topology t2 supported by nuclear markers in coalescent-based approaches (c, d). In contrast to Fig. 6, this approach assumes that Anomalopteryx didiformis, Emeus crassus and Dinornis giganteus already lost mitogenomic duplication. The area of colors at nodes in b and d corresponds to the probability of the given state, single or duplicated region. Two-colored branches correspond to the equal probability of two states, single or duplicated region. The probability value for a more likely state was also given at these nodes. Mk1 model was applied for ML approach. Fig. S9. Annotated and reannotated duplicated gene orders in mitogenomes of Neoaves and Galloanserae deposited in GenBank for: Alca torda (CM018102.1), Antrostomus carolinensis (MN356120.1), Calypte anna (CM016612.1), Cariama cristata (CM020379.1), Catharus ustulatus (CM020378.1), Colius striatus (MN356125.1), Ketupa blackistoni (LC099104), Ketupa flavipes (LC099100.1), Leptosomus discolour (MN356135.1), Merops nubicus (CM020464.1), Melanereps aurifrons (CM022134.1), Morus serrator (GU071056.1), Nyctibius grandis (CM023771.1), Pterocles burchelli (MN356340.1), Pterocles gutturalis (CM020177.1, MN356147.1), Puffinus lherminieri (MH206162.1), Rynchops niger (MN356248.1), Sterna hirundo (CM020500.1), Thalassarche chlororhynchos (MN356342.1), Uria aalge (MN356418.1) and Urocolius indicus (MN356373.1). Genes for NADH dehydrogenase subunit 5 are marked in magenta coloring, cytochrome b in red coloring, tRNA-Thr in green coloring, tRNA-Pro in yellow coloring, NADH dehydrogenase subunit 6 in cyan coloring, tRNA-Glu in grey shading, 12S rRNA in blue fonts, tRNA-Phe in bold fonts, and control regions in red fonts. Pseudogenes are in italics and underlined. Microsatellites present within Pterocles gutturalis mitogenomes between tRNA-Thr gene and tRNA-Pro pseudogene or tRNA-Thr gene and ND6 pseudogene are marked in brown. Fig. S10. Sequences of mitogenomic fragments in Neoaves and Galloanserae mitogenomes amplified and sequenced for: Alca torda (MK263222), Apus apus (MW151827), Cathartes aura (MN629891), Chrysolophus pictus (MW151829), Ciconia nigra (MH264509), Corythaixoides personatus (MW082596), Gavia arctica (MK263210), Gavia stellata (MK263209), Eurypyga helias (MW208859), Podiceps cristatus (MN629890), Podiceps grisegena (MK263194), Scopus umbretta (MW151828), Spheniscus demersus (MH264510) and Uria aalge (MK263188). Genes for NADH dehydrogenase subunit 5 are marked in magenta coloring, cytochrome b in red coloring, tRNA-Thr in green coloring, tRNA-Pro in yellow coloring, NADH dehydrogenase subunit 6 in cyan coloring, tRNA-Glu in grey shading, 12S rRNA in blue fonts, tRNA-Phe in bold fonts, and control regions in red fonts. Pseudogenes are in italics and underlined. Fig. S11. Representative cropped agarose gels of the following amplicons: cytb/12S rRNA obtained for Rhynochetos jubatus (a), Eurypyga helias (d) and Trogon collaris (f); CR1/CR2 obtained for Rhynochetos jubatus (c) and Trogon collaris (g); CR/12S rRNA obtained for Rhynochetos jubatus (b) and Eurypyga helias (e); cytb/CR obtained for Eurypyga helias (e). The numbering of amplicons separated in the agarose gel corresponds to the reaction numbered in Table S9. Lane M - GeneRuler 1 kb DNA Ladder (Thermo Scientific). Fig. S12. Reconstruction of ancestral states and mapping of mitogenomic duplications onto the Aves three tree topologies under a more liberal assumption for assignment of duplication states. The colored thick lines of the trees correspond to maximum parsimony reconstruction, whereas values at nodes mean the probability of a more likely ancestral state (proportional likelihood) provided by maximum likelihood reconstruction. The trees obtained by Prum, et al. [45] are in a, by Kimball, et al. [47] in b, whereas by Kuhl, et al. [113] in c. The duplication state was assumed for a given order if a duplication was reported in at least one species from this order. FD indicates the presence of the most fully duplicated avian region. In ML approach, AsymmMk model was applied for a and, c whereas Mk1 for b. Fig. S13. Original uncropped and unprocessed agarose gels of: CR1/CR2 amplicons obtained for Struthio camelus (a), which were also shown in the Fig. S1a; CR1/CR2 amplicons obtained for Rhea pennata (b) which were also shown in the Fig. S1b; CR1/CR2 amplicons obtained for Rhea americana (7 lanes to the left of DNA Ladder) (c), which were also shown in the Fig. S1c; CR1/CR2 amplicons obtained for Crypturellus tataupa (6 lanes to the right of the DNA Ladder) (d), which were also shown in the Fig. S1d; potential NUMTs amplified for Dromaius novaehollandiae (e), which were also shown in the Fig. S4a; potential NUMTs amplified for Casuarius casuarius (f), which were also shown in the Fig. S4b; cytb/12S rRNA amplicons obtained for Rhynochetos jubatus (g), which were shown in the Fig. S11a; CR/12S rRNA amplicons obtained for Rhynochetos jubatus (h), which were shown in the Fig. S11b; CR1/CR2 amplicons (potential NUMTs) obtained for Rhynochetos jubatus (i), which were shown in the Fig. S11c; cytb/12S rRNA amplicon obtained for Eurypyga helias (j), which was shown in the Fig. S11d; cytb/CR and CR/12S rRNA amplicons obtained for Eurypyga helias (k), which were shown in the Fig. S11e; cytb/12S rRNA amplicon obtained for Trogon collaris (the first lane prior the first DNA Ladder) (l), which was shown in the Fig. S11f; CR1/CR2 amplicons obtained for Trogon collaris (m), which were shown in the Fig. S11g.