Taxon Sampling

We selected 153 ingroup species belonging to 56 crematogastrine genera for our molecular study. These included the 100 Crematogastrini species used in Ward et al. (2015); in all but one case the remaining DNA extracts from that study were used. Fifty-three new taxa were added to achieve a nearly complete sampling of the 64 genera of the tribe. All but four described genera (Formosimyrma Terayama, Goaligongidris Xu, Peronomyrmex Viehmeyer, Secostruma Bolton) are included in the present analysis. We further increased our taxon sampling at the species level within some genera to mitigate long branches resulting from single-species representation in the previous study of the subfamily (Ward et al. 2015). Eight outgroup taxa were included: five species representing members of other tribes of the subfamily Myrmicinae (Manica bradleyi Wheeler, Pogonomyrmex imberbiculus Wheeler, Aphaenogaster swammerdami Forel, Monomorium hanneli Forel, Pheidole longispinosa Forel), as well as members of three other closely related subfamilies (Ectatomminae: Rhytidoponera chalybaea Emery, Heteroponerinae: Acanthoponera minor Forel, Formicinae: Formica moki Wheeler).  Supp Table 1 (ixy013_suppl_supplementary_table_s1_s7.pdf) (online only) contains voucher information for all taxa and further specimen data are available on AntWeb ( www.antweb.org). Voucher specimens have been deposited at the University of California, Davis; at the National Museum of Natural History, Smithsonian Institution; at the California Academy of Sciences; and at the North Carolina State University Insect Museum. Specimens included in this study were collected in accordance with local regulations and all necessary permits were obtained.

UCE Data Collection

We extracted DNA destructively or nondestructively (specimen retained after extraction) from worker ants or pupae using the DNeasy Blood and Tissue Kit (Qiagen, Valencia, CA). We quantified genomic DNA for each sample using a Qubit fluorometer (High sensitivity kit, Life Technologies, Inc., Carlsbad, CA) and sheared between <5 and 300 ng for 10–60 s (amp = 25, pulse = 10) to a target size of approximately 500–600 bp by sonication (Q800, Qsonica Inc., Newtown, CT). This sheared DNA was used as the input for a modified genomic DNA library preparation protocol (Kapa Hyper Prep Library Kit, Kapa Biosystems, Wilmington, MA) that incorporated ‘with-bead’ cleanup steps (Fisher et al. 2011) and a generic SPRI substitute (Rohland and Reich 2012, ‘speedbeads’ hereafter), as described by Faircloth et al. (2015). We used TruSeqstyle adapters during adapter ligation (Faircloth and Glenn 2012) and combined groups of 9–10 libraries at equimolar ratios into enrichment pools having final concentrations of 84.7–160 ng/µl. We enriched each pool using a set of custom-designed probes (MYcroarray, Inc., now ArborBiosciences, Ann Arbor, MI) targeting 2,590 UCE loci in ants (Branstetter et al. 2017a). We followed target enrichment procedures for the MYcroarray MYBaits kit (Blumenstiel et al. 2010) except we used a 0.2× concentration of the standard MYBaits concentration and added 0.7 µl of 500 µM custom blocking oligos designed against our custom sequence tags. We ran the hybridization reaction for 24 h at 65°C, subsequently bound all pools to streptavidin beads (Dynabeads MyOne Streptavidin T1, Life Technologies, Inc., Carlsbad, CA) and washed bound libraries according to a standard target capture protocol (Blumenstiel et al. 2010). We used the with-bead approach for PCR recovery of enriched libraries as described in Faircloth et al. (2015). We performed qPCR library quantification and combined libraries at equimolar concentrations into two final pools based on the estimated size-adjusted concentrations. We size-selected these final pools for 250–800 bp with a BluePippin (SageScience, Beverly, MA). For a more detailed description of library and enrichment protocols refer to Blaimer et al. (2016a,b). All of the UCE laboratory work was conducted in and with support of the Laboratories of Analytical Biology (L.A.B.) facilities of the National Museum of Natural History. The pooled libraries were sequenced using two lanes of a 125-bp paired-end Illumina HiSeq 2500 run at the University of Utah’s Huntsman Cancer Institute. Quality-trimmed sequence reads generated as part of this study are available from the NCBI Sequence Read Archive under accession SRP149549; trinity assembled contig files are deposited in the Dryad data repository, together with data matrices and tree files (doi:10.5061/ dryad.121b453).

Phylogenetic Inference

Dating Analyses

We inferred divergence dates within Crematogastrini from the subset of ‘50-best’ UCE loci (as described above), as well as from the two subsets of ‘50-random’ loci using BEAST2 v2.3.2 (Bouckaert et al. 2014) run via the SI/HPC. Analyses incorporated data partitioning as estimated by the greedy algorithm in PartitionFinder2. We used two alternative starting trees: 1) the best tree resulting from the concatenated RAxML analysis with 1,763 partitions and 2) the species tree estimated by ASTRAL-II (used in combination with 50-best data set only). Due to computational constraints, topology was fixed during analyses to these starting trees by setting weight = 0 for all operators acting on the tree model. Four runs each for the 50-best data and two runs each for the 50-random data sets were employed with a chain length of 220–300 million generations. We imposed a diffuse gamma distribution on the mean branch lengths (ucld.mean; alpha = 0.001, beta = 1,000) and set the starting value for the clock rate = 1e-9, but otherwise left priors at their default values. We chose a Yule tree prior and defined lognormal calibration priors based on seven calibrations used by Ward et al. (2015) (see  Supp Table 2 [online only] (ixy013_suppl_supplementary_table_s1_s7.pdf)); however, we incorporated an updated younger estimate for Baltic amber (Aleksandrova and Zaporozhets 2008) and Sicilian amber (LaPolla et al. 2013) into these calibration priors. An important consideration when using lognormal calibration priors for divergence dating is how the user-specified calibration priors interact with the tree prior to form the marginal prior (or joint prior) on the calibrated node, and further how this marginal prior influences the posterior distribution in the light of the data (Brown and Smith 2017). We therefore also performed an analysis sampling only from the marginal prior (running an ‘empty’ analysis without data) in order to confirm that 1) the marginal prior distributions were similar to our user-specified calibration prior and 2) the data were informative enough to return posterior distributions different from the marginal priors. We did this by visually inspecting distributions in Tracer v1.6. We analyzed trace files in Tracer v1.6 to determine chain convergence and burn-in, and after discarding a burn-in of at least 20%, summarized 36,300–84,300 trees as maximum clade credibility (MCC) trees using LogCombiner v2.4.3 and TreeAnnotator v2.4.3.

Biogeographic Analyses

We inferred the biogeographic history of Crematogastrini using BioGeoBEARS (Matzke 2013). We focused hereby exclusively on the genus level to avoid any bias caused by species ranges being either overrepresented or absent from our phylogeny. While this is a rather coarse approach to biogeographic reconstruction, we believe this is appropriately conservative given our taxon sampling. We thus coded geographic distributions for each genus based on the distribution of all of its described members (see  Supp Table 3 [online only] (ixy013_suppl_supplementary_table_s1_s7.pdf)). Given the disagreements between classification and phylogeny for some genera, we assigned Ancyridris polyrhachioides Wheeler, Propodilobus pingorum DuBois, and the three Lasiomyrma Terayama & Yamane species to a more inclusive Lordomyrma Emery s.l. clade. Strongylognathus testaceus Schenck was further included within Tetramorium Mayr. We designated seven biogeographic areas (Neotropical, Nearctic, Palearctic, Afrotropical, Malagasy, Indomalayan, and Australasian [after Cox 2001]), and defined dispersal constraints based on the level of connectivity between these landmasses using slightly different rates for two distinct time periods from 0–50 Ma and 50–80 Ma (see  Supp Table 4 [online only] (ixy013_suppl_supplementary_table_s1_s7.pdf)). We prepared a modified chronogram from each of the 50-best analyses based on the concatenated (50-best-CP) and the species-tree (50-best-ST) topologies by pruning all but one representative species per genus, as well as the outgroups from the tree. We used these phylogenies, the genus-level distribution matrix, and the constraint matrix as input for two sets of BioGeoBEARS analyses (50-best-CP and 50-best-ST) following guidelines and tutorials available on the BioGeoBEARS PhyloWiki ( http://phylo.wikidot.com/biogeobears). For both sets of analyses, we implemented the standard two-parameter dispersal and extinction cladogenesis (DEC) model (Ree and Smith 2008), as well as the DEC+J version of the former model that incorporates founder-event speciation by assigning a separate probability parameter ‘j’ (Matzke 2014). We defined max_range_size = 7, since a few genera occupy all seven biogeographic areas. The starting value for the founder-event speciation parameter ‘j’ was specified with jstart = 0.0001. We compared the fit of these models with a log-likelihood ratio test and repeated each analysis twice, obtaining identical results each time.

Diversification Rate Estimation

We estimated whether shifts in diversification rates occurred over time in the evolution of Crematogastrini using BAMM v2.5 (Rabosky et al. 2013, Rabosky 2014) and the associated R package BAMMtools v2.1.5 (Rabosky et al. 2014). We accounted for incomplete sampling by using clade-specific sampling probabilities. For these we assembled a species-richness matrix listing 1) the diversity of described species and 2) the diversity of described species and subspecies for each genus (see  Supp Table 5 [online only] (ixy013_suppl_supplementary_table_s1_s7.pdf)). We then accordingly created two clade-specific sampling probability matrices by assigning each tip (i.e., a species) in our phylogeny to a clade and assigned this clade a sampling fraction calculated as the number of sampled species divided by number of described species (or species + subspecies) in that clade (see  Supp Table 5 [online only] (ixy013_suppl_supplementary_table_s1_s7.pdf)). Most species were assigned to a clade equivalent to their respective genus, except we again used a more inclusive definition for the Lordomyrma s.l. clade (see above) and Tetramorium. We used the two UCE-50-best chronograms based on the concatenated and coalescent species-tree topologies plus these two sets of sampling fractions to perform four different sets of BAMM analyses. We configured our analyses using the guidelines in the BAMM documentation ( http://bamm-project.org/). We used the function ‘setBAMMpriors’ within BAMMtools to obtain appropriate priors for speciation-extinction analyses, and the expected number of shifts was left at the default value (=1). Our runs included four MCMC chains with a length of 50 million generations, sampling every 10,000 generations. We discarded a burn-in of 10% and confirmed that effective sample size values were appropriate (1,400–2,400). Results were analyzed and plotted with various functions in BAMMtools as follows. We used computeBayesFactors to identify the best-supported model of rate shifts in our data. We plotted the maximum a posteriori probability (MAP) shift configuration (getBestShiftConfiguration), the 95% credible set of distinct rate shift configurations (credibleShiftSet), and a ‘phylorate’ graph showing mean marginal posterior density of speciation rates (plot.bammdata).

Nesting Preference

To estimate the ancestral state of nesting preference in the Crematogastrini, we scored each genus by its members' combined nesting habits as either ground-dwelling, arboreal, or generalist (i.e., occupying both arboreal and ground habitats). Sixteen genera were thus scored as being exclusively arboreal, 35 as ground-dwelling, and five as generalists ( Supp Table 3 [online only] (ixy013_suppl_supplementary_table_s1_s7.pdf)), based on information we gathered from the literature, specimen records on AntWeb, and personal collecting experience. We used the rayDISC function in the R package corHMM ( https://www.R-project.org/) to reconstruct ancestral states for these trait categories. We performed reconstructions under the ‘equal rates’ (ER) model and the ‘all rates different’ (ARD) model and compared the fit of these models by performing a likelihood ratio test on the resulting –lnL scores as (1-pchisq(Δ lnL, df)). All reconstructions were performed using both the 50-best-CP and the 50-best-ST chronogram.

UCE Data Characteristics

Due to the large age range of our specimens, the concentrations of our DNA extractions ranged from <0.05 to 27.1 ng/µl (average 2.1 ng/µl). Our post-library preparation DNA concentrations were equally variable, ranging from 0.3 to 105 ng/µl (average 32.5 ng/µl). From the enriched libraries, we were able to generate 13,667–5,845,581 reads per taxon (average: 2,279,011), which Trinity assembled into 136–228,374 contigs (average: 59,518) of average length 235–606 bp per taxon (average across taxa: 347 bp) and 4–28.1× coverage (average: 9.4×). The assembled contigs matched 103–2399 UCE loci (average: 2200 loci after removal of duplicates in Phyluce) with an average length of 229–1049 bp per taxon (average across taxa: 654 bp) and 4.7–123.2× coverage (average: 42.1×). Between 0.13 and 0.94 reads (average: 0.42) were on target. The overall lowest results were obtained from a Lordomyrma epinotalis specimen that had already undergone a nondestructive DNA extraction procedure for a different study (Lucky and Sarnat 2010); it was re-extracted by us due to the lack of other available material. Details on the library concentrations and UCE capture values can be found in  Supp Table 6 (online only) (ixy013_suppl_supplementary_table_s1_s7.pdf). We calculated several descriptive statistics from our complete alignment (90% matrix, 1,763 loci), as well as for the data subsets used in the dating analysis, which are listed in  Supp Table 7 (online only) (ixy013_suppl_supplementary_table_s1_s7.pdf). In brief, our complete data matrix had a length of 913,718 bp with 18% missing data, 59% variable sites, and 49% parsimony-informative characters. The GC heterogeneity score across taxa in the alignment is relatively low at 1.2%, and mean GC content is moderate at 44%.

Phylogeny of Crematogastrini

We inferred phylogenetic relationships within the tribe Crematogastrini using a concatenated ML approach (both partitioned by locus and unpartitioned), as well as by reconstructing a species tree estimated from 1,763 UCE gene trees (binned into 883 ‘supergenes’). Both methods returned tree topologies that were overall well supported and largely congruent (Figs. 1 and 2), with the few exceptions described below. In all analyses, Rostromyrmex pasohensis Rosciszewski is estimated as the sister lineage to all remaining Crematogastrini. Based on the concatenated ML tree, we further identify 10 main clades within Crematogastrini that have good or full support (BS: 79–100), here treated as informal genus-groups. Eight of these groups are recovered, with maximum support, in the species-tree topology but two of them become paraphyletic in that analysis (Fig. 2).

Fig. 1.

Phylogeny of Crematogastrini, estimated by concatenated ML analysis. We performed best-tree and bootstrap searches (N = 100) in RAxML v8.2.7. Both panels show the best ML tree resulting from analyses of a concatenated data matrix divided per UCE locus into 1,763 partitions. (A) Phylogram showing branch lengths and emphasizing generic relationships. (B) Cladogram with BS values indicated: white squares represent 100% BS; red squares indicate the support for that node. The 10 major genus-groups are indicated; the compositions of the two genus-groups highlighted in red conflict with those resulting from the ASTRAL-II analysis presented in Fig. 2.

Fig. 2.

Cladogram estimated for Crematogastrini, by ASTRAL-II species-tree analysis. We reconstructed a species tree from 1,763 UCE gene trees based on weighted statistical binning (Bayzid et al. 2015) and 883 supergenes. Only local posterior probabilities (LPP) > 0.7 are shown. Note that ASTRAL-II support values are branch support values that measure the support for a quadripartition, not a bipartition.The 10 major genus-groups are indicated; clade names in red highlight conflicts regarding the composition of these clades with respect to Fig. 1.

Aside from the alternative topologies for the Carebara genus-group and the Mayriella genus-group, concatenated and species-tree approaches show only one other major disagreement: in the concatenated topology, the Crematogaster genus-group is sister to the Podomyrma genus-group with fairly high support (BS: 97; Fig. 1B), whereas the species-tree analysis estimates a ladderized topology and does not support a sister-group relationship between these two clades (Fig. 2). We also performed an unpartitioned analysis on the concatenated data set to test whether our partitioning strategy had an influence on the results. This unpartitioned ML analysis returned a tree very similar to the above partitioned analyses, with the notable exception that Cardiocondyla and Oxymyrmex jump outside the Carebara genus-group as a separate lineage, sister to the remaining Crematogastrini (but exclusive of Rostromyrmex) ( Supp Fig. 1 [online only] (ixy013_suppl_supplementary_figure_s1_s12.pdf)).

Since Rostromyrmex, Cardiocondyla, and Ocymyrmex are subtended by relatively long branches, we explored the possibility that LBA affects the placement of these taxa, as well as the influence of these taxa on tree topology in general, by sequentially excluding them from the analysis. The results of these exclusion tests are depicted in schematic format in  Supp Fig. 2 (online only) (ixy013_suppl_supplementary_figure_s1_s12.pdf). Excluding Rostromyrmex from the analysis had no effect on topology or support with regard to the remaining ingroup taxa ( Supp Fig. 2A [online only] (ixy013_suppl_supplementary_figure_s1_s12.pdf)), nor did the exclusion of Cardiocondyla ( Supp Fig. 2C [online only] (ixy013_suppl_supplementary_figure_s1_s12.pdf)) or Ocymyrmex and Cardiocondyla at the same time ( Supp Fig. 2D [online only] (ixy013_suppl_supplementary_figure_s1_s12.pdf)). In contrast, excluding Ocymyrmex from the analysis eroded support for the Cataulacus genus-group clade, resulting in a position of Cardiocondyla outside of this clade, grouping as sister to the remaining Crematogastrini (exclusive of Rostromyrmex). This alternative topology, however, received little support (bootstrap proportion: 44;  Supp Fig. 2B [online only] (ixy013_suppl_supplementary_figure_s1_s12.pdf)).

Timing of the Evolution of Crematogastrini

Using BEAST2 and a reduced data set of 50 UCE loci, we estimated the timescale of the evolution of Crematogastrini. We estimated ages from 50 loci with the best average bootstrap score (50-best) for both the concatenated and species-tree topologies (concatenated-partitioned: 50-best-CP; species tree: 50-best-ST). We also compared these age estimates to divergence ages estimated from two sets of randomly chosen loci (50-random-1, 50-random-2) to assure that the 50-best data set did not result in biased divergence ages. Table 1 compares divergence estimates from these four analyses for the main clades within Crematogastrini, as well as for all genera for which crown-group estimates were obtained. Complete chronograms with age estimates and 95% highest posterior density intervals can be found in the supplementary material ( Supp Figs. 3–6 [online only] (ixy013_suppl_supplementary_figure_s1_s12.pdf)). Overall, we received highly congruent age estimates from the four analyses, differing between <1 and 5.7 Ma (Table 1). The only exception is the crown-group age of Dacetinops, which was estimated ca. 20 Ma younger from the 50-random-1 and 50-random-2 data sets than from the 50-best-CP and 50-best-ST data sets. Differences between the concatenated and species-tree topologies did not appear to affect age estimates significantly. For example, the Lordomyrma genus-group, while estimated in a different position in the two analyses, is estimated with a crown age of 38.3 and 37 Ma in the 50-best-CP and 50-best-ST analysis, respectively.

Crown-group Crematogastrini are estimated in our analyses to have evolved ∼75.3–79.4 Ma in the late Cretaceous. This is remarkably close to the estimate of 70.8 Ma obtained in the larger analysis of the subfamily Myrmicinae by Ward et al. (2015)—espe-cially considering that study did not include Rostromyrmex, the earliest-branching taxon in our analysis. Diversification into the 10 main clades that we recognize here took place mainly in the late Cretaceous to early Eocene between 50 and 70 Ma, with the exception of the Lordomyrma genus-group, which is somewhat younger with an estimated age of 37–40 Ma. The Oligocene and Miocene were important settings in generating the presently recognized generic diversity of Crematogastrini, as the majority of genera appear to have evolved between 8 and 31 Ma. Notably older, with Eocene crown-age estimates, are Carebara (46.8 Ma), Nesomyrmex (40.5 Ma), and Crematogaster (39.3 Ma). A few genera are estimated as being of very young crown-group age, such as Ocymyrmex (5.8 Ma), Atopomyrmex (3.1 Ma), Lophomyrmex (4 Ma), Rotastruma (4.4 Ma), and Romblonella (1.8 Ma). These estimates must be viewed with caution as they could be influenced by incomplete taxon sampling.

Biogeography

We used BioGeoBEARS (Matzke 2013) to infer the early biogeographic history of Crematogastrini from the 50-best-CP and 50-best-ST chronograms. We compared the fit of the standard two-parameter DEC model (Ree and Smith 2008) with the DEC+J model incorporating founder-event speciation. The DEC+J model did not provide a significantly better fit in either of the two alternative analyses; therefore, we focus on discussing results estimated from the standard DEC model. Biogeographic reconstructions are highly congruent between the two alternative topologies (Fig. 3A;  Supp Fig. 7A [online only] (ixy013_suppl_supplementary_figure_s1_s12.pdf)). The ancestral range of Crematogastrini is estimated in the Afrotropical and Indomalayan regions, which also remain the two dominant ancestral distributions during much of the early evolution of the tribe. However, in the very early evolution of Crematogastrini the affinity with the Indomalayan region is indicated only by the phylogenetic position of the enigmatic Malaysian endemic and monotypic genus Rostromyrmex. The Cataulacus genus-group and the Carebara genus-group are inferred to have evolved mainly in the Afrotropical region, with a possible dispersal to the Neotropics by the ancestor of Nesomyrmex and Xenomyrmex between 62 and 68 Ma, and a dispersal to the Indomalayan region by Rhopalomastix sometime in the Eocene (after 45 Ma). This scenario would imply dispersal back to the Afrotropics by one Nesomyrmex lineage. All remaining eight clades are estimated to share a common ancestor restricted to the Indomalayan region. Further lineage diversification appears to have taken place exclusively in this region until ca. 65 Ma or later, when an ancestor of Adlerzia, Recurvidris, and Lophomyrmex dispersed to Australasia. Other notable dispersal events in the early history of Crematogastrini are inferred within the Podomyrma genus-group. Here, one dispersal of an ancestor of the four Malagasy endemic genera Eutetramorium, Royidris, Malagidris, and Vitsika to the Malagasy region is supported between 50 and 60 Ma, as is another dispersal event around the same time in the sister lineage comprised of Liomyrmex, Metapone, Huberia, and Podomyrma. For the latter, while the ancestral reconstructions remain ambiguous, an initial dispersal from the Indomalayan to the Australian region is probably most likely: Podomyrma and Huberia are Australasian taxa (Huberia is a New Zealand endemic; Podomyrma occurs in Australia and New Guinea), Liomyrmex occurs in southeast Asia and New Guinea, and only Metapone extends to the African and Malagasy regions. Again prior to ca. 50 Ma, a dispersal to the Paleotropical region is estimated within the Formicoxenus genus-group, followed by a subsequent dispersal to the Nearctic. The early evolution of the remaining clades (Vollenhovia genus-group, Crematogaster genus-group, Mayriella genus-group, Lordomyrma genus-group, Myrmecina genus-group, and Paratopula genus-group) was found to be heavily centered in the Indomalayan region. Further dispersal events explaining the broader present-day distributions of some crematogastrine genera are therefore likely associated with species-level diversification events not captured by our biogeographic study.

Table 1.

Divergence age estimates for major clades (bolded) and genera within Crematogastrini

Evolution of Nesting Preference

The ER model was estimated to provide the best fit to our nesting preference data, and tree topology had only a minor effect on state reconstructions. We found overwhelming support for a ground-dwelling ancestor of Crematogastrini, which is probably unsurprising given that the majority of the genera within the tribe are ground-dwelling. A more interesting result is the good support for four lineages of exclusively arboreal genera: 1) the subclade within the Cataulacus genus-group including Xenomyrmex, Nesomyrmex, Cataulacus, Terataner, and Atopomyrmex; 2) Rhopalomastix and Melissotarsus (part of the Carebara genus-group); 3) the Paratopula genus-group; and 4) Vombisidris and Gauromyrmex (part of the Formicoxenus genus-group). There is further weak support (PP = 0.4–0.5, Fig. 3B;  Supp Fig. 7B [online only] (ixy013_suppl_supplementary_figure_s1_s12.pdf)) for an arboreal nesting preference in the most recent common ancestor (MRCA) of the Formicoxenus and Paratopula genus-groups. Otherwise, exclusive arboreality has evolved independently only in three other genera, Dilobocondyla, Podomyrma, and Vitsika. Four or five genera have adapted to a generalist lifestyle from a ground-dwelling ancestor: Crematogaster, Tetramorium, Trichomyrmex, Leptothorax, and Temnothorax, with the ancestral state of the latter being somewhat ambiguous.

Fig. 3.

The evolution of biogeography (A) and nesting preference (B) in Crematogastrini. Ancestral reconstructions based on a modified chronogram from the 50-best concatenated-partitioned analysis in BEAST2, from which all but one representative species per genus, as well as the outgroups, have been pruned from the tree. (A) Biogeographic reconstructions with BioGeoBEARS under the DEC model. N = Nearctic, T = Neotropical, P = Palearctic, E = Afrotropical, M = Malagasy, O = Indomalayan, and A = Australasian. (B) Nesting preference reconstructed under the ER model with rayDISC in corHMM. Black = arboreal nesting; white = ground nesting. See SuppTable 3 (online only) for geographic distributions and trait data for each genus.

Genus-Level Diversification of Crematogastrini

We estimated rates of diversification in Crematogastrini and investigated whether significant rate shifts have occurred over time, using the 50-best-CP and 50-best-ST topologies and two sets of sampling fractions based on species estimates that either excluded or included subspecies. We plotted the mean marginal posterior density of speciation rates for the four different analyses in Fig. 4A and B and  Supp Fig. 8A and B (online only) (ixy013_suppl_supplementary_figure_s1_s12.pdf). These phylorate graphs show an increase in estimated rate in all analyses along the branches leading to five genera: Crematogaster, Meranoplus, Tetramorium, Temnothorax, and Ocymyrmex. The highest rate increase occurs along the branch leading to Tetramorium. The two analyses including subspecies further estimate a higher rate increase along the branch leading to Crematogaster than analyses excluding subspecies. The best-supported rate shift configuration estimated three shifts occurring throughout the diversification of Crematogastrini, regardless of which topology and species estimate was used (Fig. 4C and D;  Supp Fig. 8C and D [online only] (ixy013_suppl_supplementary_figure_s1_s12.pdf)). The four analyses all agree on the position of two of these rate shifts: one is placed on the branch leading to Tetramorium and one on the branch leading to Temnothorax. The third shift is placed on either the branch leading to Crematogaster (50-best-CP analysis including subspecies, 50-best-ST analyses excluding/including subspecies; Fig. 4D;  Supp Fig. 8C and D [online only] (ixy013_suppl_supplementary_figure_s1_s12.pdf)) or the branch leading to the Crematogaster–Meranoplus clade (50-best-CP analysis excluding subspecies, Fig. 4C). However, the frequency (=posterior probability support) for the best shift configuration was low in all four cases (f = 0.2–0.27), and several other slightly different credible shift configurations were found ( Supp Figs. 9–12 [online only] (ixy013_suppl_supplementary_figure_s1_s12.pdf)).

Fig. 4.

Diversification of Crematogastrini. We performed BAMM analyses on the chronogram resulting from the 50-best concatenated-partitioned BEAST2 analysis using clade-specific sampling probabilities to account for incomplete sampling based either on species estimates only or including species and subspecies. Panels show (A) mean phylorate plots based on species only; (B) mean phylorate plots based on species and subspecies; (C) best shift configuration based on species only; (D) best shift configuration based on species and subspecies. SuppTable 5 (online only) lists diversity estimates per genus including and excluding subspecies.

Biogeography

Divergence dating and biogeographic inferences suggest a Paleotropical origin for the Crematogastrini at approximately 75 Ma. The 10 major lineages recognized here arose over a relatively limited time interval (50–70 Ma) in the Afrotropical and Indomalayan regions. This was a period of warm, equable climate (Zachos et al. 2008, Archibald et al. 2013) when many of the major angiosperm clades also originated, especially among eudicots (Tank et al. 2015). Most extant genera of Crematogastrini have estimated crown ages of 10–40 Ma, i.e., they appear to have originated after the early Eocene climatic optimum (Zachos et al. 2001) under generally cooler and drier conditions, and again predominantly in the Old World. A few genera, such as Carebara, Crematogaster (Blaimer 2012), and Temnothorax (Prebus 2017), have infiltrated the New World, but it is striking how the main theater of crematogastrine evolution has been the Paleotropics, especially the Indomalayan region. Of the 10 genus-group lineages we recognize in this study, the majority (eight) are inferred to have an Indomalayan origin, and only two, the Cataulacus and Carebara genus-groups, are reconstructed having an ancestral African distribution. Interestingly, the latter two are the earliest-branching lineages. However, given the enigmatic position of the Indomalayan genus Rostromyrmex as sister to the remaining Crematogastrini, the biogeographic origin of the tribe remains somewhat uncertain.

Within the Podomyrma genus-group, several dispersal events are inferred away from the Indomalayan region. First, a clade of four Malagasy endemic genera is inferred to have dispersed to that island before 50 Ma. This is an interesting result, since the only other time-calibrated study investigating dispersal of ants to Madagascar had suggested a much more recent origin for the genus Crematogaster (several dispersals to Madagascar after ∼25 Ma; Blaimer 2012). Could the ancestor of these endemic groups have been the first ant to arrive in Madagascar? Neither of these four genera are currently very diverse (3–16 species). One might expect to see a burst of diversification coinciding with increased ecological opportunity if most ecological niches for ants had been freely available at the time of arrival. However, competition by more recently arriving lineages followed by extinction could have erased traces of an initial radiation of this Malagasy endemic ant clade. Some other reconstructions in the Podomyrma genus-group appear less biologically plausible given the present distribution of taxa and are probably a result of the genus-level coding used for distributions. For example, the MRCA of the clade consisting of Podomyrma, Huberia, Liomyrmex, and Metapone most likely had an Indomalayan and Australasian distribution (not Malagasy and Indomalayan), given the present-day distribution of these genera, and a dispersal to the Malagasy and Afrotropical region probably only happened much later within Metapone than our reconstructions suggest. Likewise, the inferred dispersal event of the MRCA of Adlerzia, Recurvidris, and Lophomyrmex to the Australasian region could be an artifact of assigning the Indomalayan and Australasian region equal weight in the distribution of Recurvidris and Lophomyrmex (by coding ‘OA’). These genera are predominantly Indomalayan in distribution and only a few species cross Wallace’s line into Australasia. Instead it is conceivable that only one dispersal event occurred along the stem lineage of Adlerzia, followed by separate recent dispersals within the other two genera.

Nesting Habitat Evolution

Ants successfully occupy a breadth of habitats and their vertical layers (canopy, leaf litter, and soil), yet these vary tremendously in microclimate and properties of nest sites. A significant distinction can usually be made between a community of arboreally nesting species and a community of ground-nesting species (living in the soil and leaf litter). A few larger ant genera have adapted to utilize both these vertical habitat layers (such as Crematogaster, Tetramorium, Temnothorax, and Leptothorax among Crematogastrini) but by and large we see fidelity in the choice for either one or the other. The question of whether the MRCA of all ants was arboreal or ground-dwelling, or in which habitat crown-group ants have evolved, has thus always piqued the curiosity of myrmecologists. Earlier hypotheses suggested the leaf litter was the cradle of ant evolution (Wilson and Hölldobler 2005, Perrichot et al. 2007). A subterranean, soil-dwelling ancestor was supported by the first formal analysis of this question, and the subfamily Myrmicinae was further reconstructed as ancestrally surface-dwelling (i.e., leaf litter) (Lucky et al. 2013). Perhaps unsurprisingly given these earlier results, we infer ground-dwelling (i.e., leaf litter and soil-dwelling) to be the ancestral condition in Crematogastrini. However, we discovered two larger clades within Crematogastrini that apparently have switched in whole to arboreal nesting: a clade within the Cataulacus genus-group comprising Atopomyrmex, Cataulacus, Terataner, Nesomyrmex, and Xenomyrmex, as well as the entire Paratopula genus-group. For the latter, it is even possible (based on our reconstructions) that the switch occurred in the ancestor of the Paratopula and Formicoxenus genus-groups. Another switch to arboreal nesting occurred in the lineage comprising the enigmatic genera Melissotarsus and Rhopalomastix, which possess specialized silk-producing glands (although for Rhopalomastix the extent of silk production and use is not yet entirely confirmed; Hölldobler et al. 2014). We see a few switches to arboreality in other phylogenetically isolated genera (Podomyrma, Vitsika, and Dilobocondyla), but clearly overall arboreality has evolved in a nonrandom fashion in Crematogastrini. No switches from an exclusively arboreal nesting habit back to ground-nesting are supported (although some uncertainty remains within the Formicoxenus-group). This may indicate benefits of an adaptation to this habitat layer that prevent the retention of a generalist or reversal to a ground-dwelling lifestyle. More detailed species-level analyses estimating ancestral nest-site evolution within the ‘generalist genera’ (Crematogaster, Tetramorium, Temnothorax, and Leptothorax) are necessary, which in turn could lead to slightly different conclusions at the genus level.

Diversification of Crematogastrini

The subfamily Myrmicinae is the most species-rich and diverse lineage of ants in the world. A majority of this diversity is represented by the tribe Crematogastrini (together with the Attini), which contains three of the largest, most species-rich genera within the Myrmicinae: Tetramorium, Crematogaster, and Temnothorax. Three diversification rate shifts were inferred by our study and these shifts are associated unsurprisingly with these same three genera, which are also among the geographically most widespread taxa within Crematogastrini. A moderate diversification rate increase is further estimated for Meranoplus and Ocymyrmex, but with no distinct shifts supported. These two genera are fairly small, with 90 described species for Meranoplus and 34 species for Ocymyrmex. One explanation is that Ocymyrmex is estimated with a very young age (5.8 Ma), an estimate that could have been influenced by incomplete taxon sampling in our dating analyses (see above). We suspect that the evolutionary rate for this lineage is overestimated by our BAMM analyses. Meranoplus is the sister lineage to Crematogaster, and in one of the analyses (50-best-CP analysis excluding subspecies; Fig. 4C), a rate shift is placed on the branch subtending both Meranoplus and Crematogaster. Possibly Meranoplus shares some (but not all) of the characteristics with Crematogaster that led to increased diversification in this clade. Interestingly, in a larger-scale analysis estimating diversification rate across the ants (Blanchard and Moreau 2016), a positive shift in diversification rate is also placed on the branch leading to the Crematogaster and Meranoplus clade (rather than just Crematogaster). Further rate shifts in this study are placed in very similar locations on the Crematogastrini subtree; for example, as in our analysis, a positive shift in diversification rate for Tetramorium is recovered (Blanchard and Moreau 2016). However, the positive shift placed in our estimations on the Temnothorax branch is located in that study on the branch leading to a more inclusive clade of Temnothorax, Leptothorax, Formicoxenus, and Harpagoxenus, and is further estimated to be negative, resulting in a decreased diversification rate for the clade. That study also recovers another negative shift leading to a clade comprised of Proatta, Dacatria, and Tetheamyrma (Blanchard and Moreau 2016). Presumably these differences in results compared to our study are influenced by taxon sampling, divergence dating, and the different methods used to estimate diversification rates (Medusa was used in Blanchard and Moreau 2016). The accuracy of BAMM in estimating diversification rates and rate shifts under certain circumstances has been recently debated. Moore et al. (2016) raised concerns about the sensitivity of BAMM analyses to the selected rate shift prior, which were convincingly countered by Rabosky et al. (2017). Additional criticism has been raised regarding a tendency to overestimate diversification rates in smaller clades, thereby resulting in a potential underestimation of rate shifts overall (Meyer and Wiens 2018). We conclude that our results are most likely not affected by potential erroneous diversification rate estimates, given that we recovered statistically significant shifts for the three most species-rich genera within Crematogastrini only—a result that essentially confirms observations based on taxonomic species diversity.

Positive shifts in diversification rate have usually been associated with either increased ecological opportunity, e.g., dispersal to and colonization of a new environment, or with the evolution of a key innovation, such as a novel trait (be it morphological, physiological, or genetic) that confers a competitive advantage to an organism or allows it to expand into a previously inaccessible environmental niche space. For example, turtle ants, genus Cephalotes, show an increase in diversification rate after colonizing a new habitat in the Chacoan region of South America (Price et al. 2014, Price et al. 2016). This positive rate shift was associated with phenotypic diversification in turtle ants (Price et al. 2016). Likewise, diving beetles (genus Exolina) were found to show a diversification rate shift after dispersing to Melanesia from Australia (Toussaint et al. 2015). Our results did not enable us to discern an association of the diversification rate shifts leading to Crematogaster and Tetramorium with a biogeographic dispersal event or an evolutionarily significant geological time period. However, potential biogeographic dispersal events within these genera may be responsible for these patterns, which would be masked by our generic-level approach to estimating diversification rates. With respect to the rate shift estimated for Temnothorax, if this shift would be more correctly located on the branch leading to the larger clade that includes Leptothorax, Formicoxenus, and Harpagoxenus (as estimated in Blanchard and Moreau 2016), a possible association with a dispersal of this clade from the Indomalayan to the Palearctic and Nearctic regions and a burst of diversification in these temperate regions is conceivable. Beyond a biogeographic dispersal, the colonization of additional nesting space could have led to increased diversification in these three groups. Crematogaster, Temnothorax, and Tetramorium have all evolved to nest both in ground and in arboreal habitats from an ancestrally ground-nesting state, and it seems plausible that the exploitation of the arboreal habitat has led to increased species diversification. For Crematogaster, morphological innovations also have been proposed to have led to the success of this genus. Synapomorphies for the genus are the dorsal attachment of the post-petiole to the fourth abdominal segment (first segment of gaster) and the absence of a dorsal petiolar node. These features give the ants the ability to flex the gaster forwards over the mesosoma while the petiole can be pressed tightly against the propodeum for protection (Buren 1959). Crematogaster workers also use this posture in defense paired with another unique feature: their spatulate sting applies venom topically instead of piercing the skin of the attacker (Marlier et al. 2004). Unfortunately, while hypotheses of diversification associated with dispersal and ecological opportunity could be tested with extended species-level phylogenetic and biological data, hypotheses regarding the evolutionary significance and influence of innovative features are difficult to test and must remain somewhat speculative.

In summary, we here present an improved genus-level phylogeny for the tribe Crematogastrini that provides the basis for recognizing 10 genus-group lineages. Our estimates show that these lineages have evolved over a relatively short time frame (50–70 Ma) from a Paleotropical ancestor. The position of the monotypic enigmatic Indomalayan genus Rostromyrmex as the sister lineage to all remaining lineages within the tribe prevents us from pinpointing a specific biogeographic region that acted as the cradle for the diversification of the tribe. However, we have strong support for an early diversification of Crematogastrini in the Afrotropical region within the Cardiocondyla and Carebara genus-groups, and then for a shift in emphasis toward the Indomalayan region as the center for crematogastrine diversification. We discerned several shifts in diversification rates related to the evolution of large, widespread genera, but were not able to associate these definitively with key evolutionary events such as dispersal or novel morphological characteristics. Arboreal habitats have been successfully colonized by only a few clades within Crematogastrini, but it appears that members of the tribe have remained faithful to their lofty nest sites once they adapted to utilize this habitat. Our genus-level study sets the scene for more detailed investigations of diversification and evolution in the delineated genus-group lineages using species-level phylogenies and biological data.