To carry out the phylogenetic assessment of the subgenus, a total of 25 species were analysed, which included all currently recognized species of Mesochila (20) plus four closely related Odontocheilina species, namely, Mesacanthina chalceola (Bates, 1872), Odontocheila nodicornis (Dejean, 1825), Pentacomia speculifera (Brullé, 1837) and Phyllodroma lutteomaculata Chaudoir, 1860. Finally, we used Opisthencentrus dentipennis Germar, 1843, to root the tree topology, because this species was never considered congeneric with the remaining species analysed, which all eventually were considered Odontocheila.

A total of 424 specimens of 14 South American species of Mesochila were examined, including six paratypes of M. proceroides Moravec, 2016, housed in the Coleção Entomológica do Instituto Oswaldo Cruz, Rio de Janeiro, Brazil (CEIOC). We also examined 14 specimens referable to the outgroup species. The Central American species of Mesochila (Paramesochila), M. (Mesochila) conformis (Dejean, 1831), M. (Mesochila) drechseli (Sawada and Wiesner, 1997) and M. (Mesochila) prepusula (Horn, 1907) were not found in the visited collections and thus not examined. However, their character states were possible to obtain from their extensive descriptions and detailed illustrations (Sawada and Wiesner 1997; Duran and Moravec 2013; Moravec and Brzoska 2013; Moravec 2016, 2018a, 2020).

The following institutions lent specimens for this study: CEIOC, Coleção Entomológica do Instituto Oswaldo Cruz, Rio de Janeiro, Brazil; DZRJ, Coleção Entomológica Professor José Alfredo Pinheiro Dutra, Departamento de Zoologia, Universidade Federal do Rio de Janeiro, RJ, Brazil; DZUP, Coleção Entomológica Pe. Jesus Santiago Moura, Universidade Federal do Paraná, PR, Brazil; MCNZ, Fundação Zoobotânica do Rio Grande do Sul, RS, Brazil; MNRJ, Museu Nacional, Universidade Federal do Rio de Janeiro, RJ, Brazil; MZSP, Museu de Zoologia, Universidade de São Paulo, SP, Brazil; PUCMG, Pontifícia Universidade Católica de Minas Gerais, MG, Brazil; UNIFEI, Universidade Federal de Itajubá, MG, Brazil.

A complete list of all examined material can be found in the Supporting information, Appendix S1.

The terminology follows Moravec (2018b, 2020). The male genitalia were left in a solution of 10% KOH during 48h in room temperature, in order to clarify the structure and make possible the observation of the sclerites in the inner sac. Habitus photographs were made using a Nikon 7000 digital camera with SIGMA 150 mm macrolens. The structures photographs were made on a stereomicroscopy Leica M205C with a coupled digital camera DFC 450 with the Application Suite CV3 montage software. The photographs were edited using Adobe Photoshop and the figure plates were designed with Adobe Illustrator (Adobe Systems).

The coding followed the logical basis presented by Sereno (2007). The matrix was built in MESQUITE v3.2 (Maddison and Maddison 2017), and can be found in the Supporting information, Table S1.

Recently, a debate was brought forth regarding the relative performance of probabilistic phylogenetic methods when dealing with morphological data (Goloboff et al. 2008a; Wright and Hillis 2014; Puttick et al. 2017; Goloboff et al. 2017; O’reilly et al. 2018). The argument lies on the whether Lewis (2001) model, as implemented in both maximum likelihood and Bayesian inference, provides accurate and precise tree topologies when compared to maximum parsimony (Brown et al. 2017; Puttick et al. 2017; Schrago et al. 2018). This prompted us to carry out a comparative evaluation of the performance of the major algorithms in our dataset.

Parsimony analyses were performed with TNT (Goloboff et al. 2008b), under equal (EW) and implied weights (IW) (Goloboff, 1993). All analyses were conducted using heuristic search with branch swapping TBR, employing 10,000 replicates and 100 trees saved for each replicate. For the IW analysis, we explored the topologies obtained under different concavity constant values (k=1, 2, 3, 5, 10, 100). The k values were not chosen with regular intervals because high k values tend to generate uniform and similar results to the ones of EW analysis (Goloboff 2008a). Therefore, analysis with regular k intervals may result in bias toward high k topologies (Mirande 2009). Support was assessed through Symmetric resampling (SR) with 5,000 replicates, because this method is not distorted by implied weights, and absolute Bremer was assessed for EW analysis only (Goloboff et al. 2003).

Character and character states were optimized in WINCLADA, ver. 1.00.08 (Nixon 2002). Following the sensibility analysis criterion (Giribet 2003), we employed the tree topology that was recovered more frequently in our analyses.

Maximum Likelihood (ML) analysis was performed on IQTREE (Nguyen et al. 2015) on its CIBIV Web Server (Trifinopoulos et al. 2016) using ultrafast bootstrap as branch support (Hoang et al. 2018). Bayesian Inference (BI) was performed on MRBAYES 3.2 (Huelsenbeck and Ronquist 2001). The posterior distribution of tree topologies was approximated by the Markov chain Monte Carlo (MCMC) algorithm, which was run for 10,000,000 generations and sampled every 1,000th cycle, with 10% of the initial chains discarded as burn-in. Chain convergence was evaluated in TRACER v1.6 (Rambaut et al. 2014). In both ML and BI, was used the MKV model of morphology evolution, modified from MK (Lewis 2001). This model implements a Markov process for discrete character evolution accounting for acquisition bias.

The geographic distribution of Mesochila species and other species on the phylogeny was arranged in five different areas, based on Morrone’s dominions (2014a): A (Pacific), B (Boreal Brazilian), C (South Brazilian), D (Chacoan) and E (Paraná). The South-eastern Amazonian dominion had no species recorded. The information about their distribution was based on recent revisions (Pearson et al. 1999; Moravec 2016, 2018; Roza and Mermudes 2017; Moravec 2020). A map was built using QUANTUM-GIS 2.14.3 (QGIS Development Team 2016) with the shapefile of Morrone’s dominions (Löwenberg-Neto 2014).

We used S-DIVA (Yu et al. 2010) and Bayesian Binary Method (BBM) (Ronquist and Huelsenbeck 2003) implemented in RASP (Yu et al. 2015) to reconstruct the possible ancestral ranges of Mesochila on the phylogenetic trees. To account for phylogenetic uncertainty, we used 20,004 trees from the MCMC output and ran both methods. The BBM was run with the fixed state frequencies model (Jukes-Cantor) with equal among-site rate variation for two million generations, ten chains each, and two parallel runs.

The maximum number of areas in the nodes was set to 5, matching the number of dominions in which the species occurs and allowing for all the possibilities of composed ancestral areas (Kodandaramaiah 2010). The closest outgroup taxa (or clade, if it is the case) was retained in the biogeographical analysis to avoid the bias towards a widespread ancestral in the root (Ronquist 1997; Kodandaramaiah op. cit.).

We derived 37 characters for investigating Mesochila phylogeny, encompassing body size, coloration and external covering, head, thorax, membranous wings and male genitalia. The characters are coded as binary (25) or multistate (12). The characters used in the analysis are as follows (include the length (L), consistent index (CI) and retention index (RI) for the EW analysis):

1. Body Size: (0) 6.0 to 9.3 mm; (1) 9.5 to 14 mm. L = 3; CI = 0.33; RI = 0.81.

2. Head and pronotum, dorsal surface, coloration: (0) green (Fig. 2A); (1) dark green (Fig. 2B); (2) purple (Fig. 2C); (3) black (Fig. 2D); (4) brownish green (Fig. 2E); (5) copper (Fig. 2F). L = 8; CI = 0.62; RI = 0.70.

3. Labrum of the male, relation of width to length: (0) at least 2 times wider than long (Fig. 3G); (1) 1.4 to 1.5 times wider than long (Fig. 3A). L = 4; CI = 0.33; RI = 0.33.

4. Labrum, base coloration in males: (0) green (Fig. 3A); (1) black (Fig. 3B); (2) yellowish brown (Fig. 3C); (3) pale yellow (Fig. 3E). L = 4; CI = 0.75; RI = 0.75.

5. Labrum, base color in females: (0) green; (1) black; (2) yellowish brown; (3) pale yellow. L = 5; CI = 0.60; RI = 0.60.

6. Labrum of the male, median longitudinal macula: (0) absent; (1) present (Fig. 3B). L = 2; CI = 0.50; RI = 0.66.

7. Labrum, margin of laterobasal region, shape: (0) smooth (Fig. 3D); (1) rhomboid (Fig. 3A);(2) long and curved teeth (Fig. 3F); (3) angular (Fig. 3B); (4) short teeth (Fig. 3G). L = 6; CI = 0.66; RI = 0.71.

8. Labrum, lateromedial margin, shape: (0) smooth (Fig. 3J); (1) rhomboid (Fig. 3A); (2) long and curved teeth (Fig. 3F); (3) rounded (Fig. 3D). L = 5; CI = 0.75; RI = 0.66.

9. Labrum, notch before the apical lobe when lateromedial margin is close to apical lobe length: (0) absent; (1) present (Fig. 3E). L = 3; CI = 0.33; RI = 0.60.

10. Labrum, apical lobe, form: (0) smooth or with vestigial teeth (Fig. 3C); (1) two pairs of lateral teeth with a longer medial (Fig. 3J); (2) very long median tooth (Fig. 3H); (3) three long sub-equal teeth (Fig. 3A); (4) a pair of long lateral teeth (Fig. 3F); (5) three short sub-equal teeth (Fig. 3E). L = 5; CI = 1.0; RI = 1.0.

11. Mandible, teeth, quantity: (0) three (Fig. 3K); (1) four (Fig. 3L). L = 1; CI = 1.0; RI = 1.0.

12. Maxillary and labial palpi, coloration: (0) brown, latest article black (Fig. 3M); (1) white, last article black (Fig. 3N); (2) all black (Fig. 3O); (3) all white (Moravec 2020: pl. 71, fig. A). L = 3; CI = 1.0; RI = 1.0.

13. Pronotum, relation of length to width: (0) distinctly wider than long (Fig. 2C); (1) as long as wide (Fig. 2A); (2) distinctly longer than wide (Fig. 2F). L = 6; CI = 0.33; RI = 0.42.

14. Elytra, punctuation: (0) finely punctuated with small and shallow grooves (Fig. 4E); (1) coarsely punctuated with large and deep grooves (Fig. 4F). L = 1; CI = 1.0; RI = 1.0.

15. Elytra, relation of length with width: (0) 3.9–4.2 times longer than wide (Fig. 4A); (1) 3.3–3.6 times longer than wide (Fig. 4D). L = 3; CI = 0.33; RI = 0.81.

16. Elytra, humeral spot, occurrence: (0) absent; (1) present (Fig. 4B). L = 4; CI = 0.25; RI = 0.50.

17. Elytra, humeral spot, position: (0) restricted to humerus (Fig. 4H); (1) extending to the end of the basal third of the lateral margin (Fig. 4G). L = 3; CI = 0.33; RI = 0.50.

18. Elytra, lateral and post-middle spot, occurrence: (0) absent; (1) present (Fig. 4A). L = 2; CI = 0.50; RI = 0.0.

19. Elytra, lateral and post-medial spot, shape: (0) subround to subtriangular (Fig. 4A); (1) distinctly triangular (Fig. 4B); (2) rectangular (Fig. 4C); (3) corrugated (Moravec and Brzoska 2014: fig. 13–15); (4) “half T”-shaped (Fig. 4D). L = 7; CI = 0.57; RI = 0.40.

20. Elytra, latero-apical spot, occurrence: (0) absent; (1) present (Fig. 4B). L = 3; CI = 0.33; RI = 0.50.

21. Elytra, latero-apical spot, dorsal length: (0) reaches approximately half the width of the elytron (Fig. 4C); (1) approximately reaches the elytral suture (Fig. 4D). L = 1; CI = 1.0; RI = 1.0.

22. Membranous wing, radial cell, pigmentation: (0) partial (Fig. 5A); (1) total (Fig. 5B). L = 1; CI = 1.0; RI = 1.0.

23. Membranous wing, cr vein, projection, occurrence: (0) absent (Fig. 5A); (1) present (Fig. 5B). L = 1; CI = 1.0; RI = 1.0.

24. Aedeagus, apex, shape: (0) rounded (Fig. 6A); (1) distinctly hook (Fig. 6B); (2) distinctly wide bent rounded beak (Fig. 6C). L = 4; CI = 0.50; RI = 0.78.

25. Aedeagus, apical third, shape: (0) subparallels (Fig. 6C); (1) gradually tapered (6A, 6B); (2) prolonged into narrow, cylindrical and rounded apex (Moravec 2020: pl. 72, figs A–G). L = 3; CI = 0.66; RI = 0.83.

26. Aedeagus, basal part, width: (0) subparallel (Fig. 6A); (1) wide (Moravec 2020: pl. 72, figs A–G). L = 1; CI = 1.0; RI = 1.0.

27. Aedeagus, middle part, width: (0) thin (1/7 of aedeagus length) (Fig. 6B); (1) moderately wide (1/6 of aedeagus length) (Fig. 6A); (2) wide (1/5 of aedeagus length) (Fig. 6E). L = 3; CI = 0.66; RI = 0.89.

28. Aedeagus, sclerites, position inside of aedeagus: (0) apical (Fig. 6D); (1) apical-medial (Fig. 6A); (2) on the entire aedeagus (Fig. 6E). L = 1; CI = 1.0; RI = 1.0.

29. Aedeagus, wavy and spatulate dorsal sclerites, occurrence: (0) absent; (1) present (Fig. 7A). L = 1; CI = 1.0; RI = 1.0.

30. Aedeagus, right dorsal sclerite rod shaped, occurrence: (0) absent; (1) present (Fig. 6B). L = 3; CI = 0.33; RI = 0.33.

31. Aedeagus, proximal sclerites formed by two pieces, occurrence: (0) absent; (1) present (Fig. 7B, 7C). L = 1; CI = 1.0; RI = 1.0.

32. Aedeagus, proximal sclerites form by two pieces, shape: (0) triangular (Fig. 7D); (1) subquadrate (Fig. 7C). L = 2; CI = 0.50; RI = 0.86.

33. Aedeagus, longitudinal ventral sclerite parallel to aedeagus (ventral spur), occurrence: (0) absent; (1) present (Fig. 6A). L = 1; CI = 1.0; RI = 1.0.

34. Aedeagus, longitudinal ventral sclerite parallel to aedeagus, base shape: (0) Stingray-like (Fig. 7C); (1) triangular (Moravec 2020: pl. 72, fig. A–G); (2) bifurcate (Fig. 7E); (3) u-shaped (Fig. 7F). L = 4; CI = 0.75; RI = 0.80.

35. Aedeagus, quadrangular central sclerite with lower sclerotized margin, occurrence: (0) absent; (1) present (Fig. 7C). L = 1; CI = 1.0; RI = 1.0.

36. Aedeagus, quadrangular central sclerite with projected lower left tip, occurrence: (0) absent; (1) present (Fig. 8A). L = 1; CI = 1.0; RI = 1.0.

37. Aedeagus, oval feebly sclerotized central sclerite, occurrence: (0) absent; (1) present (Fig. 8B). L = 1; CI = 1.0; RI = 1.0.

Equally weighted parsimony analysis of the 37 parsimony informative characters resulted in two equally parsimonious trees (L = 106, CI = 0.61, RI = 0.78), both exhibiting identical relationships between Mesochila species (Consensus – Fig. 9A, 9B). Implied weights analyses under different concavity constant values (k = 1, 2, 3, 5, 10, 100) yielded exactly the same topology of equally weighted analysis, with the exception of K (1), which recovered a different topology for three ingroup taxa (which will be discussed below).

Mesochila was recovered as monophyletic in all analyses, supported by both Symmetric resampling (64–73) and Bremer (4) supports. The genus is supported by the following synapomorphies: labrum, margin of laterobasal region, shape: angular (7:3, non-homoplastic); aedeagus, sclerites, internal sac positioning: apical-medial (28:1, non-homoplastic); aedeagus, proximal sclerites form by two pieces, occurrence: present (31:1, non-homoplastic); aedeagus, longitudinal ventral sclerite parallel to aedeagus (ventral spur), occurrence: present (33:1, non-homoplastic).

The earliest branching of the Mesochila genus separated Mesochila (Paramesochila) from the remaining species. A subsequent split formed two clades containing Mesochila (Mesochila) and Mesochila (Eumesochila) species, which consisted of reciprocally paraphyletic lineages.

The first clade, henceforth referred to as clade A, which contained all species of the subgenus M. (Paramesochila), was supported by both Symmetric resampling (55–71) and Bremer (3) supports. The synapomorphies for the subgenus were: mandible, teeth, quantity: three (11:0, non-homoplastic); elytra, punctuation, shape: coarsely punctuated with large and deep grooves (14:1, non-homoplastic); aedeagus, quadrangular central sclerite with projected lower left tip, occurrence: present (36:1, non-homoplastic). Mesochila (P.) horni Schilder, 1953 was found as sister to all other M. (Paramesochila) species.

The phylogenetic relationships of the remaining species – M. (P.) wappesi (Moravec and Brzoska, 2013), M. (P.) skrabali (Duran and Moravec, 2013) and M. (P.) tayutica Moravec, 2018 – were recovered as a polytomy in all analyses. The evolutionary affinity of these three species was the better-supported within the genus, with high values of both Symmetric Resampling (99) and Bremer (12) supports. The synapomorphies of the Central American species were: head and pronotum, dorsal surface, coloration: dark green (2:1, non-homoplastic); maxillary and labial palpi, coloration: all white (12:3, non-homoplastic); elytra, humeral spot, position: (0) restricted to humerus (17:0, homoplastic); elytra, latero-apical spot, dorsal length: approximately reaches the elytral suture (21:1, non-homoplastic); aedeagus, apical third, shape: prolonged into narrow, cylindrical and rounded apex (25:2, homoplastic); aedeagus, basal part, width: wide (26:1, non-homoplastic); aedeagus, middle part, width: wide (27:2, homoplastic); aedeagus, right dorsal sclerite rod shaped, occurrence: present (30:1, homoplastic); aedeagus, longitudinal ventral sclerite parallel to aedeagus, base shape: triangular (35:1, non-homoplastic).

A second clade, named clade B, was formed by M. (Mesochila) and M. (Eumesochila) species, although with low support values. It was characterized by a single synapomorphy: aedeagus, apex, shape: thin hook (24:1, non-homoplastic). Clade B was divided into two subclades: clade C, composed of M. (E.) distigma (Dejean, 1825) + M. (E.) discrepans (Horn, 1893) as sister group to M. (M.) drechseli (Sawada and Wiesner, 1997) + a polytomy containing M. (M.) brevipennis (W. Horn, 1907), M. (M.) prepusula (W. Horn, 1907) and M. (M.) moraveci Roza and Mermudes, 2019; and the other (Clade D) composed of the remaining species of the genus.

The clade C has low support values with two synapomorphies: elytra, humeral spot, occurrence: absent (16:0, homoplastic); aedeagus, middle part, width: thin (27:0, non-homoplastic). The clade M. (E.) distigma + M. (E.) discrepans was supported by both Symmetric resampling (50–65), at least in most analyses, and Bremer (3) support. The clade was sustained by two synapomorphies: body size: 9.5 to 14 mm (1:1, homoplastic); pronotum, relation of length to width: distinctly wider than long (13:0, homoplastic). Its sister clade, M. (M.) drechseli + a polytomy of M. (M.) brevipennis, M. (M.) prepusula and M. (M.) moraveci (herein designed as brevipennis species-group) is supported by both Symmetric Resampling (80–89) and Bremer (8) supports, and presents the following synapomorphies: aedeagus, apex, shape: distinctly wide bent rounded beak (24:2, non-homoplastic); aedeagus, apical third, shape: subparallels (25:0, non-homoplastic); aedeagus, longitudinal ventral sclerite parallel to aedeagus, base shape: stingray-like (34:0, homoplastic); aedeagus, oval feebly sclerotized central sclerite, occurrence: present (37:1, non-homoplastic). The polytomy implicate M. (M.) brevipennis, M. (M.) prepusula and M. (M.) moraveci has low support values, and only one synapomorphy: elytra, lateral and post-medial spot, shape: distinctly triangular (19:1, homoplastic).

The clade D had low Bremer support and Symmetric Resampling values of above 50 only in two IW analyses. The group has three synapomorphies: elytra, relation of length with width: 3.9–4.2 times longer than wide (15:0, homoplastic); aedeagus, longitudinal ventral sclerite parallel to aedeagus, base shape: u-shaped (34:3, non-homoplastic); aedeagus, quadrangular central sclerite with lower sclerotized margin, occurrence: present (35:1, non-homoplastic). M. (M.) brasiliensis (Dejean, 1825) was found as the sister group of the rest of the species. Next, a polytomy of M. (M.) conformis (Dejean, 1831) and M. (E.) distincta (Dejean, 1831) was recovered as the sister group to the rest of Clade D. This topology has low support value in all analyses for both supports and has only three synapomorphies: body size: 9.5 to 14 mm (1:1, homoplastic); head and pronotum, dorsal surface, coloration: brownish green (2:4, homoplastic); elytra, humeral spot, position: restricted to humerus (17:0, homoplastic).

The only different tree recovered a clade formed by (M. (M.) conformis + (M. (M.) brasiliensis + M. (E.) distincta)) as sister to the rest of clade D in the K1 analysis. This topology, however, presented low support value.

Next, the lade E formed by M. (M.) procera (Chaudoir, 1860) and M. (M.) proceroides Moravec, 2016 recovered as sister of M. (M.) biguttata (Dejean, 1825) and all the smaragdula species-group sensu Moravec (2020). This position is supported by Symmetrical Resampling (60–74) and Bremer (2) supports. They share the following synapomorphies: membranous wing, radial cell, pigmentation: total (22:1, non-homoplastic); membranous wing, cr vein, projection, occurrence: present (23:1, non-homoplastic); aedeagus, wavy and spatulated dorsal sclerites, occurrence: present (29:1, homoplastic); aedeagus, longitudinal ventral sclerite parallel to aedeagus, base shape: Stingray-like (34:0, homoplastic). The clade M. (M.) procera + M. (M.) proceroides was well supported by Symmetrical resampling (63–65) and Bremer (4) supports, but shares only one synapomorphy: pronotum, relation of length to width: distinctly longer than wide (13:2, homoplastic).

The clade formed by M. (M.) biguttata and all the smaragdula species-group was highly supported by Symmetrical Resampling (83–86) and Bremer (11) supports, and share five synapomorphies: head and pronotum, dorsal surface, coloration: green (2:0, non-homoplastic); labrum, base coloration in males: black (4:1, non-homoplastic); labrum, base color in females: black (5:1, homoplastic); labrum, median longitudinal macula in the male, occurrence: present (6:1, non-homoplastic); elytra, latero-apical spot, occurrence: absent (20:0, homoplastic).

Inside smaragdula species-group, it was found a polytomy of M. (M.) smaragdula (Dejean, 1825), M. (M.) completemaculata (Horn, 1922) and a clade formed by M. (M.) viridis (Dejean, 1831) and M. (M.) cyaneomarginata (Horn, 1920). It is supported by Bremer (2) support, but has low support from Symmetrical Resampling, except in K10 and K100 analyses (51). Presents four synapomorphies: labrum, margin of laterobasal region, shape: (1) rhomboid (7:1, homoplastic); labrum, notch before the apical lobe when lateromedial margin is close to apical lobe length, occurrence: absent (9:0, homoplastic); labrum, apical lobe, form: three long sub-equal teeth (10:3, non-homoplastic); aedeagus, apex, shape: rounded (24:0, homoplastic). The clade formed by M. (M.) viridis and M. (M.) cyaneomarginata is supported by Bremer support (2) but has low values of Symmetrical resampling in all analyses. It shares two synapomorphies: labrum, base coloration in males: green (4:0); labrum, base coloration in females: green (5:0).

ML analysis recovered a tree with log-likelihood score of -ln 471.111 (Fig. 10A), with topological relationships very similar to that of maximum parsimony. Mesochila was also recovered as monophyletic, with high support of Ultrafast bootstrap (95). The main differences between ML and MP trees concerned the resolution of the polytomies recovered by MP. The ML-resolved groups consisted of: (M. (P.) tayutica + (M. (P.) wappesi + M. (P.) skrabali)) in clade A, (M. (M.) moraveci + (M. (M.) drechseli + (M. (M.) brevipennis + M. (M.) prepusula))) in clade C, (M. (M.) brasiliensis + (M. (E.) distincta + (M. (M.) conformis + rest of clade D))) in clade D, and (M. (M.) smaragdula + (M. (M.) completemaculata + (M. (M.) viridis + M. (M.) cyaneomarginata))) in smaragdula species group. However, most of these groups were poorly supported by ultrafast bootstrap values below 95 (Minh et al. 2013). The exceptions were clade A, the Central American species clade, Clade E and M. (M.) biguttata + smaragdula species group clade.

BI also recovered a tree topology (Fig. 10B) similar to the MP tree, but the occurrence of polytomies was higher. Apart from Mesochila, which was consistently monophyletic with a high posterior clade probability (PP = 0.99), the Central American species clade, M. (M.) brevipennis clade, clade E and M. (M.) biguttata + smaragdula species group clade (all with PP > 95), the statistical support for clades in the Bayesian tree was overall low. BI–MP differences were as follows: clade A was recovered in a polytomy with clade C and clade D; The Central American species were also recovered as a polytomy; M. (M.) brevipennis clade was also found as a polytomy; M. (M.) brasiliensis and M. (E.) distincta were found in a polytomy at the base of clade D.

The map of the distribution (Fig. 11) showed that most species are restrained to one or two dominions, with the exception of M. (E.) discrepans and M. (E.) distigma. Both these species are widely amplilocated. We performed the S-DIVA and BMM analyses for the topology yielded by the ML analysis (Fig. 12, 13), since both methods require a polytomy-free tree topology as input. Because topological variation between phylogenetic methods was not significant, and the ML-resolved polytomies concerned mainly lineages with the same dominion distribution, we expect that the adoption of a fully resolved tree had little impact in the biogeographic analysis.

The S-DIVA analysis inferred 18 dispersal and three vicariant events (Fig. 12). The most likely ancestral area of Mesochila (node 42) was recovered as area E with 51% marginal probability. Other possibilities were area D (29% probability) and area DE (18% probability). The analysis estimated a 5% probability of two dispersion events from this original area E to areas E and EAD separately. The ancestral distribution of the clade A (node 26) was assigned to ADE (35% probability), AE (34%) or AD (31%). There is 34% probability of a vicariant event from area ADE to area A and DE, and A is the ancestral area of the Central American species (node 25) with 100% probability. The following node on the Central American clade (node 24) has ancestral area A (100% probability) with no dispersal or vicariant events.

The ancestral area of the clade B (node 41) was recovered as E (35% probability), followed by DE (33%) or D (31%). The analysis recovered a 9% chance of one vicariant event separating area DE into areas D and E. Clade C (node 31) has the ancestral D (41% probability) or DE (36%), with alternative areas E (18%), with 15% chance of one dispersal event, D to DE, and one vicariant event separating area DE into area D and E. Inside Clade C, M. (E.) distigma + M. (E.) discrepans clade (node 27) has the ancestral area of E (50% probability) or D (49%). There is 50% probability of six dispersal events, which rendered the current amplilocate distribution of M. (E.) distigma and M. (E.) discrepans. In the other clade, brevipennis species-group (node 30) had the ancestral area D (73% probability), or alternatively DE (27%). The subsequent ancestral areas of the group were also most likely D (100% probability). There is one dispersion in node 29 (for M. drechseli – DE) and other in node 28 (for M. brevipennis – CD).

The clade D has the ancestral area (node 40) of E (66% probability), or alternatively DE (33%). There is 70% possibility of one dispersal event from area E to area DE. All subsequent nodes have most likely the ancestral area of E (100% probability), with dispersal events happening on node 36 (100% chance for M. biguttata – DE).

The BMM analysis resulted in the possibility of 22 dispersal and one vicariance events (Fig. 13). The most likely ancestral area of Mesochila (node 42) was recovered as a composed area of two different areas (DE) with 69% marginal probability. Other possibilities were area D (17% probability) and area E (8% probability). There are 32% probability of two events of dispersion from this original DE area to areas D and E separately, and posterior colonization of area DE independently in the two following ancestral nodes. The ancestral distribution of the clade A (node 26) was found as most probable as area DE (with 57% probability), and alternative areas being D (27%) or E (10%). There is 50% probability of a dispersal to area DEA, with a subsequent vicariant event separating area A, which is ancestral area of the Central American species (node 25) with 88% probability (alternative area AD presented 6% probability). The following node on the Central American clade (node 24) has ancestral area A (99% probability) with no dispersal or vicariant events.

The ancestral area of the clade B (node 41) was recovered as DE (81% probability), with alternative area D (8%) or CDE (6%). There is 39% chance of two dispersal events in similar way of node 42. Clade C (node 31) has the ancestral DE (55% probability), with alternative areas CDE (26%) and D (10%). There is 18% chance of two dispersal events, one from are DE to area D, and other from are DE to area CDE. Inside Clade C, M. (E.) distigma + M. (E.) discrepans clade (node 27) has the ancestral area of CDE (63% probability), with two alternative ancestral areas: BCDE (17%) and ACDE (8%). There is 62% probability of five dispersal events, which rendered the current amplilocate distribution of M. (E.) distigma and M. (E.) discrepans. In the other clade, brevipennis species-group (node 30) had the ancestral area D (53% probability), or alternatively DE (35%) or CD (7%). The subsequent ancestral areas of the group were also most likely D. There is one dispersion in node 29 (for M. drechseli – DE) and other in node 28 (for M. brevipennis – CD).

Clade D has the ancestral area (node 40) of DE (87% probability), or alternatively E (10%). There is 77% possibility of one dispersal event from area DE to area E. All subsequent nodes have most likely the ancestral area of E (at least 70% probability), with dispersal events happening on node 36 (67% chance for M. biguttata – DE).

Here we presented the first phylogeny of an Odontocheilina lineage. Tree topologies obtained by MP and ML were well resolved, with the exception of two internal clades, in clade A and C, and four taxa inside clade D. BI yielded considerably more polytomies, including the earliest branching of Mesochila. Overall, homoplasy likely had little impact on the results, since the topology found by the implied weight analysis was identical to the equal weight parsimony analysis. Although most analyses generated groups with low support values, the main clades inside Mesochila (clade A, M. (E.) discrepans + M. (E.) distigma clade, brevipennis species-group and clade D) were usually well supported. Besides that, topologies were consistent across different analyses, which corroborates employing this phylogeny to evaluate the systematics of this group.

The topologies generated by MP and ML exhibited better resolutions when compared with the Bayesian tree. In theory, previous studies have reported that BI of morphological data is accurate, meaning the true tree frequently lies within the 95% credibility interval of topologies. However, higher accuracy generally comes at the expense of lower precision, resulting in trees with lower resolution, especially in small datasets, as the present study (O’reilly et al. 2016; Puttick et al. 2017; Schrago et al. 2018). This prompted some authors to advocate in favour of the implied weight analysis for morphological characters (Goloboff et al. 2017). The debate around this subject is ongoing and generally there is a lack of studies employing empirical analyses instead of simulations. The only exception is the recent study of Schrago et al. (2018), who carried out an extensive comparative evaluation of MP and BI using more than 100 published matrices from MorphoBank. The authors confirmed the results from simulations, showing that BI under the MKv model generates wide credibility intervals and, consequently, lower resolution.

Based on the present study, Mesochila is a well-supported monophyletic group, including all 20 species previously assigned to the genera. We obtained, however, a different configuration for the subgenera. M. (Paramesochila) was recovered as a well-supported phylogenetic group, with M. (P.) horni as sister group to the remaining species. The other subgenera, M. (Mesochila) and M. (Eumesochila) were recovered as reciprocally paraphyletic.

We identified three main clades, in clade B, which were recovered in all analyses, usually well-supported: M. (E.) discrepans + M. (E.) distigma (which corresponds with M. (Eumesochila) with the exclusion of M. (E.) distincta; brevipennis species-group, which was placed outside of the main M. (Mesochila) clade; and clade C, which includes the remaining M. (Mesochila) species with the inclusion of M. (E.) distincta. The position of three taxa at clade C, M. (M.) brasiliensis, M. (M.) conformis and M. (E.) distincta, is weakly supported and varied between analyses. Their morphological resemblance with the rest of clade D, however, makes their inclusion in this group more likely. Due to the fact that this four previous mentioned clades were stable in all analyses carried out, frequently highly supported, we were prompted to propose a new infrageneric configuration for the genus:

The S-DIVA analysis presented competing ancestral areas for Mesochila and its more internal nodes, many of them with very close probabilities. The results of BMM were more consistent, with most of the more probable areas with high probabilities. This lack of resolution is a common issue in S-DIVA (Kodandaramaiah 2010). We tried to attenuate this bias including the sister group of the focus taxa of the analysis (Ronquist 1997; Kodandaramaiah op. cit.).

This measure, however, seems to have had little effect on the results of this analysis. Even if the results have been more resolved, the sister group of Mesochila is unknown. Pentacomia speculifera and Mesacanthina chalceola were retained in the analysis because they were found as the sister clade of the ingroup in most phylogenetic results. But, as the phylogeny did not have the objective of testing the relationship of Odontocheilina genera, the relationship between Pentacomia, Mesacanthina and Mesochila cannot be trusted. Therefore, any resolved result for the ancestral area of the genus should be taken with caution, at least. However, when comparing both S-DIVA and BMM results, it was clear that the ancestral area was Chacoan dominion, Parana dominion or a composed area of the two. Those dominions encompass Atlantic rainforest, Cerrado and Caatinga biomes.

The recovery of the ancestral split in M. (Paramesochila), with an ancestral area ADE or with posterior dispersal to this area shows an inedited relationship between the south of Central America and the Atlantic rainforest (In this case, M. horni occurs in the Atlantic rainforest biome, sometimes in the border with Caatinga biome). The possible vicariant event here suggests that those two regions were in contact sometime prior to the speciation event between these two lineages. Pleistocenic shifts of vegetation are postulated as a probable cause for the link between Amazon and Atlantic rainforest (Costa 2003; Sobral-Souza 2015), and these shifts may have also connected the northwest Amazon and the Caatinga (De Oliveira et al. 2015). Future studies in the group, with the descriptions of new species with a median distribution between these two areas, may or may not clarify this question.

The recovery of an ancestral area of D, E or DE for the clade B suggests that M. (Mesochila), M. (Eumesochila) and M. (Neomesochila) subgen. nov. diversified in these dominions (mainly in Cerrado and Atlantic rainforest), and only recently some species have dispersed to other areas (like Pacific, Boreal Brazilian and South Brazilian dominions, which are mainly covered by tropical wet forests or drier lowland vegetation). Another interesting point is the absence of any species in the South-eastern Amazonian dominion, and almost complete absence in the Caatinga biome (in Chacoan dominion). The amplilocated distribution of P. (M.) discrepans and P. (M.) distigma due to dispersion events may correspond to the hypotheses of Pleistocenic shifts of vegetation that cause a posterior link between Amazon and Atlantic rainforest (Costa 2003; Sobral-Souza 2015). These species could be unable to establish themselves in the above-mentioned areas, or the absence in South-eastern Amazonian dominion and Caatinga may be an artefact of the lack of sampling in those locations.

Finally, as we do not have a hypothesis for the divergence time of the lineages, any cause for these speciations postulated here is merely speculative. Future studies using molecular data can give more insights about the biogeographical events occurring in the history of this group.