Around the world in 10 million years: Rapid dispersal of a kleptoparasitoid spider wasp (Pompilidae: Ceropales)

Cosmopolitan distributions have classically been explained by Pangaean vicariance. However, evidence of recently diverged cosmopolitan groups has re‐opened consideration on the processes involved. Our aim is to estimate the processes leading to the worldwide distribution of the kleptoparasitoid genus Ceropales.


| INTRODUC TI ON
The processes leading to worldwide distributions have been the subject of long-standing discussion among biogeographers. Early explanations are attributed to Darwin who proposed that widespread taxa dispersed long distance from a single center of origin and geographical barriers fuelled natural selection and induced speciation (Darwin, 1859). This dispersal paradigm was further refined and continued to be supported (Bremer, 1992;Brundin, 1966;Darlington, 1957;Erwin, 1981;Fleagle, 1999;Hennig, 1966;Matthew, 1914;Wallace, 1876). It has nevertheless been widely criticized due to the lack of falsifiable hypotheses (Ball, 1975;Morrone & Crisci, 1995;) and lack of evidence for possible ways of dispersal such as land bridges (Fichman, 1977;Simpson, 1940;Stankiewicz et al., 2006). The recognition of Earth systems leading to geographical instability (Wegener, 1915;Wilson, 1963) gave rise to an alternative model that incorporated continental drift and sea-level fluctuations to explain the distribution of fossils and living organisms due to dynamism of the planet (Charig, 1971;Colbert, 1973;Cracraft, 1974;Croizat, 1968): the vicariance theory (Croizat, 1958;Schram, 1977). Strictly speaking, vicariance occurs when the formation of a barrier divides an ancestral species' range and fragments its population, leading to genetic divergence and speciation.
Dispersal across such a barrier was assumed to be highly improbable and deemed to be a rare chance event (Ball, 1975;Croizat et al., 1974;Humphries, 2000;Nelson, 1978).
While vicariance biogeography appears to undermine the traditional dispersal paradigm, vicariant explanations are not dependent upon the rejection of dispersal patterns per se (Rosen, 1975;. Platnick and Nelson (1978) proposed that certain kinds of dispersal (i.e., range expansion or 'diffusion') actually imply a "hidden" vicariance when a species' range is enlarged and subsequently divided. These areas can have a secondary overlap of biotas by the collapse of the barrier. The strict definition of dispersal, therefore, is a long-distance dispersal, often over inhospitable habitat (Rosen, 1978). Due to the unlikeliness of dispersal across substantial barriers, many taxa -in particular, cosmopolitan taxa -had been assumed to be of Pangean, Gondwanan, or Laurasian origin (e.g., Craw, 1985;Jaeger & Martin, 1984;Patterson, 1981a;1981b;Rosen, 1975;Tiffney, 1985;Wiley, 1988). The outcome of such a purely vicariant history underlying tectonic shifts may be represented by sister relationships between taxa from areas that were adjacent before the continental break-up ( Figure 1).
Although vicariance analysis provides testable hypotheses, it is limited by the necessity for detailed phylogenies (Weitzman & Weitzman, 1982). With the advent of molecular phylogenetics and divergence time estimation, we can now muster robust phylogenies, but the question arises whether we can falsify vicariance hypotheses using a range of possible ages for processes underlying cladogenic events. Even if the branching pattern of a lineage is congruent with the vicariant fragmentation of an area (i.e., the geological history of continents), the age of the lineage also must be congruent with geology. This is known as pseudo-congruence (Moore & Donoghue, 2003).
While on the one hand recent divergence time estimation analyses have validated vicariance hypotheses (Bakkes et al., 2018;Jerome et al., 2014;Joshi & Edgecombe, 2019;Toussaint et al., 2016), on the other hand it could be shown that they have falsified vicariance explanations in cases showing that dispersals across even large barriers like oceans have taken place (McDowall, 2002;Raxworthy et al., 2002;Schrago & Russo, 2003;Yoder et al., 2003). Furthermore, longdistance dispersal has been demonstrated for certain cosmopolitan taxa with disjunct distributions that diverged only in the Tertiary (Ceccarelli et al., 2016) or that are considered less mobile (Rota et al., 2016;Trewick, 2000;Ward, 2014). Perhaps recent widespread distributions may be more dependent on a species' ability to colonize new niches than on dispersal ability. For example, generalist species have better dispersal ability because a plastic diet lets the them exploit a variety of food resources in new environments (Dennis et al., 2000;Jønsson et al., 2016). The evolution of a generalist parasitic life style could perhaps be a strategy that may facilitate long-distance dispersal (Grégoir et al., 2015). Pompilidae, or spider wasps, are parasitoids of spiders, which usually paralyze, oviposit, and transport a single host to a nesting site.
Pompilidae are generally not obligate parasites of a specific species, but instead parasitize a range of spiders with a particular ecological niche. Being linked to the host's ecology implies certain limitations to their diversification and their ability to disperse over large distances. For these wasps it has been shown that host switching may provide new niches and is correlated with higher diversification rates (Rodriguez et al., 2016). A few genera within the family have developed the ability to parasitize other spider wasps' hosts as so-called kleptoparasitoids (Evans, 1953). Within those, Ceropales (Latreille) are unique compared to other pompilid kleptoparasitoids in that they lay their egg in another spider wasp's host before it is placed in a nesting site, while its being transported (O'Neill, 2001). Compared to other pompilids, kleptoparasitoids may have a higher chance of colonizing new niches and therefore dispersing more rapidly because of a broader niche that does not require them to adapt to the spider's defense.
For recently diverged taxa such as cosmopolitan Ceropales, the Pangaean vicariance hypothesis can be falsified and long-distance dispersal may be hypothesized. It is intriguing, however, how such small organisms have conquered the globe in such a short period of time. Here, we study the historical biogeography of Ceropales spider wasps using molecular data and divergence time estimation methods. Our aim is to resolve the timing of Ceropales' diversification, as well as its ancestral range. Most importantly, we aim to reconstruct the dispersal events that led to its current worldwide distribution.

| Molecular methods
DNA was extracted from pinned specimens using the Roche DNA Isolation Kit for Cells and Tissues (Roche Molecular Systems, Inc.).
The specimen pin was removed and specimens placed in a proteinase K and lysis buffer solution from 24 to 48 h. Further steps of the extraction followed the instruction manual. Amplification and sequencing of long-wavelength rhodopsin (LWRh), and the D2-D3 regions of the 28S ribosomal RNA (28S) followed the methods outlined in Pilgrim and Pitts (2006). Amplification of cytochrome c oxidase I followed the methods in Rodriguez et al. (2014) (for primers used see Table S1 in Appendix S1). Raw sequencing reads were processed and assembled with Sequencher 4.1 (Gene Codes Corp.).

| Phylogenetic and dating analyses
Alignment of each marker was performed with ClustalW (Thompson et al., 1994) implemented in Geneious 6.1 (created by Biomatters, available at http://www.genei ous.com) and manually refined. Introns were removed from the LWRh alignment. All markers were concatenated into a supermatrix and the model of molecular evolution for each marker and codon position (the latter for LWRh and COI) was determined using PartitionFinder 1.01 (Lanfear et al., 2012; http://www.rober tlanf ear.com/parti tionf inder/). The combined data phylogeny was inferred through a Bayesian analysis in MrBayes 3.2 (Huelsenbeck & Ronquist, 2001) by running two independent MCMC (Markov Chain Monte Carlo) for 100,000,000 generations. Chain convergence and ESS (Estimated Sample Size) were assessed in tracer 1.5 (Rambaut et al., 2003). The first 10% generations were removed as burn-in. Single-gene topologies were also estimated under the best-fit model and visually inspected to detect possible topological incongruencies. An uncorrelated lognormal relaxed-clock model (Drummond et al., 2006;Drummond & Rambaut, 2007) was used for divergence time estimation in Beast 1.7.5 (Drummond et al., 2012) using the partitioned supermatrix. Substitution models were unlinked among partitions (see Table S2 in Appendix S1), with the underlying clock and trees linked. Because of the lack of fossils in Ceropalinae, we performed a secondary calibration using previously published ages for Ceropalinae and Ceropales . Hard minimum ages were placed as priors for the age of crown-group Ceropalinae (lognormal distribution, mean in real space = 8.5, stdev = 1, offset = 13.3) and Ceropales (lognormal distribution, mean in real space = 8.3, SD = 1, offset = 5.3) by setting the offset of the lognormal distribution to the lower bound of the 95% HPD (Highest Posterior Density) reported by Waichert, Rodriguez, et al. (2015). This approach allowed for the highest density of the prior probability to be placed close to the minimum age encountered by Waichert, Rodriguez, et al. (2015), allowing the MCMC to search for older dates at a lower density. The MCMC search was run for 100,000,000 generations.
Chain convergence and ESS were assessed with tracer 1.5. Independent runs were assembled with LogcoMBiner 1.7.5 (Drummond et al., 2012).
The first 10% of the generations were discarded as burn-in. We performed an event-based likelihood ancestral-range estimation implementing three models of range evolution in BioGeoBEARS (Matzke, 2014). and versions of these models incorporating jump dispersal for a total of six models. The models evaluated were as follows: (1) (Ronquist, 1996); (4) DIVALIKE+J, including jump dispersal (Matzke, 2014); (5) BAYAREAlike, a likelihood implementation of the process assumptions of BayArea, originally Bayesian (Landis et al., 2013); and (6) BAYAREAlike+J, allowing jump-dispersal (Matzke, 2014). All six models were evaluated under a constrained analysis, where dispersal to nonadjacent areas was not allowed and an unconstrained analysis allowing all possible dispersal events; we therefore evaluated 12 possible scenarios. These constraints were placed mainly because (1) there is a very low probability of interoceanic dispersal, (2) the small ranges of spider wasps are dictated by their limited dispersal ability, and (3) the age of Pompilidae suggests that younger areas do not include disjunct distributions because they have been isolated since the origin of the group (see Rodriguez et al., 2015). The log-likelihood of each pair of nested models (i.e., DEC and DEC+j, DIVAlike and DIVAlike+j, BAYAREAlike and BAYAREAlike+j) was compared using a likelihood-ratio test. In addition, the AIC (Akaike Information Criterion, Burnham & Anderson, 1998) was calculated for all models. Using the model with the best AIC score, we performed Biogeographic Stochastic Mapping (BSM) in BioGeoBEARS (Dupin et al., 2017) to simulate possible biogeographic processes. We simulated 200 BSM using the DEC+J model and the consensus tree for the BEAST analysis. Simulated processes were extracted to determine the process with the highest support.

| Phylogenetic and dating analyses
The best partitioning strategy for the combined dataset included six partitions (see Table S2 in Appendix S1). Overall, we recovered high posterior probabilities (PP) for the majority of nodes (PP > 0.95, Figure 2). Subgenera and species groups currently used in Ceropales's classification (Móczár, 1994) were mostly paraphyletic except for Priesnerius, which was represented by a single species (Figure 1). The taxonomic implications of these results will be discussed elsewhere.
Our results suggested that Ceropales had its origin in the Tortonian, early Miocene (10.6 M, 15.7-6.5 95% HPD), and recovered five major clades including taxa found in shared areas of endemism.
These clades originated in the Late Miocene or Pliocene. The Oriental region was the only area not represented by a monophyletic group as the inclusion of the Australian taxa renders it paraphyletic. The mean ages of the nodes and their 95% confidence intervals (CI) for the highest posterior density (HPD) are shown in Figure S1 of Appendix S1.

| Historical biogeography
The constrained DEC+J analysis produced the highest likelihood scores and the lowest AIC and AICc scores (Table 1, Figure 2). As a conservative approach, and to account for error and missing lineages (extinct or species not sampled), ages for events are reported as an interval of time from the stem-group median age to the crown-group median age, and the 95% HPD are reported as the upper bound 95% HPD from the stem group and the lower bound 95% HPD from the crown group.

| DISCUSS ION
With an origin in the Miocene, Ceropales constitutes a young group, whose worldwide distribution cannot be explained by vicariance events caused by continental drift (Figure 1), and its range expansion F I G U R E 2 Estimate of the biogeography history of Ceropales using the best-fit DEC+J (Dispersal Extinction Cladogenesis with Jump dispersal) model. Nodes with Bayesian posterior probability (BPP) between 0.85 and 0.94 are denoted with an asterisk BPP between 0.75 and 0.84 are denoted by black circles. All remaining nodes had BPP = 1. The legend on the bottom left corner represents the five areas shown in the map (Aitoff projection), and an additional area, composed by North America and Eurasia. For each node, pie chart slices represent ancestral area probabilities resulting from the DEC+J analysis. The time scale is in millions of years. The bar graph shows the relative proportion of different cladogenetic events at each 1-million-year time slot.  (Waichert, Rodriguez, et al., 2015). Suitable hosts were already widespread throughout the continents when Ceropales began to diversify.
Colonization was probably favored by the already diversified hosts, which reduced limiting factors such as food resource and nest construction.

| The origin of Old World-New World disjunct distributions
The

Mesoamerica and South America
Four independent dispersal events to South America were obtained through our analyses. These events are inferred to be younger than had been obtained for other pompilids Waichert, Rodriguez, et al., 2015), where dispersal occurred between

| Dispersal to the Ethiopian region
The presence of worldwide-distributed taxa in the Ethiopian region has classically been explained by vicariance events that split a widespread Gondwanan population. The expected phylogenetic pattern for such event is a sister relationship between Ethiopian and Neotropical clades

TA B L E 1
Relative probability of each model implemented in BIOGEOBEARS using Akaike Information Criterion (AIC) and AICc, which controls for sample size.  (Riccieri et al., 2017), the Pleistocene has been a highly reported time for dispersal between the Ethiopian and Palaearctic regions (Frajman & Schönswetter, 2017) for a number of terrestrial invertebrates, such as spiders, land snails, and hoverflies (Opatova et al., 2016;Sherpa et al., 2018;Ståhls et al., 2016).

| Dispersal to the Australian region
Our results show a sister relationship between taxa from the Eurasian and Australian regions, with a recent divergence in the Pleistocene ( Figure 2). The observed pattern is incongruent with the expected.
The Eurasian region has classically been observed as sister to North America as a result of the break-up of Laurasia ( Figure 1). Under this same scenario, the Australian insect fauna is a mix between relictual and disjunct elements (Yeates & Cassis, 2017), with origins from Gondwanan break-up, dispersal through temporary land bridges or long-distance dispersal ( Figure 1). Recent divergences as a result of dispersal in the Pleistocene are consistent with the drying of extensive areas, such as the Sunda Plains, during glacial cycles (Lohman et al., 2011). These large, exposed areas may have served as land bridges for dispersal between the Oriental and Australian regions. This route has been proposed as the main dispersal path used by humans ca.
65,000 years ago, thousands of years before the last glacial maximum (De Deckker et al., 2019). In the Pliocene-Pleistocene, this route of dispersal has been documented for invertebrate groups such as land snails (Köhler & Criscione, 2013) and moths (Spitsyn et al., 2018).

| CON CLUS IONS
This is the first study to address the biogeographical processes that produced the current distributions of a cosmopolitan pompilid lineage.
In general, our results fit previously suggested hypotheses for the processes underlying these distributions. Ceropales is the youngest widespread genus within Pompilidae (Waichert, Rodriguez, et al., 2015). Females exhibit kleptoparasitoid behavior, where adults oviposit on spiders already paralyzed by other pompilids. This behavior could be advantageous by reducing limiting colonization factors, for example, adaption to the spider's defenses. This could result in an expanded niche and allowed the group to colonize, which may have facilitated dispersal through all continents in a short period of time. and an anonymous reviewer for significantly improving this manuscript. We also thank Jaime Florez (ANIC, CSIRO) for assistance with R script development. Specimens used for this study were sampled from the existing EMUS collection. Information on specimen date and place of collection as well as who collected them can be found with the Dryad dataset at https://doi.org/10.5061/dryad.c59zw 3r5f.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data that support the findings of this study are openly available in Dryad, https://doi.org/10.5061/dryad.c59zw 3r5f.