Access the full text.
Sign up today, get DeepDyve free for 14 days.
Loading next page...
/lp/ou_press/when-darwin-s-special-difficulty-promotes-diversification-in-insects-Dt3DDKLIWi
- Publisher
- Oxford University Press
- Copyright
- © The Author(s) 2018. Published by Oxford University Press, on behalf of the Society of Systematic Biologists. All rights reserved. For Permissions, please email: journals.permissions@oup.com
- ISSN
- 1063-5157
- eISSN
- 1076-836X
- DOI
- 10.1093/sysbio/syy014
- pmid
- 29481653
- Publisher site
- See Article on Publisher Site
Abstract
Abstract Eusociality, Darwin’s special difficulty, has been widely investigated but remains a topic of great debate in organismal biology. Eusocial species challenge existing theories, and the impact of highly integrated societies on diversification dynamics is controversial with opposing assertions and hypotheses in the literature. Here, using phylogenetic approaches in termites—the first group that has evolved eusociality—we assessed the fundamental prediction that eusocial lineages have higher diversification rates than non-eusocial clades. We found multiple lines of evidence that eusociality provided higher diversification as compared to non-eusociality. This is particularly exacerbated for eusocial species with “true” workers as compared to species with “false” workers. Because most species with “true” workers have an entirely prokaryotic microbiota, the latter feature is also related to higher diversification rates, but it should be investigated further, notably in relation to angiosperm diversification. Overall, this study suggests that societies with “true” workers are not only more successful at ecological timescales but also over millions of years, which further implies that both organism- and species-level traits act on species selection. Aggregate trait, emergent trait, eusocial insects, macroevolution, species selection, termites Eusociality is the most integrated organizational level of animal sociality and is defined by cooperative brood care, overlapping generations within a colony of adults, and division of labor into reproductive and non-reproductive groups (Wilson and Hölldobler 2005). Although the origin of eusociality is still widely investigated and debated (e.g., Nowak et al. 2010; Nonacs 2011; Rousset and Lion 2011; Johnstone et al. 2012), this rare but widespread behavioral characteristic has independently appeared in insects (at least 7–20 times depending on the definition of eusociality), crustaceans (at least thrice), and mammals (at least twice) (Fig. 1; Bourke 2011). Within insects, eusocial species have colonies with caste differences: queens and reproductive males take the roles of sole reproducers, while soldiers defending the colony and workers foraging and maintaining resources contribute to creating a living situation favorable for the colony (Wilson 1971). In his theory of natural selection, Darwin (1859) famously portrayed the eusocial insects as a “special difficulty, which at first appeared to [him] insuperable and actually fatal to [his] theory”. Indeed, eusociality presents an apparent paradox: if adaptive evolution is mediated by differential reproduction of individuals on which natural selection acts, how can individuals incapable of passing on their genes evolve and persist? Figure 1. View largeDownload slide Origin and timeline of eusocial lineages for mammals, crustaceans, and insects (bees, wasps, ants, and termites). Bars represent temporal ranges with solid lines based on the fossil record, whereas broken lines indicate possible eusociality supported by molecular dating studies (ants: Moreau et al. 2006; bees: Cardinal and Danforth 2011, wasps: Hines et al. 2007, termites: Ware et al. 2010; Legendre et al. 2015b) but not supported by fossil evidence. Figure 1. View largeDownload slide Origin and timeline of eusocial lineages for mammals, crustaceans, and insects (bees, wasps, ants, and termites). Bars represent temporal ranges with solid lines based on the fossil record, whereas broken lines indicate possible eusociality supported by molecular dating studies (ants: Moreau et al. 2006; bees: Cardinal and Danforth 2011, wasps: Hines et al. 2007, termites: Ware et al. 2010; Legendre et al. 2015b) but not supported by fossil evidence. The multilevel selection theory could settle this paradox. This theory has a long and troubled history (Okasha 2006) but is increasingly accepted and supported, especially at the species level (Jablonski 2008). Species selection results from differences in speciation and extinction rates among species within a clade, and these differences are not merely related to the fitness at the level of individuals (Jablonski 2008; Rabosky and McCune 2010). Data and theory suggest that both organism-level and species-level traits (called aggregate and emergent traits, respectively) act on species selection. Aggregate traits are characteristics of individual organisms (e.g., body size), whereas emergent traits are species properties that are not reducible to organismal traits such as eusociality. Extending this vocabulary, authors often refer to the emergent fitness of species within clades when discussing species selection (Jablonski 2008). At the macroevolutionary timescale, the origin and evolution of eusociality are increasingly well-understood due to recent phylogenetic analyzes, along with estimates of divergence times and diversification rates notably for arthropod groups (ants: Moreau et al. 2006; Pie and Feitosa 2016, bees: Cardinal and Danforth 2011, beetles: Farrell et al. 2001, wasps: Hines et al. 2007, spiders: Agnarsson et al. 2006, termites: Ware et al. 2010). These studies have investigated the number of origins of eusociality and the putative factors behind the extraordinary diversification of eusocial groups (e.g., the Cretaceous radiation of Angiosperms: Moreau et al. 2006; Cardinal and Danforth 2013). But another challenge surrounds eusocial organism evolution: the impact of highly integrated societies on the dynamics of diversification in these organisms (Davis et al. 2009; Ware et al. 2010). This topic has received less attention probably due to a lack of methods to handle large phylogenies that include both non-eusocial and eusocial lineages. Although the role of eusociality has not been formally studied with state-of-the-art diversification methods, the literature is filled with opposing assertions and hypotheses. In terms of ecological and evolutionary success, modern eusocial insects can either be seen as species-rich groups such as ants, bees, or termites (Wilson 1990; Chapman and Bourke 2001) or as a group with lower species diversity than non-eusocial relatives such as ambrosia beetles, or poplar spiral gall aphids (Wilson 1992). Although the success of eusocial groups could naively lead us to think that eusociality was a trigger of diversification, Wilson (1992) argued that social insects (contrary to social vertebrates) show lower diversification (because they are “less speciose”) than non-eusocial insects. Based on studies of population genetics and molecular evolution, sociality reduces effective population size and possibly reduces levels of genetic variation, which would in turn increase the risk of extinction of social organisms (Lanfear et al. 2014). Hence, eusociality could indeed be seen as an impediment to diversification. The fossil record shows that taxa went extinct within eusocial groups (ants: LaPolla et al. 2013, termites: Krishna et al. 2010) but, so far, eusocial insect groups have been successful at the macroevolutionary scale. For the most part, they have persisted since the origins of these highly integrated behaviors, leading Wilson (1990) to state that “eusociality conveys evolutionary long life to social insects”. In reality, several emergent traits influence both extinction and speciation rates, resulting in unpredictable trade-offs (Jablonski 2008). The ensuing questions “does eusociality confer higher diversification rates than non-eusocial species?” and “what is the dynamic of diversification of social insects and what factor could explain their success and potential shift in diversification rates?” are long-standing research topics (Darwin 1859; Wilson 1992; Wilson and Hölldobler 2005). Biotic and abiotic factors are important in diversification processes (Benton 2009) and range from climatic/geologic events (e.g., Condamine et al. 2013), to competition (e.g., Liow et al. 2015; Silvestro et al. 2015a), and key innovations (e.g., Rainford et al. 2014; Sánchez-García and Matheny 2016). Eusociality can be classified as a key innovation but has surprisingly been understudied when one considers the ecological success that eusocial insect species have experienced (Wilson 1971, 1992; Chapman and Bourke 2001). This ecological dominance is mainly explained by one property of eusociality: the reproductive division of labor. In other words, some individuals do not reproduce (i.e., they belong to sterile castes), and this has been qualified as the most important feature of eusociality (Wilson 1990). Within insects, termites (Dictyoptera, order Blattodea, infraorder Isoptera) contain only eusocial species, and it is thought that they were the first modern eusocial animals to evolve, sometime in the Late Jurassic (ca. 155 Ma) (Legendre et al. 2015b; Engel et al. 2016). A termite colony is usually differentiated into reproductive, worker and soldier castes (Fig. 2a). The latter two castes, which are sterile and perform highly specialized tasks, derive from different stages depending on the species (Noirot and Pasteels 1987). Whereas soldiers likely have a single origin, the evolution of workers is more controversial (Watson and Sewell 1981; Noirot and Pasteels 1987; Thompson et al. 2000; Grandcolas and D’Haese 2002; Inward 2007, Legendre et al. 2008, 2013; Roisin and Korb 2011). Importantly, two types of workers were distinguished early: Some species have a “true” worker caste, while others have a “false” worker caste (also called pseudergates; Fig. 2a). “True” workers differ from pseudergates notably because they diverge early and irreversibly from the imaginal line (Noirot and Pasteels 1987), while pseudergates do not. Societies with “true” workers have been qualified as more “socially” highly integrated, which might have consequences for their ecological and evolutionary success given that $$>80$$% of all termites have “true” workers (Korb 2009; Brune 2014). Figure 2. View largeDownload slide Eusociality in termites: life cycle (a) and b) evolutionary origin associated to putative key innovations. The life cycle differs for species with pseudergates and species with “true workers” (highlighted by the dashed box). The phylogeny of termites (Isoptera) illustrates that the origin of eusociality is nested within the radiation of cockroaches (termites and cockroaches together forming the Blattodea). The phylogeny also highlights important events in the evolution of the digestive symbiosis: the presence of cellulolytic flagellates, acquired by a common ancestor of termites and Cryptocercidae, is associated to their wood-feeding lifestyle; the loss of flagellates in Termitidae gave rise to an enormous dietary diversification (chronogram adapted from Legendre et al. 2015b). Dashed lines for the common ancestor of cockroaches denote that the group is older than 201 Ma. The asterisk (*) specifies that some families or subfamilies are paraphyletic according to the latest studies. Figure 2. View largeDownload slide Eusociality in termites: life cycle (a) and b) evolutionary origin associated to putative key innovations. The life cycle differs for species with pseudergates and species with “true workers” (highlighted by the dashed box). The phylogeny of termites (Isoptera) illustrates that the origin of eusociality is nested within the radiation of cockroaches (termites and cockroaches together forming the Blattodea). The phylogeny also highlights important events in the evolution of the digestive symbiosis: the presence of cellulolytic flagellates, acquired by a common ancestor of termites and Cryptocercidae, is associated to their wood-feeding lifestyle; the loss of flagellates in Termitidae gave rise to an enormous dietary diversification (chronogram adapted from Legendre et al. 2015b). Dashed lines for the common ancestor of cockroaches denote that the group is older than 201 Ma. The asterisk (*) specifies that some families or subfamilies are paraphyletic according to the latest studies. The study of termite digestive systems might also bring important insights on the ecological and evolutionary successes of these insects (Fig. 2b). Although digestive strategies and the gut microbiota of termites remain to be fully deciphered, important progress has been made and two groups can be distinguished by their primary cellulolytic partners and their hindgut compartmentation (Brune 2014). The first group comprises termite species that harbor cellulolytic flagellates in their gut; they group all termite species except those belonging to Termitidae. In the second group (i.e., Termitidae), these flagellates are absent; instead Termitidae have an entirely prokaryotic microbiota. They also show a higher compartmentation of the hindgut (except for Macrotermitinae) and a large dietary diversification (Brune and Dietrich 2015). Because most termite species belong to the family Termitidae (Beccaloni and Eggleton 2013), their peculiar digestive system might have been a trigger of diversification. A parallel can be made with cockroach species; even though the gut microbiota of most cockroach species has not been investigated, flagellates have been found in some wood-eating species, including, but not restricted to (e.g., Pellens et al. 2002), Cryptocercus, the alleged sister group of termites (Cleveland et al. 1934; Klass et al. 2008). Because of characteristics shared with “lower” termites (Klass et al. 2008), it is then generally assumed that Cryptocercus belongs to the group with cellulolytic flagellates in their gut (Brune and Dietrich 2015), contrary to other supposedly or asserted wood-eating species (but see Pellens et al. 2007 and the section Trait-dependent diversification for more details and alternative codings). Our knowledge on the role of eusociality in the dynamic of diversification is limited. Only one previous study has investigated the diversification of termites (Davis et al. 2009), suggesting increases in diversification early in termite evolution. However, this study was based on a family-level supertree, and methods to estimate diversification were still in their infancy. Besides, this study did not differentiate the potential role of “true” workers or that of gut microbiota composition in termite diversification. To better understand the triggers of diversification in Dictyoptera (mantises, cockroaches, and termites), we analyze the most comprehensive time-calibrated phylogeny of the group. We investigate the diversification dynamic of termites compared to their dictyopteran counterparts using a battery of birth–death methods to corroborate and strengthen these results. More specifically, we tested the competing hypotheses in which eusociality is associated with lower (Wilson 1992) or higher (Davis et al. 2009) net diversification rate compared to non-eusociality. We also tested whether the presence of “true” workers in termite societies can be considered a key innovation relative to societies with pseudergates, by assessing the putative role on diversification of these castes on termite evolution. Similarly, we tested whether an entirely prokaryotic gut microbiota could also be a key innovation for termite diversification. Materials and Methods Analytical Pipeline We used a recently published time-calibrated phylogeny of Dictyoptera including 762 species representatives of living dictyopteran diversity (ca. 10,000 known species, Legendre et al. 2015b). This phylogeny resulted from a molecular data set comprising four mitochondrial and two nuclear markers. Divergence times were estimated using a relaxed-clock approach coupled with 17 fossils, which were used as minimum age constraints on different nodes and a maximum age for the tree root (Legendre et al. 2015b). Three putatively key innovations were investigated: eusociality, societies with “true” workers and entirely prokaryotic gut microbiota. To test the role of these innovations, we used three approaches: i) the maximum likelihood (ML) approach of time-dependent diversification (Morlon et al. 2011), implemented in the R-package RPANDA v.1.2 (Morlon et al. 2016), ii) the Bayesian analysis of macroevolutionary mixture (BAMM v.2.5, Rabosky et al. 2013), implemented in the R-package BAMMtools v.2.1.4 (Rabosky 2014), and iii) the ML approach of trait-dependent diversification (Maddison et al. 2007; FitzJohn et al. 2009), implemented in the R-package diversitree v.0.9-8 (FitzJohn 2012). Each method is designed to estimate speciation and extinction rates, and the three methods are used to cross-test hypotheses and corroborate results. Nonetheless, it is worth mentioning that each method differs at several points in the way speciation and extinction rates are estimated. For instance, trait-dependent birth–death models (diversitree) estimate constant speciation and extinction rates, while time-dependent birth–death models (BAMM and RPANDA) estimate the speciation and extinction rates and their variation through time. Therefore, we expect some differences in the values of the estimated diversification rates that are inherent to each approach. Our diversification analyzes should be seen as complementary for the inferred diversification trend rather than corroborative on the values and magnitude of the speciation and extinction rates. Across-Clade and Time-Variation Diversification The ML approach of Morlon et al. (2011) is a birth–death method that extends previous birth–death methods such that speciation and/or extinction rates may change continuously through time, and subclades may have different speciation and extinction rates. This method has the advantage of not assuming constant extinction rate over time (unlike BAMM, Rabosky et al. 2013), and allows clades to have declining diversity because extinction can exceed speciation, meaning that diversification rates can be negative (Morlon et al. 2011). We designed six nested diversification models to test with this approach: i) a Yule model, where speciation is constant and extinction is null; ii) a constant birth–death model, where speciation and extinction rates are constant; iii) a variable speciation rate model without extinction; iv) a variable speciation rate model with constant extinction; v) a rate-constant speciation and variable extinction rate model; and vi) a model in which both speciation and extinction rates vary. Models were compared by computing the ML score of each model and the resulting Akaike information criterion corrected by sample size (AICc). First, the Dictyoptera tree was analyzed as a whole using this approach and taking into account missing species ($$f = 0.08$$). To test the hypotheses that eusociality, societies with “true” workers and gut microbiota composition act on diversification rates, we compared the likelihoods of models that allow for different patterns of rate variation in different clades. More specifically, we first allowed for a rate shift along the branch corresponding to the apparition of eusociality (the termite crown node), and we sequentially tested for a rate shift between the crown node of termites and the node corresponding to the family Termitidae (prokaryotic microbiota). Therefore, we partitioned the Dictyoptera phylogeny into seven evolutionary scenarios: from a backbone without the termites, and the termites as a subtree; to a backbone without the Termitidae, and the Termitidae as a subtree (with all five possibilities in between; for a schematic illustration see Supplementary Fig. S1 available on Dryad at http://dx.doi.org/10.5061/dryad.2tg34). After creating these backbones and corresponding subtrees, and accounting for the missing species in each of them, we fitted the same diversification models as detailed above but independently on each subtree and its corresponding backbone tree (Supplementary Fig. S1 available on Dryad). The global log-likelihood was calculated as the sum of both subclade and backbone log-likelihoods as determined by the corresponding best-fitting models. This approach allowed us to compare the global log-likelihood (under the assumption of a single speciation rate and extinction rate) to the log-likelihood of one shift at the termite crown (eusociality), and to the log-likelihood of one shift at the crown of Termitidae (prokaryotic microbiota). We directly compare the seven scenarios using AICc. We also used BAMM to estimate speciation and extinction rates through time and among/within clades (Rabosky et al. 2013). BAMM was constructed to study complex evolutionary processes on phylogenetic trees, potentially shaped by a heterogeneous mixture of distinct evolutionary dynamics of speciation and extinction across clades. BAMM can automatically detect rate shifts and sample distinct evolutionary dynamics that explain the diversification dynamics of a clade without a priori hypotheses on how many and where these shifts might occur. Evolutionary dynamics can involve time-variable diversification rates; in BAMM, speciation is allowed to vary exponentially through time while extinction is maintained constant: subclades in a tree may diversify faster (or slower) than others. This Bayesian approach can be useful in detecting shifts of diversification potentially associated with key innovations (Rabosky 2014). We ran BAMM by setting four Markov chain Monte Carlo running (MCMC) for 20 million generations and sampled every 2000 generations. A compound Poisson process is implemented for the prior probability of a rate shift along any branch. We used a gradient of prior values ranging from 0.1 to 50 to test the sensitivity to the prior, because it has been shown that BAMM can be affected by the prior (Moore et al. 2016 but see Rabosky et al. 2017). We accounted for nonrandom incomplete taxon sampling using the implemented analytical correction; we set a sampling fraction for mantises ($$f = 0.118$$), cockroaches ($$f = 0.042$$), and termites ($$f = 0.094$$). We performed independent runs (with a 15% burn-in) using different seeds to assess the convergence of the runs with effective sample size. We processed the output data using the BAMMtools (Rabosky 2014) by estimating i) the mean global rates of diversification through time, ii) the configuration of the diversification rate shifts evaluating alternative diversification models as compared by posterior probabilities, and iii) the clade-specific rates through time when a distinct macroevolutionary regime is identified. Trait-Dependent Diversification We used trait-dependent diversification models to simultaneously model trait evolution and its impact on diversification (Maddison et al. 2007). In those models, species are characterized by an evolving trait, and their diversification follows a birth–death process in which speciation and extinction rates may depend on the trait state. We created three data sets of traits. First, we categorized all dictyopteran species as being non-eusocial or eusocial. This two-trait scheme distinguishes the termites (all species are eusocial) from the rest of dictyopterans (all non-eusocial). Second, we made a three-trait scheme by keeping the non-eusocial category, and dividing the eusocial trait into two: eusocial with pseudergates (“false” workers) and eusocial with “true” workers. Third, we made another three-trait scheme to test the role of the primary cellulolytic partners in diversification. We followed Brune and Dietrich (2015) to distinguish Termitidae (entirely prokaryotic microbiota), Cryptocercus$$+$$ other termites (cellulolytic flagellates), and the remaining Dictyoptera (no specialized microbiota for lignocellulose digestion). All these coding schemes are imperfect because, like any broad categorization process, they might oversimplify the reality found in nature (Legendre et al. 2015a; Goutte et al. 2016). The category “remaining Dictyoptera” in the latter coding scheme, for instance, is imperfect because it is delineated by default and because of a lack of information for some cockroach species. However, this three-state coding reflects our current understanding of the primary cellulolytic partners in Dictyoptera as reviewed in Brune and Dietrich (2015). Anyway, to take into account our lack of knowledge on cockroach microbiota and the existence of species that, arguably, could not fit well in these three categories – while limiting the number of categories to save statistical power, we created two supplementary trait data sets. In the first alternative coding, the wood-feeding non-eusocial lineages were coded with cellulolytic flagellates (same as “lower” termites) even though the existence and role of cellulolytic flagellates have not always been investigated and proved. These lineages represent several genera of the subfamily Panesthiinae, and the genera Cyrtotria, Lauraesilpha, Paramuzoa, Parasphaeria, and Colapteroblatta (Grandcolas 1993; Grandcolas et al. 2002; Pellens et al. 2002). In the second alternative coding, the Macrotermitinae were not coded with the Termitidae (i.e., initial coding supported by their entirely prokaryotic gut microbiota) but with the other termites (coding supported by their relatively simple hindgut structure). We applied the Binary State Speciation and Extinction model (BiSSE, Maddison et al. 2007) with constant diversification rate for the two-trait data set. We then performed the Multi-State Speciation and Extinction model (MuSSE, FitzJohn et al. 2009) on the three-trait data sets. Both BiSSE and MuSSE models account for incomplete taxon sampling, which is informed as a sampling fraction of species at present having a given trait (FitzJohn et al. 2009). Simulation studies have shown that a large, well-sampled tree is required by these methods, whereas trees containing fewer than 300 species may lack sufficient phylogenetic signal to produce enough statistical power (Davis et al. 2013). The BiSSE model has six distinct parameters: two speciation rates without character change (i.e., in situ speciation) associated with non-eusocial species ($$\lambda_{\mathrm{N}})$$ and eusocial species ($$\lambda_{\rm E})$$, two extinction rates associated with non-eusocial ($$\mu_{\mathrm{N}})$$ and eusocial species ($$\mu_{\rm E})$$, and two transition rates (i.e., anagenetic change) with one from non-eusocial to eusocial ($$q_{\mathrm{N-E}})$$, and from eusocial to non-eusocial ($$q_{\mathrm{E-N}})$$. Using the same rationale, but applied to three traits, the MuSSE model has 12 distinct parameters (three for speciation, three for extinction, and six for transitions). Here we used 100 dictyopteran phylogenies (sampled from the Bayesian analysis in Legendre et al. 2015b) and combined them with our three data sets of traits. Analyzes were performed using the R-package diversitree (FitzJohn 2012) using the functions make.bisse and make.musse to construct the likelihood functions for each model based on the data, and the functions constrain and find.mle to apply different diversification scenarios. For each model, we computed the AICc based on the log-likelihood and the number of parameters. We checked support for the selected model against all models using the difference between AICc ($$\Delta$$AIC) and the Akaike weight (AIC$$\omega$$). The scenario supported with the lowest AICc was considered the best when $$\Delta$$AIC $$>2$$ against other models and with AIC$$\omega > 0.5$$ (otherwise the model with less parameter was instead considered the best). Finally, we used the consensus tree and MCMC simulations with the best model to examine the confidence interval of the parameter estimates. We used an exponential prior and started the chain with the parameters obtained by ML (FitzJohn 2012). We ran 20,000 MCMC steps and applied a 10% burn-in. Net diversification rates were then computed. SSE models have recently been criticized due to type I error (Rabosky and Goldberg 2015). To test whether our diversification results were biased, we estimated the difference of fit ($$\Delta$$AIC) between the best model and a null model (i.e., with no state dependence) and compared this with the difference between the same models as estimated from 100 simulated trait data sets for each SSE model. Results Diversification Across-Clade and Time The BAMM analyzes supported a model with six evolutionary regimes (i.e., five rate shifts; Supplementary Fig. S2 and Table S1 available on Dryad). The post burn-in posterior distribution indicated that five shifts occur with a posterior probability (PP) of 0.54, and PP $$=$$ 0.33 and 0.10 for four and six shifts, respectively (less or more shifts occur with PP $$<$$ 0.05; Supplementary Table S1 available on Dryad). The 95% credible set of shift configurations and marginal shift probabilities strongly favored a model that includes shifts along branches within termites (Supplementary Figs. S3 and S4 available on Dryad). The shift configurations within this credible set contained also shifts at internal nodes in mantises and at a more derived position within cockroaches (Blaberidae: Panesthiinae). The different shift configurations in the credible set allowed quantifying uncertainty in placement of a termite shift (Supplementary Fig. S4 available on Dryad). The shift configurations sampled at the highest frequency always contained a shift within termites (cumulative posterior probability of 0.62). The best configuration shift retained five core shifts, one located within the termites (crown Termitidae), two within Mantodea and two others within cockroaches (Supplementary Fig. S5 available on Dryad). We repeated the analyzes by changing the Poisson process for the prior probability of a rate shift along any branch: these analyzes supported a similar diversification pattern, although we found an additional clade-specific macroevolutionary regime when a higher prior probability of rate shift was used (Supplementary Fig. S6 and Table S1 available on Dryad). We found evidence for low and constant-rate diversification through time for Dictyoptera punctuated by few shifts. When diversification changed, speciation rates increased for the termites in particular at the base of Termitidae (Fig. 3). Although other increases of speciation were found, the most elevated speciation rate (and also net diversification rate) occurred within termites. Interestingly, when termite macroevolutionary regime was studied in isolation, the net diversification rate of the clade increased towards the present (Fig. 3). On the contrary, the other macroevolutionary regimes showed a more common pattern of diversification slowdown over time. Figure 3. View largeDownload slide Diversification and eusociality in Dictyoptera. a) Time-calibrated phylogeny of Dictyoptera (adapted from Legendre et al. 2015b) showing the origin of eusociality (all termites), the origin of an entirely prokaryotic microbiota (Termitidae), and the shifts in diversification rate identified with BAMM. b) Results of BiSSE show that lineages with eusocial species have higher net diversification rates than non-eusocial species. c) Results of MuSSE show that lineages with “true” workers have higher net diversification rates than other dictyopteran lineages, including termite lineages with pseudergates. d) Results of MuSSE show that lineages with an entirely prokaryotic microbiota (and a complex hindgut compartmentation) have a higher diversification rate than other dictyopteran lineages. For BiSSE and MuSSE plots, Bayesian posterior distributions represent the 95% credibility interval of each estimated parameter. e) BAMM analyzes suggest a strong increase in net diversification rate of termites towards the present, especially in termites with “true” workers (mostly Termitidae and Rhinotermitidae). The shaded areas represent the 95% credibility interval of each estimated parameter. Figure 3. View largeDownload slide Diversification and eusociality in Dictyoptera. a) Time-calibrated phylogeny of Dictyoptera (adapted from Legendre et al. 2015b) showing the origin of eusociality (all termites), the origin of an entirely prokaryotic microbiota (Termitidae), and the shifts in diversification rate identified with BAMM. b) Results of BiSSE show that lineages with eusocial species have higher net diversification rates than non-eusocial species. c) Results of MuSSE show that lineages with “true” workers have higher net diversification rates than other dictyopteran lineages, including termite lineages with pseudergates. d) Results of MuSSE show that lineages with an entirely prokaryotic microbiota (and a complex hindgut compartmentation) have a higher diversification rate than other dictyopteran lineages. For BiSSE and MuSSE plots, Bayesian posterior distributions represent the 95% credibility interval of each estimated parameter. e) BAMM analyzes suggest a strong increase in net diversification rate of termites towards the present, especially in termites with “true” workers (mostly Termitidae and Rhinotermitidae). The shaded areas represent the 95% credibility interval of each estimated parameter. Termite Key Innovations and Diversification The RPANDA analyzes refined the location of the rate shift within termites, and allowed testing the effect of three termite key innovations: eusociality, societies with “true” workers, and entirely prokaryotic gut microbiota, which have evolved sequentially from the origin of termites to the origin of Termitidae. Applying a series of time-dependent birth–death models (including a rate shift at the crown node of a termite subclade), we compared the fit of alternative scenarios to locate the best change of diversification rates. We first found that any scenario including a rate shift within termites is better than a scenario having no shift at all (Table 1). Comparing all diversification scenarios including a shift in the termite tree, we found that the best-fitting scenario includes a shift along a branch supporting termites that only contain species with “true” workers in their societies (Scenario 6, Table 1). There is strong support for this scenario (lowest $$\Delta$$AIC $$=$$ 12.7) as compared to competing hypotheses of other locations in the termite tree such as the evolution of eusociality (Scenario 1), or the change in gut microbiota composition (Scenario 7). Table 1. Support for an increase of diversification rates for termites, especially for lineages including “true” workers Scenario with a rate shift located at: Models Whole phylogeny The MRCA of Isoptera The MRCA of Euisoptera The MRCA of Kalotermitidae and Neoisoptera The MRCA of Neoisoptera The MRCA of Neoisoptera minus Rhinotermitinae The MRCA of Termitidae, Coptotermitinae and Heterotermitinae The MRCA of Termitidae NP 3 4 4 5 4 5 5 4 logL –3889.247 –3709.158 –3706.857 –3693.107 –3658.84 –3659.554 –3651.496 –3691.284 AICc 7784.526 7426.385 7421.781 7396.308 7325.754 7329.206 7313.092 7390.653 $$\Delta$$AIC 471.434 113.293 108.689 83.216 12.662 16.114 0 77.561 AIC$$\omega$$ 0 0 0 0 $$\approx$$ 0.002 $$<$$ 0.0001 0.998 0 Subclade speciation rate 0.191 0.177 0.178 0.185 0.186 0.191 0.193 0.197 Subclade extinction rate 0.175 0 0 0 0 0 0 0 Subclade net div. rate 0.016 0.177 0.178 0.185 0.186 0.191 0.193 0.197 Backbone speciation rate — 0.044 0.044 0.057 0.044 0.058 0.058 0.084 Backbone extinction rate — 0 0 0.024 0 0.024 0.023 0.06 Backbone net. div. rate — 0.044 0.044 0.033 0.044 0.034 0.035 0.024 Scenario with a rate shift located at: Models Whole phylogeny The MRCA of Isoptera The MRCA of Euisoptera The MRCA of Kalotermitidae and Neoisoptera The MRCA of Neoisoptera The MRCA of Neoisoptera minus Rhinotermitinae The MRCA of Termitidae, Coptotermitinae and Heterotermitinae The MRCA of Termitidae NP 3 4 4 5 4 5 5 4 logL –3889.247 –3709.158 –3706.857 –3693.107 –3658.84 –3659.554 –3651.496 –3691.284 AICc 7784.526 7426.385 7421.781 7396.308 7325.754 7329.206 7313.092 7390.653 $$\Delta$$AIC 471.434 113.293 108.689 83.216 12.662 16.114 0 77.561 AIC$$\omega$$ 0 0 0 0 $$\approx$$ 0.002 $$<$$ 0.0001 0.998 0 Subclade speciation rate 0.191 0.177 0.178 0.185 0.186 0.191 0.193 0.197 Subclade extinction rate 0.175 0 0 0 0 0 0 0 Subclade net div. rate 0.016 0.177 0.178 0.185 0.186 0.191 0.193 0.197 Backbone speciation rate — 0.044 0.044 0.057 0.044 0.058 0.058 0.084 Backbone extinction rate — 0 0 0.024 0 0.024 0.023 0.06 Backbone net. div. rate — 0.044 0.044 0.033 0.044 0.034 0.035 0.024 Notes: The table reports the log-likelihood (logL) and the corrected Akaike Information Criterion (AICc) calculated for the whole phylogeny of Dictyoptera (assuming no shift of diversification) and seven evolutionary scenarios including a shift of diversification (as specified in the table). For each model, the parameter values of the best-fitting model are reported (rates are in events/Myr). The “backbone” rates represent the estimated diversification parameters for the tree containing all Dictyoptera species except the termite subclades. The termite subclades are analyzed aside. The logL, the number of parameters (NP), and AICc are thus the sum of the logL/parameters/AICc of the backbone tree plus the logL/parameters/AICc of the termite subclade. The best-fitting scenario is determined by $$\Delta $$AIC and AIC$$\omega $$ and includes a shift within the termites that only include species with “true” workers in their societies (Scenario 7, highlighted in bold). Table 1. Support for an increase of diversification rates for termites, especially for lineages including “true” workers Scenario with a rate shift located at: Models Whole phylogeny The MRCA of Isoptera The MRCA of Euisoptera The MRCA of Kalotermitidae and Neoisoptera The MRCA of Neoisoptera The MRCA of Neoisoptera minus Rhinotermitinae The MRCA of Termitidae, Coptotermitinae and Heterotermitinae The MRCA of Termitidae NP 3 4 4 5 4 5 5 4 logL –3889.247 –3709.158 –3706.857 –3693.107 –3658.84 –3659.554 –3651.496 –3691.284 AICc 7784.526 7426.385 7421.781 7396.308 7325.754 7329.206 7313.092 7390.653 $$\Delta$$AIC 471.434 113.293 108.689 83.216 12.662 16.114 0 77.561 AIC$$\omega$$ 0 0 0 0 $$\approx$$ 0.002 $$<$$ 0.0001 0.998 0 Subclade speciation rate 0.191 0.177 0.178 0.185 0.186 0.191 0.193 0.197 Subclade extinction rate 0.175 0 0 0 0 0 0 0 Subclade net div. rate 0.016 0.177 0.178 0.185 0.186 0.191 0.193 0.197 Backbone speciation rate — 0.044 0.044 0.057 0.044 0.058 0.058 0.084 Backbone extinction rate — 0 0 0.024 0 0.024 0.023 0.06 Backbone net. div. rate — 0.044 0.044 0.033 0.044 0.034 0.035 0.024 Scenario with a rate shift located at: Models Whole phylogeny The MRCA of Isoptera The MRCA of Euisoptera The MRCA of Kalotermitidae and Neoisoptera The MRCA of Neoisoptera The MRCA of Neoisoptera minus Rhinotermitinae The MRCA of Termitidae, Coptotermitinae and Heterotermitinae The MRCA of Termitidae NP 3 4 4 5 4 5 5 4 logL –3889.247 –3709.158 –3706.857 –3693.107 –3658.84 –3659.554 –3651.496 –3691.284 AICc 7784.526 7426.385 7421.781 7396.308 7325.754 7329.206 7313.092 7390.653 $$\Delta$$AIC 471.434 113.293 108.689 83.216 12.662 16.114 0 77.561 AIC$$\omega$$ 0 0 0 0 $$\approx$$ 0.002 $$<$$ 0.0001 0.998 0 Subclade speciation rate 0.191 0.177 0.178 0.185 0.186 0.191 0.193 0.197 Subclade extinction rate 0.175 0 0 0 0 0 0 0 Subclade net div. rate 0.016 0.177 0.178 0.185 0.186 0.191 0.193 0.197 Backbone speciation rate — 0.044 0.044 0.057 0.044 0.058 0.058 0.084 Backbone extinction rate — 0 0 0.024 0 0.024 0.023 0.06 Backbone net. div. rate — 0.044 0.044 0.033 0.044 0.034 0.035 0.024 Notes: The table reports the log-likelihood (logL) and the corrected Akaike Information Criterion (AICc) calculated for the whole phylogeny of Dictyoptera (assuming no shift of diversification) and seven evolutionary scenarios including a shift of diversification (as specified in the table). For each model, the parameter values of the best-fitting model are reported (rates are in events/Myr). The “backbone” rates represent the estimated diversification parameters for the tree containing all Dictyoptera species except the termite subclades. The termite subclades are analyzed aside. The logL, the number of parameters (NP), and AICc are thus the sum of the logL/parameters/AICc of the backbone tree plus the logL/parameters/AICc of the termite subclade. The best-fitting scenario is determined by $$\Delta $$AIC and AIC$$\omega $$ and includes a shift within the termites that only include species with “true” workers in their societies (Scenario 7, highlighted in bold). BAMM identified an important rate shift within the mantises (166.8 Ma, Fig. 3), which corresponds to an increase of diversification. Using RPANDA, we tested whether a model including this shift performs better than the model with the best-fitting shift in termites. We reproduced the approach explained above but we isolated the mantis clade as a subtree, from the rest of the Dictyoptera (this time including termites). We analyzed all the same models and found that the shift within the Mantodea does not improve the likelihood as compared with the model including the shift within termites (AICc for the shift in mantises $$=$$ 7640.83; resulting in a $$\Delta$$AIC $$=$$ 327 between the best scenario including a shift within termites and this scenario, Table 1). RPANDA analyzes endorsed BAMM results: The best-fit scenario indicated a speciation rate increasing through time for both the termite tree and the backbone. Nonetheless, we estimated an elevated speciation rate of 0.193 lineage/Myr for termites, strikingly contrasting with the speciation rate of the remaining dictyopteran lineages (0.058 lineage/Myr). Moreover, we found the termite clade diversified with low or zero extinction rates, again contrasting with the extinction rate for the rest of the tree (0.023 lineage/Myr). Consequently, the net diversification rate of the termite tree is 5.5 times higher than the remaining part of the dictyopteran tree (Scenario 7, Table 1). We corroborated the BAMM and RPANDA analyzes with BiSSE and MuSSE applied to a two-trait data set and two three-trait data sets, respectively. We found that the best-fitting BiSSE model was the one with different speciation and extinction, but with equal transition rates. Simpler or more complex models were not supported ($$\Delta$$AIC $$=$$ 4.13 with the second best-fit BiSSE model, which has one extra parameter, Table 2, Supplementary Table S2 available on Dryad). In the best-fit model, both speciation and extinction rates were higher for eusocial species (Table 2), and the net diversification rate was 2-fold higher for eusocial lineages (Fig. 3). Bayesian MCMC analysis showed that speciation, extinction, and net diversification rates were significantly different between the two traits (Supplementary Fig. S7 available on Dryad). Table 2. Supports for eusociality as a driver of diversification Models NP logL AICc $$\Delta $$AICc AIC$$\omega $$ $$\lambda 0$$ $$\lambda 1$$ $$\mu$$0 $$\mu$$1 $$q$$01 $$q$$10 Null model 3 $$-3831$$.072 7668.175 218.2023 0 0.1980 0.1786 2.01E$$-$$05 $$\lambda $$0 $$\ne \lambda 1$$, $$\mu 0 = \mu 1$$, $$q 01 = q10$$ 4 $$-3752$$.206 7512.465 62.492 0 0.1335 0.1867 0.1129 1.64E$$-$$05 $$\lambda 0 = \lambda 1 $$, $$\mu 0 \ne \mu 1$$, $$ q01 = q 10$$ 4 $$-3773$$.267 7554.586 104.6133 0 0.1626 0.1441 0.0921 1.24E$$-$$05 $$\lambda 0 = \lambda 1 $$, $$\mu 0 = \mu 1$$, $$q 01 \ne 10$$ 4 $$-383$$0.967 7669.988 220.0149 0 0.1979 0.1784 2.43E$$-$$05 4.97E$$-$$06 $$\lambda \textbf{0} \ne \lambda 1$$, $$\mu \textbf{0} \ne \mu 1$$, $$q01 = q \textbf{10}$$ 5 $$-3719$$.947 7449.973 0 0.887 0.1004 0.4145 0.0766 0.3721 0.3721 $$\lambda 0 \ne \lambda 1$$, $$\mu 0 = \mu 1$$, $$ q 01 \ne q\ne $$10 5 $$-3752$$.346 7514.772 64.7993 0 0.1280 0.1815 0.1071 1.81E$$-$$05 1.20E$$-$$05 $$\lambda $$0 $$= \lambda 1$$, $$\mu $$0 $$\ne \mu 1$$,$$ q$$01 $$\ne q\ne $$10 5 $$-3774$$.607 7559.294 109.3214 0 0.1541 0.1349 0.0810 1.27E$$-$$05 7.40E$$-$$06 $$\lambda $$0 $$\ne \lambda 1$$, $$\mu $$0 $$\ne \mu 1$$,$$ q$$01 $$\ne q\ne $$10 6 $$-372$$0.997 7454.106 4.1336 0.113 0.0910 0.4099 0.0667 0.3674 4.35E$$-$$05 5.91E$$-$$06 Models NP logL AICc $$\Delta $$AICc AIC$$\omega $$ $$\lambda 0$$ $$\lambda 1$$ $$\mu$$0 $$\mu$$1 $$q$$01 $$q$$10 Null model 3 $$-3831$$.072 7668.175 218.2023 0 0.1980 0.1786 2.01E$$-$$05 $$\lambda $$0 $$\ne \lambda 1$$, $$\mu 0 = \mu 1$$, $$q 01 = q10$$ 4 $$-3752$$.206 7512.465 62.492 0 0.1335 0.1867 0.1129 1.64E$$-$$05 $$\lambda 0 = \lambda 1 $$, $$\mu 0 \ne \mu 1$$, $$ q01 = q 10$$ 4 $$-3773$$.267 7554.586 104.6133 0 0.1626 0.1441 0.0921 1.24E$$-$$05 $$\lambda 0 = \lambda 1 $$, $$\mu 0 = \mu 1$$, $$q 01 \ne 10$$ 4 $$-383$$0.967 7669.988 220.0149 0 0.1979 0.1784 2.43E$$-$$05 4.97E$$-$$06 $$\lambda \textbf{0} \ne \lambda 1$$, $$\mu \textbf{0} \ne \mu 1$$, $$q01 = q \textbf{10}$$ 5 $$-3719$$.947 7449.973 0 0.887 0.1004 0.4145 0.0766 0.3721 0.3721 $$\lambda 0 \ne \lambda 1$$, $$\mu 0 = \mu 1$$, $$ q 01 \ne q\ne $$10 5 $$-3752$$.346 7514.772 64.7993 0 0.1280 0.1815 0.1071 1.81E$$-$$05 1.20E$$-$$05 $$\lambda $$0 $$= \lambda 1$$, $$\mu $$0 $$\ne \mu 1$$,$$ q$$01 $$\ne q\ne $$10 5 $$-3774$$.607 7559.294 109.3214 0 0.1541 0.1349 0.0810 1.27E$$-$$05 7.40E$$-$$06 $$\lambda $$0 $$\ne \lambda 1$$, $$\mu $$0 $$\ne \mu 1$$,$$ q$$01 $$\ne q\ne $$10 6 $$-372$$0.997 7454.106 4.1336 0.113 0.0910 0.4099 0.0667 0.3674 4.35E$$-$$05 5.91E$$-$$06 Notes: The table reports the models applied, their number of parameters (NP), the means for the log-likelihood (logL) and the corrected Akaike Information Criterion (AICc), the difference of AICc between the best model (lowest AIC) and a given model ($$\Delta $$AICc), the Akaike weight (AIC$$\omega )$$, and the values for each parameter of the corresponding model based on 100 dated trees. The best model is highlighted in bold. Standard errors for each parameter estimate are reported in Supplementary Table S2 available on Dryad. 0$$=$$ non-eusocial and 1$$=$$ eusocial; $$\lambda =$$ speciation rate (one for each character state, $$\lambda 1$$ is the speciation rate for eusocial species); $$\mu$$$$=$$ extinction rate (one for each character state, $$\mu 1$$ is the extinction rate for eusocial species); and $$q =$$ transition rate between character states. Table 2. Supports for eusociality as a driver of diversification Models NP logL AICc $$\Delta $$AICc AIC$$\omega $$ $$\lambda 0$$ $$\lambda 1$$ $$\mu$$0 $$\mu$$1 $$q$$01 $$q$$10 Null model 3 $$-3831$$.072 7668.175 218.2023 0 0.1980 0.1786 2.01E$$-$$05 $$\lambda $$0 $$\ne \lambda 1$$, $$\mu 0 = \mu 1$$, $$q 01 = q10$$ 4 $$-3752$$.206 7512.465 62.492 0 0.1335 0.1867 0.1129 1.64E$$-$$05 $$\lambda 0 = \lambda 1 $$, $$\mu 0 \ne \mu 1$$, $$ q01 = q 10$$ 4 $$-3773$$.267 7554.586 104.6133 0 0.1626 0.1441 0.0921 1.24E$$-$$05 $$\lambda 0 = \lambda 1 $$, $$\mu 0 = \mu 1$$, $$q 01 \ne 10$$ 4 $$-383$$0.967 7669.988 220.0149 0 0.1979 0.1784 2.43E$$-$$05 4.97E$$-$$06 $$\lambda \textbf{0} \ne \lambda 1$$, $$\mu \textbf{0} \ne \mu 1$$, $$q01 = q \textbf{10}$$ 5 $$-3719$$.947 7449.973 0 0.887 0.1004 0.4145 0.0766 0.3721 0.3721 $$\lambda 0 \ne \lambda 1$$, $$\mu 0 = \mu 1$$, $$ q 01 \ne q\ne $$10 5 $$-3752$$.346 7514.772 64.7993 0 0.1280 0.1815 0.1071 1.81E$$-$$05 1.20E$$-$$05 $$\lambda $$0 $$= \lambda 1$$, $$\mu $$0 $$\ne \mu 1$$,$$ q$$01 $$\ne q\ne $$10 5 $$-3774$$.607 7559.294 109.3214 0 0.1541 0.1349 0.0810 1.27E$$-$$05 7.40E$$-$$06 $$\lambda $$0 $$\ne \lambda 1$$, $$\mu $$0 $$\ne \mu 1$$,$$ q$$01 $$\ne q\ne $$10 6 $$-372$$0.997 7454.106 4.1336 0.113 0.0910 0.4099 0.0667 0.3674 4.35E$$-$$05 5.91E$$-$$06 Models NP logL AICc $$\Delta $$AICc AIC$$\omega $$ $$\lambda 0$$ $$\lambda 1$$ $$\mu$$0 $$\mu$$1 $$q$$01 $$q$$10 Null model 3 $$-3831$$.072 7668.175 218.2023 0 0.1980 0.1786 2.01E$$-$$05 $$\lambda $$0 $$\ne \lambda 1$$, $$\mu 0 = \mu 1$$, $$q 01 = q10$$ 4 $$-3752$$.206 7512.465 62.492 0 0.1335 0.1867 0.1129 1.64E$$-$$05 $$\lambda 0 = \lambda 1 $$, $$\mu 0 \ne \mu 1$$, $$ q01 = q 10$$ 4 $$-3773$$.267 7554.586 104.6133 0 0.1626 0.1441 0.0921 1.24E$$-$$05 $$\lambda 0 = \lambda 1 $$, $$\mu 0 = \mu 1$$, $$q 01 \ne 10$$ 4 $$-383$$0.967 7669.988 220.0149 0 0.1979 0.1784 2.43E$$-$$05 4.97E$$-$$06 $$\lambda \textbf{0} \ne \lambda 1$$, $$\mu \textbf{0} \ne \mu 1$$, $$q01 = q \textbf{10}$$ 5 $$-3719$$.947 7449.973 0 0.887 0.1004 0.4145 0.0766 0.3721 0.3721 $$\lambda 0 \ne \lambda 1$$, $$\mu 0 = \mu 1$$, $$ q 01 \ne q\ne $$10 5 $$-3752$$.346 7514.772 64.7993 0 0.1280 0.1815 0.1071 1.81E$$-$$05 1.20E$$-$$05 $$\lambda $$0 $$= \lambda 1$$, $$\mu $$0 $$\ne \mu 1$$,$$ q$$01 $$\ne q\ne $$10 5 $$-3774$$.607 7559.294 109.3214 0 0.1541 0.1349 0.0810 1.27E$$-$$05 7.40E$$-$$06 $$\lambda $$0 $$\ne \lambda 1$$, $$\mu $$0 $$\ne \mu 1$$,$$ q$$01 $$\ne q\ne $$10 6 $$-372$$0.997 7454.106 4.1336 0.113 0.0910 0.4099 0.0667 0.3674 4.35E$$-$$05 5.91E$$-$$06 Notes: The table reports the models applied, their number of parameters (NP), the means for the log-likelihood (logL) and the corrected Akaike Information Criterion (AICc), the difference of AICc between the best model (lowest AIC) and a given model ($$\Delta $$AICc), the Akaike weight (AIC$$\omega )$$, and the values for each parameter of the corresponding model based on 100 dated trees. The best model is highlighted in bold. Standard errors for each parameter estimate are reported in Supplementary Table S2 available on Dryad. 0$$=$$ non-eusocial and 1$$=$$ eusocial; $$\lambda =$$ speciation rate (one for each character state, $$\lambda 1$$ is the speciation rate for eusocial species); $$\mu$$$$=$$ extinction rate (one for each character state, $$\mu 1$$ is the extinction rate for eusocial species); and $$q =$$ transition rate between character states. We further investigated whether difference of diversification rates occurred within the termites by splitting the eusociality trait into two (societies with pseudergates or societies with “true” workers). We found that the best-fitting MuSSE model supported the most complex model, in which all parameters are free. Simpler models were not supported ($$\Delta $$AICc $$=$$ 15.8 with the second best-fit MuSSE model, which has one less parameter, Table 3, Supplementary Table S2 available on Dryad). In the best-fit model, speciation and extinction rates were higher for eusocial lineages with “true” workers. Speciation and extinction rates were both low for the two other traits (non-eusocial, and eusocial with pseudergates). The net diversification rate was higher for eusocial lineages with “true” workers (Fig. 3). The trait-dependent rates of speciation, extinction, and net diversification for the eusocial species with “true” workers were significantly different from the rates of the two other traits, and the two other traits were not significantly different from each other (Supplementary Fig. S8 available on Dryad). Table 3. Supports for eusociality with true workers as drivers of diversification Models NP logL AICc $$\Delta $$AIC AI$$\omega $$ $$\lambda 1$$ $$\lambda 2$$ $$\lambda 3$$ $$\mu$$1 $$\mu$$2 $$\mu$$3 $$q$$12 $$q$$13 $$q$$21 $$q$$23 $$q$$31 $$q$$32 Null model 3 $$-3889$$.534 7785.100 373.94 0 0.1719 0.1525 4.21E$$-$$05 All $$\lambda$$ are varying 5 $$-3717$$.017 7444.114 32.95 0 0.0683 0.0727 0.1699 0.0438 2.43E$$-$$05 All $$\mu$$ are varying 5 $$-3751$$.789 7513.657 102.49 0 0.1214 0.1020 0.0974 0.0000 1.43E$$-$$05 All $$q$$ are varying 8 $$-3879$$.999 7776.191 365.03 0 0.1713 0.1519 1.29E$$-$$05 6.67E$$-$$07 1.67E$$-$$06 0.00161 4.64E$$-$$06 0.00023 All $$\lambda$$ and $$\mu$$ are varying 7 $$-37$$06.089 7426.327 15.16 $$\approx $$0.005 0.0637 0.0580 0.3198 0.0384 0.0266 0.2267 3.85E$$-$$05 All $$\lambda$$ and $$ q$$ are varying 10 $$-37$$03.355 7427.005 15.84 $$\approx $$0.005 0.0645 0.0693 0.1736 0.0394 1.55E$$-$$05 1.33E$$-$$07 7.90E$$-$$07 0.00155 9.81E$$-$$06 2.73E$$-$$06 All $$\mu $$ and $$q$$ are varying 10 $$-3738$$.133 7496.560 85.40 0 0.1258 0.1067 0.1030 0.0000 1.87E$$-$$05 7.14E$$-$$08 2.63E$$-$$07 0.00043 2.73E$$-$$06 0.00044 All rates are free 12 $${\bf -3693}$$.372 7411.163 0 0.99 0.0649 0.0378 0.3196 0.0398 0.0000 0.2244 1.09E$$-$$05 2.43E$$-$$07 6.45E$$-$$07 0.00224 1.32E$$-$$06 0.00027 Models NP logL AICc $$\Delta $$AIC AI$$\omega $$ $$\lambda 1$$ $$\lambda 2$$ $$\lambda 3$$ $$\mu$$1 $$\mu$$2 $$\mu$$3 $$q$$12 $$q$$13 $$q$$21 $$q$$23 $$q$$31 $$q$$32 Null model 3 $$-3889$$.534 7785.100 373.94 0 0.1719 0.1525 4.21E$$-$$05 All $$\lambda$$ are varying 5 $$-3717$$.017 7444.114 32.95 0 0.0683 0.0727 0.1699 0.0438 2.43E$$-$$05 All $$\mu$$ are varying 5 $$-3751$$.789 7513.657 102.49 0 0.1214 0.1020 0.0974 0.0000 1.43E$$-$$05 All $$q$$ are varying 8 $$-3879$$.999 7776.191 365.03 0 0.1713 0.1519 1.29E$$-$$05 6.67E$$-$$07 1.67E$$-$$06 0.00161 4.64E$$-$$06 0.00023 All $$\lambda$$ and $$\mu$$ are varying 7 $$-37$$06.089 7426.327 15.16 $$\approx $$0.005 0.0637 0.0580 0.3198 0.0384 0.0266 0.2267 3.85E$$-$$05 All $$\lambda$$ and $$ q$$ are varying 10 $$-37$$03.355 7427.005 15.84 $$\approx $$0.005 0.0645 0.0693 0.1736 0.0394 1.55E$$-$$05 1.33E$$-$$07 7.90E$$-$$07 0.00155 9.81E$$-$$06 2.73E$$-$$06 All $$\mu $$ and $$q$$ are varying 10 $$-3738$$.133 7496.560 85.40 0 0.1258 0.1067 0.1030 0.0000 1.87E$$-$$05 7.14E$$-$$08 2.63E$$-$$07 0.00043 2.73E$$-$$06 0.00044 All rates are free 12 $${\bf -3693}$$.372 7411.163 0 0.99 0.0649 0.0378 0.3196 0.0398 0.0000 0.2244 1.09E$$-$$05 2.43E$$-$$07 6.45E$$-$$07 0.00224 1.32E$$-$$06 0.00027 Notes:1$$=$$ non-eusocial or societies without worker, 2$$=$$ societies with pseudergates, and 3$$=$$ societies with true workers; $$\lambda =$$ speciation rate (one for each character state, $$\lambda $$1 is the speciation rate for species without worker); $$\mu$$$$=$$ extinction rate (one for each character state, $$\mu 1$$ is the extinction rate for species without worker); and $$q =$$ transition rate between character states. Table 3. Supports for eusociality with true workers as drivers of diversification Models NP logL AICc $$\Delta $$AIC AI$$\omega $$ $$\lambda 1$$ $$\lambda 2$$ $$\lambda 3$$ $$\mu$$1 $$\mu$$2 $$\mu$$3 $$q$$12 $$q$$13 $$q$$21 $$q$$23 $$q$$31 $$q$$32 Null model 3 $$-3889$$.534 7785.100 373.94 0 0.1719 0.1525 4.21E$$-$$05 All $$\lambda$$ are varying 5 $$-3717$$.017 7444.114 32.95 0 0.0683 0.0727 0.1699 0.0438 2.43E$$-$$05 All $$\mu$$ are varying 5 $$-3751$$.789 7513.657 102.49 0 0.1214 0.1020 0.0974 0.0000 1.43E$$-$$05 All $$q$$ are varying 8 $$-3879$$.999 7776.191 365.03 0 0.1713 0.1519 1.29E$$-$$05 6.67E$$-$$07 1.67E$$-$$06 0.00161 4.64E$$-$$06 0.00023 All $$\lambda$$ and $$\mu$$ are varying 7 $$-37$$06.089 7426.327 15.16 $$\approx $$0.005 0.0637 0.0580 0.3198 0.0384 0.0266 0.2267 3.85E$$-$$05 All $$\lambda$$ and $$ q$$ are varying 10 $$-37$$03.355 7427.005 15.84 $$\approx $$0.005 0.0645 0.0693 0.1736 0.0394 1.55E$$-$$05 1.33E$$-$$07 7.90E$$-$$07 0.00155 9.81E$$-$$06 2.73E$$-$$06 All $$\mu $$ and $$q$$ are varying 10 $$-3738$$.133 7496.560 85.40 0 0.1258 0.1067 0.1030 0.0000 1.87E$$-$$05 7.14E$$-$$08 2.63E$$-$$07 0.00043 2.73E$$-$$06 0.00044 All rates are free 12 $${\bf -3693}$$.372 7411.163 0 0.99 0.0649 0.0378 0.3196 0.0398 0.0000 0.2244 1.09E$$-$$05 2.43E$$-$$07 6.45E$$-$$07 0.00224 1.32E$$-$$06 0.00027 Models NP logL AICc $$\Delta $$AIC AI$$\omega $$ $$\lambda 1$$ $$\lambda 2$$ $$\lambda 3$$ $$\mu$$1 $$\mu$$2 $$\mu$$3 $$q$$12 $$q$$13 $$q$$21 $$q$$23 $$q$$31 $$q$$32 Null model 3 $$-3889$$.534 7785.100 373.94 0 0.1719 0.1525 4.21E$$-$$05 All $$\lambda$$ are varying 5 $$-3717$$.017 7444.114 32.95 0 0.0683 0.0727 0.1699 0.0438 2.43E$$-$$05 All $$\mu$$ are varying 5 $$-3751$$.789 7513.657 102.49 0 0.1214 0.1020 0.0974 0.0000 1.43E$$-$$05 All $$q$$ are varying 8 $$-3879$$.999 7776.191 365.03 0 0.1713 0.1519 1.29E$$-$$05 6.67E$$-$$07 1.67E$$-$$06 0.00161 4.64E$$-$$06 0.00023 All $$\lambda$$ and $$\mu$$ are varying 7 $$-37$$06.089 7426.327 15.16 $$\approx $$0.005 0.0637 0.0580 0.3198 0.0384 0.0266 0.2267 3.85E$$-$$05 All $$\lambda$$ and $$ q$$ are varying 10 $$-37$$03.355 7427.005 15.84 $$\approx $$0.005 0.0645 0.0693 0.1736 0.0394 1.55E$$-$$05 1.33E$$-$$07 7.90E$$-$$07 0.00155 9.81E$$-$$06 2.73E$$-$$06 All $$\mu $$ and $$q$$ are varying 10 $$-3738$$.133 7496.560 85.40 0 0.1258 0.1067 0.1030 0.0000 1.87E$$-$$05 7.14E$$-$$08 2.63E$$-$$07 0.00043 2.73E$$-$$06 0.00044 All rates are free 12 $${\bf -3693}$$.372 7411.163 0 0.99 0.0649 0.0378 0.3196 0.0398 0.0000 0.2244 1.09E$$-$$05 2.43E$$-$$07 6.45E$$-$$07 0.00224 1.32E$$-$$06 0.00027 Notes:1$$=$$ non-eusocial or societies without worker, 2$$=$$ societies with pseudergates, and 3$$=$$ societies with true workers; $$\lambda =$$ speciation rate (one for each character state, $$\lambda $$1 is the speciation rate for species without worker); $$\mu$$$$=$$ extinction rate (one for each character state, $$\mu 1$$ is the extinction rate for species without worker); and $$q =$$ transition rate between character states. As for the primary cellulolytic partners, we found that the net diversification rate was higher for the Termitidae, the lineage with an entirely prokaryotic gut microbiota (Fig. 3; Table 4; Supplementary Fig. S9 and Table S2 available on Dryad). The two other traits were not significantly different in terms of diversification even though the lineages with cellulolytic flagellates showed higher speciation and extinction rates than mantises and cockroaches (except Cryptocercus that harbors flagellates). Table 4. Supports for gut microbiota as drivers of diversification Models NP logL AICc $$\Delta $$AIC AIC$$\omega $$ $$\lambda 1$$ $$\lambda 2$$ $$\lambda 3$$ $$\mu$$1 $$\mu$$2 $$\mu$$3 $$q$$12 $$q$$13 $$q$$21 $$q$$23 $$q$$31 $$q$$32 Null model 3 $$-384$$0.575 7687.181 274.01 0 0.1959 0.1764 3.03E$$-$$05 All $$\lambda$$ are varying 5 $$-371$$0.247 7430.573 17.41 0 0.1166 0.1270 0.2157 0.0941 1.72E$$-$$05 All $$\mu$$ are varying 5 $$-3731$$.518 7473.115 59.95 0 0.1642 0.1456 0.1373 0.0447 1.21E$$-$$05 All $$q$$ are varying 8 $$-3835$$.418 7687.028 273.86 0 0.1957 0.1763 2.40E$$-$$05 6.75E$$-$$07 4.36E$$-$$06 0.00046 1.88E$$-$$06 3.10E$$-$$06 All$$\lambda $$and$$\mu $$are varying 7 $${\bf -3699}$$.510 7413.168 0 0.703 0.1070 0.1894 0.3022 0.0832 0.1652 0.1979 2.39E$$-$$05 All $$\lambda $$ and $$q$$ are varying 10 $$-37$$03.777 7427.847 14.68 0 0.1147 0.1254 0.2157 0.0921 1.98E$$-$$05 2.74E$$-$$07 3.45E$$-$$06 0.00025 3.21E$$-$$06 2.81E$$-$$06 All $$\mu$$ and $$q$$ are varying 10 $$-3735$$.700 7491.694 78.53 0 0.1683 0.1500 0.1304 0.0555 2.01E$$-$$05 7.33E$$-$$07 3.67E$$-$$06 0.01504 2.54E$$-$$06 0.00220 All rates are free 12 $$-3695$$.237 7414.892 1.72 0.297 0.1037 0.1730 0.3033 0.0795 0.1472 0.1975 3.11E$$-$$05 4.02E$$-$$07 4.87E$$-$$06 0.00029 4.52E$$-$$06 2.54E$$-$$06 Models NP logL AICc $$\Delta $$AIC AIC$$\omega $$ $$\lambda 1$$ $$\lambda 2$$ $$\lambda 3$$ $$\mu$$1 $$\mu$$2 $$\mu$$3 $$q$$12 $$q$$13 $$q$$21 $$q$$23 $$q$$31 $$q$$32 Null model 3 $$-384$$0.575 7687.181 274.01 0 0.1959 0.1764 3.03E$$-$$05 All $$\lambda$$ are varying 5 $$-371$$0.247 7430.573 17.41 0 0.1166 0.1270 0.2157 0.0941 1.72E$$-$$05 All $$\mu$$ are varying 5 $$-3731$$.518 7473.115 59.95 0 0.1642 0.1456 0.1373 0.0447 1.21E$$-$$05 All $$q$$ are varying 8 $$-3835$$.418 7687.028 273.86 0 0.1957 0.1763 2.40E$$-$$05 6.75E$$-$$07 4.36E$$-$$06 0.00046 1.88E$$-$$06 3.10E$$-$$06 All$$\lambda $$and$$\mu $$are varying 7 $${\bf -3699}$$.510 7413.168 0 0.703 0.1070 0.1894 0.3022 0.0832 0.1652 0.1979 2.39E$$-$$05 All $$\lambda $$ and $$q$$ are varying 10 $$-37$$03.777 7427.847 14.68 0 0.1147 0.1254 0.2157 0.0921 1.98E$$-$$05 2.74E$$-$$07 3.45E$$-$$06 0.00025 3.21E$$-$$06 2.81E$$-$$06 All $$\mu$$ and $$q$$ are varying 10 $$-3735$$.700 7491.694 78.53 0 0.1683 0.1500 0.1304 0.0555 2.01E$$-$$05 7.33E$$-$$07 3.67E$$-$$06 0.01504 2.54E$$-$$06 0.00220 All rates are free 12 $$-3695$$.237 7414.892 1.72 0.297 0.1037 0.1730 0.3033 0.0795 0.1472 0.1975 3.11E$$-$$05 4.02E$$-$$07 4.87E$$-$$06 0.00029 4.52E$$-$$06 2.54E$$-$$06 Notes: 1$$=$$ gut with no specialized microbiota for lignocellulose digestion, 2*$$=$$ gut with cellulolytic flagellates, and 3*$$=$$ gut with entirely prokaryotic microbiota; $$\lambda =$$ speciation rate (one for each character state, $$\lambda 1$$ is the speciation rate for species without microbiota); $$\mu$$$$=$$ extinction rate (one for each character state, $$\mu 1$$ is the extinction rate for species without microbiota); and $$q =$$ transition rate between character states. *Results of alternative coding strategies (see main text) are provided in Supplementary Figures S11 and S12 available on Dryad. Table 4. Supports for gut microbiota as drivers of diversification Models NP logL AICc $$\Delta $$AIC AIC$$\omega $$ $$\lambda 1$$ $$\lambda 2$$ $$\lambda 3$$ $$\mu$$1 $$\mu$$2 $$\mu$$3 $$q$$12 $$q$$13 $$q$$21 $$q$$23 $$q$$31 $$q$$32 Null model 3 $$-384$$0.575 7687.181 274.01 0 0.1959 0.1764 3.03E$$-$$05 All $$\lambda$$ are varying 5 $$-371$$0.247 7430.573 17.41 0 0.1166 0.1270 0.2157 0.0941 1.72E$$-$$05 All $$\mu$$ are varying 5 $$-3731$$.518 7473.115 59.95 0 0.1642 0.1456 0.1373 0.0447 1.21E$$-$$05 All $$q$$ are varying 8 $$-3835$$.418 7687.028 273.86 0 0.1957 0.1763 2.40E$$-$$05 6.75E$$-$$07 4.36E$$-$$06 0.00046 1.88E$$-$$06 3.10E$$-$$06 All$$\lambda $$and$$\mu $$are varying 7 $${\bf -3699}$$.510 7413.168 0 0.703 0.1070 0.1894 0.3022 0.0832 0.1652 0.1979 2.39E$$-$$05 All $$\lambda $$ and $$q$$ are varying 10 $$-37$$03.777 7427.847 14.68 0 0.1147 0.1254 0.2157 0.0921 1.98E$$-$$05 2.74E$$-$$07 3.45E$$-$$06 0.00025 3.21E$$-$$06 2.81E$$-$$06 All $$\mu$$ and $$q$$ are varying 10 $$-3735$$.700 7491.694 78.53 0 0.1683 0.1500 0.1304 0.0555 2.01E$$-$$05 7.33E$$-$$07 3.67E$$-$$06 0.01504 2.54E$$-$$06 0.00220 All rates are free 12 $$-3695$$.237 7414.892 1.72 0.297 0.1037 0.1730 0.3033 0.0795 0.1472 0.1975 3.11E$$-$$05 4.02E$$-$$07 4.87E$$-$$06 0.00029 4.52E$$-$$06 2.54E$$-$$06 Models NP logL AICc $$\Delta $$AIC AIC$$\omega $$ $$\lambda 1$$ $$\lambda 2$$ $$\lambda 3$$ $$\mu$$1 $$\mu$$2 $$\mu$$3 $$q$$12 $$q$$13 $$q$$21 $$q$$23 $$q$$31 $$q$$32 Null model 3 $$-384$$0.575 7687.181 274.01 0 0.1959 0.1764 3.03E$$-$$05 All $$\lambda$$ are varying 5 $$-371$$0.247 7430.573 17.41 0 0.1166 0.1270 0.2157 0.0941 1.72E$$-$$05 All $$\mu$$ are varying 5 $$-3731$$.518 7473.115 59.95 0 0.1642 0.1456 0.1373 0.0447 1.21E$$-$$05 All $$q$$ are varying 8 $$-3835$$.418 7687.028 273.86 0 0.1957 0.1763 2.40E$$-$$05 6.75E$$-$$07 4.36E$$-$$06 0.00046 1.88E$$-$$06 3.10E$$-$$06 All$$\lambda $$and$$\mu $$are varying 7 $${\bf -3699}$$.510 7413.168 0 0.703 0.1070 0.1894 0.3022 0.0832 0.1652 0.1979 2.39E$$-$$05 All $$\lambda $$ and $$q$$ are varying 10 $$-37$$03.777 7427.847 14.68 0 0.1147 0.1254 0.2157 0.0921 1.98E$$-$$05 2.74E$$-$$07 3.45E$$-$$06 0.00025 3.21E$$-$$06 2.81E$$-$$06 All $$\mu$$ and $$q$$ are varying 10 $$-3735$$.700 7491.694 78.53 0 0.1683 0.1500 0.1304 0.0555 2.01E$$-$$05 7.33E$$-$$07 3.67E$$-$$06 0.01504 2.54E$$-$$06 0.00220 All rates are free 12 $$-3695$$.237 7414.892 1.72 0.297 0.1037 0.1730 0.3033 0.0795 0.1472 0.1975 3.11E$$-$$05 4.02E$$-$$07 4.87E$$-$$06 0.00029 4.52E$$-$$06 2.54E$$-$$06 Notes: 1$$=$$ gut with no specialized microbiota for lignocellulose digestion, 2*$$=$$ gut with cellulolytic flagellates, and 3*$$=$$ gut with entirely prokaryotic microbiota; $$\lambda =$$ speciation rate (one for each character state, $$\lambda 1$$ is the speciation rate for species without microbiota); $$\mu$$$$=$$ extinction rate (one for each character state, $$\mu 1$$ is the extinction rate for species without microbiota); and $$q =$$ transition rate between character states. *Results of alternative coding strategies (see main text) are provided in Supplementary Figures S11 and S12 available on Dryad. Randomization analyzes of the three data sets of traits showed that all the SSE-based results were robust to type-I error (Supplementary Fig. S10 available on Dryad) indicating that our inferences are not biased. We performed additional MuSSE analyzes by recoding either the wood-feeding non-eusocial lineages or Macrotermitinae with Cryptocercus spp. and all the termites but the other Termitidae, and we found that the results were not sensitive to this effect: eusocial lineages with a prokaryotic microbiota or with a complex hindgut compartmentation diversified significantly faster than their counterparts (Supplementary Figs. S11 and S12 available on Dryad). Discussion The origin of eusociality has been widely investigated at the organismal level. For instance, Hall and Goodisman (2012) explored the differences in rate of substitution for queen- and worker-selected loci, whereas Rehan and Toth (2015) reviewed several hypotheses about molecular evolution of eusociality and tried to integrate them in a synthetic framework. Above the species level, the origin and evolution of eusociality has been also addressed with molecular dated phylogenies of eusocial groups (Farrell et al. 2001; Agnarsson et al. 2006; Moreau et al. 2006; Hines et al. 2007; Cardinal and Danforth 2011, 2013; Pie and Feitosa 2016). These studies have shown that eusociality is a labile trait that has appeared independently in many clades. However, the consequences of eusociality as a driver of species diversification have received less attention (Davis et al. 2009; Ware et al. 2010). What role does eusociality play at this macroevolutionary scale? In social insects, only two studies mentioned this aspect and suggested contradictory hypotheses. Edward O. Wilson (1992) proposed that social insects would have a lower diversification rate than non-eusocial insects, while Davis et al. (2009) suggested that termites would have a higher diversification rate than their close relatives; the study of Ware et al. (2010) did not compare termite diversification with the non-eusocial clades. These opposite predictions are in fact not so surprising because, for numerous aggregate and emergent traits, both low speciation and extinction rates are predicted (Jablonski 2008). These traits would then lead either to increase or decrease in diversification depending on the intensity of speciation and extinction. Speciation rate supposedly increases with molecular rate (Webster et al. 2003; Lanfear et al. 2010; but see Rabosky and Matute 2013), which in turn depends on several features including population size, population abundance, and genetic population. Harvey et al. (2017) showed, for instance, a link between population differentiation within species and speciation rates inferred from phylogenies in New World birds. All the aforementioned features are emergent traits connected to eusociality. In termites, the archetype society mode involves large and abundant populations headed by pairs of reproductives, resulting in colonies with strong genetic structures and genetically isolated by distance (Goodisman and Crozier 2002; Thompson et al. 2007; Dupont et al. 2009; Eggleton 2011; Vargo and Husseneder 2011). According to classical predictions, high abundance would imply low extinction and low speciation rates, while highly structured populations would imply both high extinction and speciation rates. As for large population size, predictions suggest low extinction rates and are contradictory about speciation rates (Jablonski 2008). For Dictyoptera, we found here higher extinction and speciation rates for eusocial lineages when compared to non-eusocial lineages (using three different approaches). Because speciation is much higher than extinction, the net diversification rate is higher for eusocial species than for non-eusocial species, resulting in a strong positive link between eusociality and net diversification. These results suggest that a high level of genetic structure could be the main factor acting on Dictyoptera diversification. However, the classical vision of termites is outdated. Population structure in termites probably shows a large spectrum of variation that should be investigated further (Vargo and Husseneder 2011) and compared to cockroach population structures that should also vary between solitary, gregarious and subsocial species (Bell et al. 2007), or flying and wingless species. Population size—and size-dependent properties—is also poorly known (Korb 2009) but range from a few hundreds to millions of individuals. Even within the family Termitidae comprising only species with “true” workers, there is a large variance with nest population of a few thousands to a few millions of individuals (Lepage and Darlington 2000). But, even though the underlying factors remain to be investigated further, eusociality seems to act as a key innovation favoring diversification in Dictyoptera. Key innovations are alluded to whenever an evolutionary novelty is thought to have favored the diversification of a clade. Three main categories of key innovations can be distinguished (Heard and Hauser 1995): i) innovations opening new adaptive zones sensuSimpson (1944); ii) innovations increasing fitness and allowing to outcompete other species; and iii) innovations increasing the propensity for reproductive or ecological specialization. These three categories are not necessarily exclusive. For termites, eusociality unlikely opened new adaptive zones (Definition 1). Termites act mainly as decomposers, a niche already existing and occupied before the emergence of this clade (Raymond et al. 2000). However, a highly integrated society, especially with “true” workers, might have allowed termites to outcompete other species and increased their propensity for reproductive and ecological specialization (Definitions 2 and 3). This scenario is classically put forward at the ecological timescales, wherein the higher dominance and the large ecological success of termites are known, and could be here translated at the macroevolutionary timescale. Interestingly, all the diversification analyzes (three different types of models) show that both speciation and extinction rates were higher for eusocial lineages with “true” workers resulting in a higher net diversification rate for these lineages. Accordingly, our results suggest that societies with “true” workers are not only more successful at ecological timescales but also over millions of years, which further implies that both organism- and species-level traits (aggregate and emergent traits, respectively) act on species selection. Other factors are also probably at play in Dictyoptera diversification and they should be examined in future research (Jablonski 2008; Benton 2009). The diversification of angiosperms in the Cretaceous is often proposed to have played an important role in insect diversification, especially social insects (Moreau et al. 2006; Cardinal and Danforth 2013). Interestingly, molecular dating phylogenies of termites concur to show that termites appeared in the very Late Jurassic—Early Cretaceous (Ware et al. 2010; Bourguignon et al. 2015; Legendre et al. 2015b). Although there is a debate on the origin of angiosperms, from the Early Cretaceous to Triassic (e.g., Beaulieu et al. 2015; Magallón et al. 2015; Silvestro et al. 2015b; Foster et al. 2017), the main period of diversification always occurred in the Cretaceous and coincides well with the rise of many extant insect clades. It is nonetheless questioned whether the angiosperm radiation has boosted insect diversification (Rainford and Mayhew 2015; Condamine et al. 2016). The increase in termite diversification in the Cretaceous can be associated with the rise to dominance of angiosperms; a similar scenario is proposed for another important eusocial insect clade, the ants (Moreau et al. 2006). Other dramatic landscape modifications that happened in the last 300 million years might have promoted shifts in diversification in termites. For instance, global climate changes or the break-up of ancient continents (Bourguignon et al. 2015) could have induced different dynamics of diversification between non-eusocial and eusocial lineages. Future studies could investigate the impact of such important environmental changes using both phylogenies and fossil data. Also, to better understand these intricate mechanisms, approaches integrating various factors such as caste and colony size but also foraging types, dispersal abilities, diets and digestive capacity (Abe 1987; Jeanson and Deneubourg 2009; Brune and Dietrich 2015), should be used to sort out the effect of these factors. Regarding diet specialization, which is strongly related to gut anatomy and microbiota composition (Bignell and Eggleton 1995; Köhler et al. 2012; Mikaelyan et al. 2015), the concomitance of termite and angiosperm diversifications brings credit to its importance in termite diversification because changes in intestinal anatomy and microbiota complexity would have opened new niches to diversify in. Rabosky (2009) underlined the possible change in species “carrying capacity” for some clades in a given area. Instead of playing on diversification rates per se, may eusociality impact ecological limits in a given environment, which would in turn increase termite diversification? This remains a possibility that could be linked to our side-result that net diversification rates within termites increase towards the present, which is unusual, but the underlying mechanisms are still to be deciphered. The apparent increased diversification rate in termites at the present can be caused by a phenomenon, called the pull of the recent, which is due to a more complete sampling of recent and still extant species than in the past (Jablonski 2008). This is slightly different from the pull of the present that is due to constant birth–death models showing an upward turn in species number towards the present (Nee 2006). Most of the phylogenies are better explained by a slowdown of speciation rates towards the present (e.g., Phillimore and Price 2008; Morlon et al. 2010). An increase in speciation through time is rarely inferred, thus the increase in termite speciation goes against the classical pattern. We do not think the increase in speciation can be attributed to the pull of the present/recent because the increase started tens of million years ago, while the effect of the pull of the present/recent is usually considered to be important in the last million years (Nee 2006). The species sampling is relatively homogeneous across the phylogenetic tree, except for the Panesthiinae which are more sampled than the rest, but more importantly the termite subtree is not better sampled than the other groups (Legendre et al. 2015b). Therefore, we think the impact of the pull of the present/recent is limited here. Finally, even though there were probably multiple origins of the “true” worker caste in termites (Legendre et al. 2008, 2013), Dictyoptera offers a single replicate to investigate the role of eusociality in clade diversification. We expect advances in our understanding of the role of eusociality at the species selection level, when similar studies would be conducted with other eusocial organisms, especially in the insect order Hymenoptera (Peters et al. 2017) where ants, bees, and wasps constitute conspicuous elements of ecosystems and have evolved eusociality independently. Supplementary Material Data available from the Dryad Digital Repository: http://dx.doi.org/10.5061/dryad.2tg34. Funding This work was financially support by a Marie Curie International Outgoing Fellowship [Project 627684 BIOMME to F.L.C.]. Acknowledgments We thank the associate editor Brian Wiegmann and two anonymous referees for the excellent comments that improved the study. We also thank Giovanni Fagua, Gaël Kergoat, Hélène Morlon, André Nel, and Felix Sperling for proofreading and commenting on this manuscript. References Abe T. 1987 . Evolution of life types in termites. In: Kawano S., Connell J.H., Hidaka T., editors. Evolution and coadaptation in biotic communities. Tokyo : University of Tokyo Press . p. 125 – 148 . Agnarsson I., Avilés L., Coddington J.A., Maddison W.P. 2006 . Sociality in theridiid spiders: repeated origins of an evolutionary dead end. Evolution 60 : 2342 – 2351 . Google Scholar CrossRef Search ADS PubMed Beaulieu J.M., O’Meara B.C., Crane P., Donoghue M.J. 2015 . Heterogeneous rates of molecular evolution and diversification could explain the Triassic age estimate for angiosperms. Syst. Biol. 64 : 869 – 878 . Google Scholar CrossRef Search ADS PubMed Beccaloni G.W., Eggleton, P. 2013 . Order Blattodea. Zootaxa 3703 : 46 – 48 . Google Scholar CrossRef Search ADS Bell W.J., Roth L.M., Nalepa C.A. 2007 . Cockroaches: ecology, behavior, and natural history. Baltimore : The Johns Hopkins University Press . Benton M.J. 2009 . The Red Queen and the Court Jester: species diversity and the role of biotic and abiotic factors through time. Science 323 : 728 – 732 . Google Scholar CrossRef Search ADS PubMed Bignell D.E., Eggleton R. 1995 . On the elevated intestinal pH of higher termites (Isoptera: Termitidae). Insectes Soc. 42 : 57 – 69 . Google Scholar CrossRef Search ADS Bourguignon T., Lo N., Cameron S.L., Šobotník J., Hayashi Y., Shigenobu S., Watanabe D., Roisin Y., Miura T., Evans T.A. 2015 . The evolutionary history of termites as inferred from 66 mitochondrial genomes. Mol. Biol. Evol. 32 : 406 – 421 . Google Scholar CrossRef Search ADS PubMed Bourke A. 2011 . Principles of social evolution. Oxford : Oxford University Press . Brune A. 2014 . Symbiotic digestion of lignocellulose in termite guts. Nat. Rev. Microbiol. 12 : 168 – 180 . Google Scholar CrossRef Search ADS PubMed Brune A., Dietrich C. 2015 . The gut microbiota of termites: digesting the diversity in the light of ecology and evolution. Annu. Rev. Microbiol. 69 : 145 – 166 . Google Scholar CrossRef Search ADS PubMed Cardinal S., Danforth B.N. 2011 . The antiquity and evolutionary history of social behavior in bees. PLoS One 6 : e21086 . Google Scholar CrossRef Search ADS PubMed Cardinal S., Danforth B.N. 2013 . Bees diversified in the age of eudicots. Proc. R. Soc. B 280 : 20122686 . Google Scholar CrossRef Search ADS Chapman R.E., Bourke A.F.G. 2001 . The influence of sociality on the conservation biology of social insects. Ecol. Lett. 4 : 650 – 662 . Google Scholar CrossRef Search ADS Cleveland L.R., Hall S.R., Sanders E.P., Collier J. 1934 . The wood feeding roach Cryptocercus, its protozoa, and the symbiosis between protozoa and roach. Mem. Am. Acad. Arts Sci. 17 : 185 – 342 . Condamine F.L., Clapham M.E., Kergoat G.J. 2016 . Global patterns of insect diversification: towards a reconciliation of fossil and molecular evidence? Sci. Rep. 6 : 19208 . Google Scholar CrossRef Search ADS PubMed Condamine F.L., Rolland J., Morlon H. 2013 . Macroevolutionary perspectives to environmental change. Ecol. Lett. 16 : 72 – 85 . Google Scholar CrossRef Search ADS PubMed Darwin C. 1859 . On the origin of species by means of natural selection, or preservation of favored races in the struggle for life. London : Murray . Davis R.B., Baldauf S.L., Mayhew P.J. 2009 . Eusociality and the success of the termites: insights from a supertree of dictyopteran families. J. Evol. Biol. 22 : 1750 – 1761 . Google Scholar CrossRef Search ADS PubMed Davis M.P., Midford P.E., Maddison W. 2013 . Exploring power and parameter estimation of the BiSSE method for analyzing species diversification. BMC Evol. Biol. 13 : 38 . Google Scholar CrossRef Search ADS PubMed Dupont L., Roy V., Bakkali A., Harry M. 2009 . Genetic variability of the soil-feeding termite Labiotermes labralis (Termitidae, Nasutitermitinae) in the Amazonian primary forest and remnant patches. Insect Conserv. Divers. 2 : 53 – 61 . Google Scholar CrossRef Search ADS Eggleton P. 2011 . An introduction to termites: biology, taxonomy and functional morphology. In: Bignell D.E., Roisin Y., Lo N., editors. Biology of termites: a modern synthesis. Dordrecht : Springer . p. 1 – 26 . Engel M.S., Barden P., Riccio M.L., Grimaldi D.A. 2016 . Morphologically specialized termite castes and advanced sociality in the Early Cretaceous. Curr. Biol. 26 : 522 – 530 . Google Scholar CrossRef Search ADS PubMed Farrell B.D., Sequeira A.S., O’Meara B.C., Normark B.B., Chung J.H., Jordal B.H. 2001 . The evolution of agriculture in beetles (Curculionidae: Scolytinae and Platypodinae). Evolution 55 : 2011 – 2027 . Google Scholar CrossRef Search ADS PubMed FitzJohn R.G. 2012 . Diversitree?: comparative phylogenetic analyses of diversification in R. Methods Ecol. Evol. 3 : 1084 – 1092 . Google Scholar CrossRef Search ADS FitzJohn R.G., Maddison W.P., Otto S.P. 2009 . Estimating trait-dependent speciation and extinction rates from incompletely resolved phylogenies. Syst. Biol. 58 : 595 – 611 . Google Scholar CrossRef Search ADS PubMed Foster C.S., Sauquet H., Van der Merwe M., McPherson H., Rossetto M., Ho S.Y. 2017 . Evaluating the impact of genomic data and priors on Bayesian estimates of the angiosperm evolutionary timescale. Syst. Biol. 66 : 338 – 351 . Google Scholar PubMed Goodisman M.A.D., Crozier R.H. 2002 . Population and colony genetic structure of the primitive termite Mastotermes darwiniensis. Evolution . 56 : 70 – 83 . Google Scholar CrossRef Search ADS PubMed Goutte S., Dubois A., Howard S.D., Marquez R., Rowley J.J.L., Dehling J.M., Grandcolas P., Xiong R., Legendre F. 2016 . Environmental constraints and call evolution in torrent dwelling frogs. Evolution 70 : 811 – 826 . Google Scholar CrossRef Search ADS PubMed Grandcolas P. 1993 . Le genre Paramuzoa Roth, 1973: sa répartition et un cas de xylophagie chez les Nyctiborinae (Dictyoptera, Blattaria). Bull. Soc. Entomol. Fr. 98 : 131 – 138 . Grandcolas P., Bellés X., Piulachs M.D., D’Haese C. 2002 . Le genre Lauraesilpha Grandcolas 1997: nouvelles espèces, endémisme, séquences d’ARN ribosomique et caractères d’appartenance aux Blattidae (Insectes, Dictyoptères, Blattidae, Tyronicinae). In: Najt J. & Grandcolas P., editors. Zoologia neocaledonica 5 , vol. 187 . Syst. Mémoires du Muséum Natl. Paris: d’Histoire Nat. p. 117 – 131 . Grandcolas P., D’Haese C. 2002 . The origin of a “true” worker caste in termites: phylogenetic evidence is not decisive. J. Evol. Biol. 15 : 885 – 888 . Google Scholar CrossRef Search ADS Hall D.W., Goodisman M.A.D. 2012 . The effects of kin selection on rates of molecular evolution in social insects. Evolution 66 : 2080 – 2093 . Google Scholar CrossRef Search ADS PubMed Harvey M.G., Seeholzer G.F., Smith B.T., Rabosky D.L., Cuervo A.M., Brumfield R.T. 2017 . Positive association between population genetic differentiation and speciation rates in New World birds. Proc. Natl Acad. Sci. USA 114 : 6328 – 6333 . Heard S.B., Hauser D.L. 1995 . Key evolutionary innovations and their ecological mechanisms. Hist. Biol. 10 : 151 – 173 . Google Scholar CrossRef Search ADS Hines H.M., Hunt J.H., O’Connor T.K., Gillespie J.J., Cameron S.A. 2007 . Multigene phylogeny reveals eusociality evolved twice in vespid wasps. Proc. Natl Acad. Sci. USA 104 : 3295 – 3299 . Inward D., Vogler A.P., Eggleton P. 2007 . A comprehensive phylogenetic analysis of termites (Isoptera) illuminates key aspects of their evolutionary biology. Mol. Phylogenet. Evol. 44 : 953 – 967 . Google Scholar CrossRef Search ADS PubMed Jablonski D. 2008 . Species selection: theory and data. Annu. Rev. Ecol. Evol. Syst. 39 : 501 – 524 . Google Scholar CrossRef Search ADS Jeanson R., Deneubourg J.-L. 2009 . Positive feedback, convergent collective patterns and social transitions in arthropods. In: Gadau J., Fewell J., editors. Organization of insect societies: from genome to sociocomplexity. Cambridge, MA : Harvard University Press . p. 460 – 482 . Johnstone R.A., Cant M.A., Field J. 2012 . Sex-biased dispersal, haplodiploidy and the evolution of helping in social insects. Proc. Biol. Sci. 279 : 787 – 93 . Google Scholar CrossRef Search ADS PubMed Klass K.-D., Nalepa C.A., Lo N. 2008 . Wood-feeding cockroaches as models for termite evolution (Insecta: Dictyoptera): Cryptocercus vs. Parasphaeria boleiriana. Mol. Phylogenet. Evol. 46 : 809 – 817 . Google Scholar CrossRef Search ADS PubMed Köhler T., Dietrich C., Scheffrahn R.H., Brune A. 2012 . High-resolution analysis of gut environment and bacterial microbiota reveals functional compartmentation of the gut in wood-feeding higher termites (Nasutitermes spp.). Appl. Environ. Microbiol. 78 : 4691 – 4701 . Google Scholar CrossRef Search ADS PubMed Korb J. 2009 . Termites: an alternative road to eusociality and the importance of group benefits in social insects. In: Gadau J., Fewell J., editors. Organization of insect societies: from genome to sociocomplexity. Cambridge, MA : Harvard University Press. p. 128 – 147 . Krishna K., Grimaldi D.A., Krishna V., Engel M.S. 2010 . Treatise on the isoptera of the world. Bull. Am. Mus. Nat. Hist. 377 : 1 – 2681 . Google Scholar CrossRef Search ADS Lanfear R., Ho S.Y.W., Love D., Bromham L. 2010 . Mutation rate is linked to diversification in birds. Proc. Natl Acad. Sci. USA 107 : 20423 – 20428 . Lanfear R., Kokko H., Eyre-Walker A. 2014 . Population size and the rate of evolution. Trends Ecol. Evol. 29 : 33 – 41 . Google Scholar CrossRef Search ADS PubMed LaPolla J.S., Dlussky G.M., Perrichot V. 2013 . Ants and the fossil record. Annu. Rev. Entomol. 58 : 609 – 630 . Google Scholar CrossRef Search ADS PubMed Legendre F., Deleporte P., Depraetere M., Gasc A., Pellens R., Grandcolas P. 2015a . Dyadic behavioural interactions in cockroaches (Blaberidae): ecomorphological and evolutionary implications. Behaviour 152 : 1229 – 1256 . Google Scholar CrossRef Search ADS Legendre F., Nel A., Svenson G.J., Robillard T., Pellens R., Grandcolas P. 2015b . Phylogeny of Dictyoptera: dating the origin of cockroaches, praying mantises and termites with molecular data and controlled fossil evidence. PLoS One 10 : e0130127 . Google Scholar CrossRef Search ADS Legendre F., Whiting M.F., Bordereau C., Cancello E.M., Evans T.A., Grandcolas P. 2008 . The phylogeny of termites (Dictyoptera: Isoptera) based on mitochondrial and nuclear markers: implications for the evolution of the worker and pseudergate castes, and foraging behaviors. Mol. Phylogenet. Evol. 48 : 615 – 627 . Google Scholar CrossRef Search ADS PubMed Legendre F., Whiting M.F., Grandcolas P. 2013 . Phylogenetic analyses of termite post-embryonic sequences illuminate caste and developmental pathway evolution. Evol. Dev. 15 : 146 – 157 . Google Scholar CrossRef Search ADS PubMed Lepage M., Darlington J.P.E.C. 2000 . Population dynamics of termites. In: Abe, T. Bignel, D.E., Higashi, M., editors. Termites: evolution, sociality, symbioses, ecology. Dordrecht : Springer. p. 333 – 361 . Liow L.H., Reitan T., Harnik P.G. 2015 . Ecological interactions on macroevolutionary time scales: clams and brachiopods are more than ships that pass in the night. Ecol. Lett. 18 : 1030 – 1039 . Google Scholar CrossRef Search ADS PubMed Maddison W.P., Midford P.E., Otto S.P. 2007 . Estimating a binary character’s effect on speciation and extinction. Syst. Biol. 56 : 701 – 710 . Google Scholar CrossRef Search ADS PubMed Magallón S., GómezAcevedo S., SánchezReyes L.L., HernándezHernández T. 2015 . A metacalibrated timetree documents the early rise of flowering plant phylogenetic diversity. New Phytol. 207 : 437 – 453 . Google Scholar CrossRef Search ADS PubMed Mikaelyan A., Dietrich C., Köhler T., Poulsen M., Sillam-Dussès D., Brune A. 2015 . Diet is the primary determinant of bacterial community structure in the guts of higher termites. Mol. Ecol. 24 : 5284 – 5295 . Google Scholar CrossRef Search ADS PubMed Moore B.R., Höhna S., May M.R., Rannala B., Huelsenbeck J.P. 2016 . Critically evaluating the theory and performance of Bayesian analysis of macroevolutionary mixtures. Proc. Natl Acad. Sci. USA 113 : 9569 – 9574 . Moreau C.S., Bell C.D., Vila R., Archibald S.B., Pierce N.E. 2006 . Phylogeny of the ants: diversification in the age of angiosperms. Science 312 : 101 – 104 . Google Scholar CrossRef Search ADS PubMed Morlon H., Lewitus E., Condamine F.L., Manceau M., Clavel J., Drury J. 2016 . RPANDA: An R package for macroevolutionary analyses on phylogenetic trees. Methods Ecol. Evol. 7 : 589 – 597 . Google Scholar CrossRef Search ADS Morlon H., Parsons T.L., Plotkin J.B. 2011 . Reconciling molecular phylogenies with the fossil record. Proc. Natl Acad. Sci. USA 108 : 16327 – 16332 . Morlon H., Potts M.D., Plotkin J.B. 2010 . Inferring the dynamics of diversification: a coalescent approach. PLoS Biol. 8 : e1000493 . Google Scholar CrossRef Search ADS PubMed Nee S. 2006 . Birth-death models in macroevolution. Annu. Rev. Ecol. Evol. Syst. 37 : 1 – 17 . Google Scholar CrossRef Search ADS Noirot C., Pasteels J.M. 1987 . Ontogenetic development and evolution of the worker caste in termites. Experientia 43 : 851 – 952 . Google Scholar CrossRef Search ADS Nonacs P. 2011 . Kinship, greenbeards, and runaway social selection in the evolution of social insect cooperation. Proc. Natl Acad. Sci. USA 108 : 10808 – 10815 . Nowak M.A., Tarnita C.E., Wilson E.O. 2010 . The evolution of eusociality. Nature 466 : 1057 – 1062 . Google Scholar CrossRef Search ADS PubMed Okasha S. 2006 . The levels of selection debate: philosophical issues. Philos. Compass 1 : 74 – 85 . Google Scholar CrossRef Search ADS Pellens R., D’Haese C.A., Bellés X., Piulachs M.D., Legendre F., Wheeler W.C., Grandcolas P. 2007 . The evolutionary transition from subsocial to eusocial behaviour in Dictyoptera: phylogenetic evidence for modification of the “shift-in-dependent-care” hypothesis with a new subsocial cockroach. Mol. Phylogenet. Evol. 43 : 616 – 626 . Google Scholar CrossRef Search ADS PubMed Pellens R., Grandcolas P., da Silva-Neto I.D. 2002 . A new and independently evolved case of xylophagy and the presence of intestinal flagellates in the cockroach Parasphaeria boleiriana (Dictyoptera, Blaberidae, Zetoborinae) from the remnants of the Brazilian Atlantic forest. Can. J. Zool. 80 : 350 – 359 . Google Scholar CrossRef Search ADS Peters R.S., Krogmann L., Mayer C., Donath A., Gunkel S., Meusemann K., Kozlov A., Podsiadlowski L., Petersen M., Lanfear R., Diez P.A., Heraty J., Kjer K.M., Klopfstein S., Meier R., Polidori C., Schmitt T., Liu S., Zhou X., Wappler T., Rust J., Misof B., Niehuis O. 2017 . Evolutionary history of the Hymenoptera. Curr. Biol. 27 : 1013 – 1018 . Google Scholar CrossRef Search ADS PubMed Phillimore A.B., Price T.D. 2008 . Density-dependent cladogenesis in birds. PLoS Biol. 6 : e71 . Google Scholar CrossRef Search ADS PubMed Pie M.R., Feitosa R.M. 2016 . Relictual ant lineages (Hymenoptera: Formicidae) and their evolutionary implications. Myrmecol. News 22 : 55 – 58 . Rabosky D.L. 2009 . Ecological limits and diversification rate: alternative paradigms to explain the variation in species richness among clades and regions. Ecol. Lett. 12 : 735 – 743 . Google Scholar CrossRef Search ADS PubMed Rabosky D.L. 2014 . Automatic detection of key innovations, rate shifts, and diversity-dependence on phylogenetic trees. PLoS One 9 : e89543 . Google Scholar CrossRef Search ADS PubMed Rabosky D.L., McCune A.R. 2010 . Reinventing species selection with molecular phylogenies. Trends Ecol. Evol. 25 : 68 – 74 . Google Scholar CrossRef Search ADS PubMed Rabosky D.L., Goldberg E.E. 2015 . Model inadequacy and mistaken inferences of trait-dependent speciation. Syst. Biol. 64 : 340 – 355 . Google Scholar CrossRef Search ADS PubMed Rabosky D.L., Matute D.R. 2013 . Macroevolutionary speciation rates are decoupled from the evolution of intrinsic reproductive isolation in Drosophila and birds. Proc. Natl Acad. Sci. USA 110 : 15354 – 15359 . Rabosky D.L., Grundler M., Anderson C., Title P., Shi J.J., Brown J.W., Huang H., Larson J.G. 2014 . BAMMtools: an R package for the analysis of evolutionary dynamics on phylogenetic trees. Methods Ecol. Evol. 5 : 701 – 707 . Google Scholar CrossRef Search ADS Rabosky D.L., Mitchell J.S., Chang J. 2017 . Is BAMM flawed? Theoretical and practical concerns in the analysis of multi-rate diversification models. Syst. Biol. 66 : 477 – 498 . Google Scholar CrossRef Search ADS PubMed Rabosky D.L., Santini F., Eastman J., Smith S.A., Sidlauskas B., Chang J., Alfaro M.E. 2013 . Rates of speciation and morphological evolution are correlated across the largest vertebrate radiation. Nat. Commun. 4 : 1958 . Google Scholar CrossRef Search ADS PubMed Rainford J.L., Mayhew P.J. 2015 . Diet evolution and clade richness in Hexapoda: a phylogenetic study of higher taxa. Am. Nat. 186 : 777 – 791 . Google Scholar CrossRef Search ADS PubMed Rainford J.L., Hofreiter M., Nicholson D.B., Mayhew P.J. 2014 . Phylogenetic distribution of extant richness suggests metamorphosis is a key innovation driving diversification in insects. PLoS One 9 : e109085 . Google Scholar CrossRef Search ADS PubMed Raymond A., Cutlip P., Sweet M. 2000 . Rates and processes of terrestrial nutrient cycling in the Paleozoic: the world before beetles, termites, and flies. In: Allmon W., Bottjer D., editors. Evolutionary paleoecology: the ecological context of macroevolutionary change. New York : Columbia University Press . p. 235 – 283 . Rehan S.M., Toth A.L. 2015 . Climbing the social ladder: the molecular evolution of sociality. Trends Ecol. Evol. 30 : 426 – 433 . Google Scholar CrossRef Search ADS PubMed Roisin Y., Korb J. 2011 . Social organisation and the status of workers in termites. In: Bignell D.E., Roisin Y., Lo N., editors. Biology of termites: a modern synthesis. Dordrecht : Springer . p. 133 – 164 . Rousset F., Lion S. 2011 . Much ado about nothing: Nowak et al.’s charge against inclusive fitness theory. J. Evol. Biol. 24 : 1386 – 1392 . Google Scholar CrossRef Search ADS PubMed Sánchez-García M., Matheny P.B. 2016 . Is the switch to an ectomycorrhizal state an evolutionary key innovation in mushroom-forming fungi? A case study in the Tricholomatineae (Agaricales). Evolution 71 : 51 – 65 . Google Scholar CrossRef Search ADS PubMed Silvestro D., Antonelli A., Salamin N., Quental T.B. 2015a . The role of clade competition in the diversification of North American canids. Proc. Natl Acad. Sci. USA 112 : 8684 – 8689 . Silvestro D., CascalesMiñana B., Bacon C.D., Antonelli A. 2015b . Revisiting the origin and diversification of vascular plants through a comprehensive Bayesian analysis of the fossil record. New Phytol. 207 : 425 – 436 . Google Scholar CrossRef Search ADS Simpson G.G. 1944 . Tempo and Mode in Evolution. New York : Columbia University . Thompson G.J., Kitade O., Lo N., Crozier R.H. 2000 . Phylogenetic evidence for a single, ancestral origin of a “true” worker caste in termites. J. Evol. Biol. 13 : 869 – 881 . Google Scholar CrossRef Search ADS Thompson G.J., Lenz M., Crozier R.H., Crespi B.J. 2007 . Molecular-genetic analyses of dispersal and breeding behaviour in the Australian termite Coptotermes lacteus: evidence for non-random mating in a swarm-dispersal mating system. Aust. J. Zool. 55 : 219 – 227 . Google Scholar CrossRef Search ADS Vargo E.L., Husseneder C. 2011 . Genetic structure of termite colonies and populations. In: Bignell D.E., Roisin Y., Lo N., editors. Biology of termites: a modern synthesis. Dordrecht : Springer . p. 321 – 347 . Ware J.L., Grimaldi D.A., Engel M.S. 2010 . The effects of fossil placement and calibration on divergence times and rates: an example from the termites (Insecta: Isoptera). Arthropod Struct. Dev. 39 : 204 – 219 . Google Scholar CrossRef Search ADS PubMed Watson J.A.L., Sewell J.J. 1985 . Caste development in Mastotermes and Kalotermes: which is primitive? In: Watson J.A.L., Okot-Kotber B.M., Noirot C., editors. Caste differentiation in social insects. Book Section. Oxford : Pergamon . p. 27 – 40 . Webster A.J., Payne R.J.H., Pagel M. 2003 . Molecular phylogenies link rates of evolution and speciation. Science 301 : 478 – 478 . Google Scholar CrossRef Search ADS PubMed Wilson E.O. 1971 . The insect societies. Cambridge (MA) : The Belknap Press of the Harvard University Press . 548 p. Wilson E.O. 1990 . Success and dominance in ecosystems: the case of the social insect. Oldendorf-L, Germany: Ecology Institute. 104 p. Wilson E.O. 1992 . The effects of complex social life on evolution and biodiversity. Oikos 63 : 13 – 18 . Google Scholar CrossRef Search ADS Wilson E.O., Hölldobler B. 2005 . Eusociality: origin and consequences. Proc. Natl Acad. Sci USA 102 : 13367 – 13371 . © The Author(s) 2018. Published by Oxford University Press, on behalf of the Society of Systematic Biologists. All rights reserved. For Permissions, please email: journals.permissions@oup.com This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://academic.oup.com/journals/pages/about_us/legal/notices)
Journal
Systematic Biology – Oxford University Press
Published: Sep 1, 2018