Hybridizing salamanders experience accelerated diversification

Whether hybridization generates or erodes species diversity has long been debated, but to date most studies have been conducted at small taxonomic scales. Salamanders (order Caudata) represent a taxonomic order in which hybridization plays a prevalent ecological and evolutionary role. We employed a recently developed model of trait-dependent diversification to test the hypothesis that hybridization impacts the diversification dynamics of species that are currently hybridizing. We find strong evidence supporting this hypothesis, showing that hybridizing salamander lineages have significantly greater net-diversification rates than non-hybridizing lineages. This pattern is driven by concurrently increased speciation rates and decreased extinction rates in hybridizing lineages. Our results support the hypothesis that hybridization can act as a generative force in macroevolutionary diversification.

A leading unresolved question in evolutionary biology is whether hybridization, defined as the interbreeding between two genetically distinct lineages 1 , acts as a creative or destructive evolutionary force [2][3][4][5][6] . The prevailing view in the animal literature is that hybridization constrains lineage diversification because hybrid lineages are often documented to be less fit than parentals 4,5,7 . Under this scenario, hybridization is predicted to increase extinction rates. Further, introgressive hybridization has the potential to "wash away" accumulating divergence among incompletely isolated lineages [8][9][10] , leading to a prediction of decreased speciation rates. In contrast, the prevailing view in the plant literature is that hybridization enhances adaptive potential by introducing novel genetic and phenotypic variation [11][12][13] . Reinforcement, or the accumulation of post-zygotic reproductive isolation through selection against hybrids 14 , has long been considered to expedite the speciation of diverging lineages 9 . Reinforcement has the potential to act even when heterosis occurs, as hybrids may be largely sterile, such as in Triturus newts 15,16 . Additionally, hybridization-mediated shuffling of old genetic variants may fuel rapid diversification, as outlined by the combinatorial view of speciation 17 . Accordingly, hybridization is predicted to increase speciation rates and/or decrease extinction rates.
Whereas evidence of hybrid speciation in plants has long been abundant, evidence for widespread hybrid speciation in animals is relatively scarce 6,10,[18][19][20][21] . In allopolypoid hybrid speciation, a mode of speciation common in plants but rare in animals, nearly complete reproductive isolation may evolve in a single generation due to a change in ploidy 11 . In contrast, homoploid hybrid speciation typically has to occur in the face of continued gene flow, which acts to homogenize the diverging hybrid lineages 10 .
Interest in hybrid-mediated speciation has recently burgeoned, but studies have typically been limited in taxonomic scope. Studies of the effect of hybridization on diversification have most commonly been conducted among closely related pairs or small clades of taxa, and results have been equivocal 1,[22][23][24][25][26][27][28] . Additionally, studies of hybridization often occur along different stages of the speciation continuum 29 , whereby hybridization can appear as a force that either facilitates or impedes speciation. Recent work highlights this uncertainty by demonstrating that the outcomes of hybridization depend on the underlying nature of selection pressures and demography 30 . We suggest that studies at broad phylogenetic and macroevolutionary timescales can help overcome these limitations by providing a phylogenetic context in which to view repeated hybridization events over evolutionary timescales and the consequent impact on lineage diversification rates.
www.nature.com/scientificreports www.nature.com/scientificreports/ Here, we conduct a taxonomically-broad test of the relationship of hybridization with macroevolutionary diversification rates. We study salamanders (order Caudata, ca. 716 spp as of October 2018 31 ), which are particularly suitable for this study because hybridization is pervasive and has been studied extensively (i.e., nearly 1/3 of N. American species hybridize: Supplementary Fig. S1, Table S1). Additionally, sufficient sequence data are available to resolve the phylogenetic relationships among most (~63%) taxa within this group 32 .
If hybridization plays a meaningful role in the diversification process, differences in diversification rates among hybridizing and non-hybridizing taxa are expected. Thus, we test the hypothesis that there is a difference in diversification rates (speciation and extinction rates) between contemporaneously hybridizing and non-hybridizing salamander lineages. Note that we are simply testing whether contemporary hybridization influences diversification rates, not whether ancient hybridization facilitated the present radiation as postulated by the hybrid swarm hypothesis 6,52 , because our experimental design cannot address this (See supplement). We replicate this test across four datasets to investigate the robustness of our results: (1) including all available data; (2) excluding species that do not exhibit sympatry (defined as <10% geographic range overlap) thereby lacking the opportunity to hybridize; (3) only the family Plethodontidae, which are the most widely hybridizing and diverse of the 10 salamander families; (4) all (nine) salamander families except the Plethodontidae.

Materials and Methods
Data collection. We used the time-calibrated phylogeny of Amphibia 32 as the source for downstream analyses. This tree 32 was constructed using nine nuclear genes and three mitochondrial genes as data using RAxMLv7.2.8 53 , and time-calibrated using treePL 54 . Using the APE package (v.3.4 55 ) in R, we extracted the subtree containing salamanders for subsequent analyses. This approach yielded a tree containing 469 of 716 extant species of salamanders ( Supplementary Fig. S1) 32 representing approximately two-thirds of the known diversity.
Each species was scored as 'non-hybridizable' (NH) or 'hybridizable' (H), based on an extensive literature review of hybridization using the search engines Google Scholar and ISI Web of Knowledge. To do so we paired each species with the terms 'hybrid' and 'introgress' as well as dialectical and structural variants thereof (e.g., 'hybridization' , 'introgressed') to search for cases of hybridization. Criteria for hybridizability under the first, "narrow" definition were the documentation of hybridization among natively distributed species, specifically the observation of heterospecific mating or hybrid offspring in the wild as detected by substantial and repeated morphological and molecular intermediacy. Criteria for hybridizability under the second, "broad" definition included those criteria described previously, as well as: 1) the observation of hybridization occurring in laboratory settings or among introduced and native species, and/or; 2) inference of historical introgression as determined using molecular lines of evidence (e.g. substantial and replicated genealogical discordance among molecular markers or detection via Approximate Bayesian Computation methods). This latter definition is less conservative than the former. We excluded cases where one or two individuals exhibited cyto-nuclear discordance between markers such as allozymes and mtDNA (as observed for members of Sirenidae 56 , and several Hynobiids 57-60 ), as we did not perceive these cases to meet the criteria of substantial and replicated genealogical discordance. Although inference of hybridization via the detection of genealogical discordance warrants caution, our narrow definition of hybridization does not recognize these species as hybridizable. In total, we retrieved 56 papers (date-range: 1957-2017: Supplementary Table S1). We confirmed hybridization for roughly 11 and 13 percent of extant Caudates for the narrow and broad datasets, respectively (78 and 92 of 716 species: 17 and 20 percent of sampled taxa; Supplementary Table S1). Documented hybridization was absent from four families (Cryptobranchidae, Sirenidae, Proteidae, Rhyacotritonidae) in the narrow dataset and two families (Rhyacotritonidae, Sirenidae) in the broad dataset.
In addition to analyzing the entire salamander subclade from the Pyron 32 amphibian phylogeny (469 total species), we compared diversification rates using three additional datasets, produced using both the narrow and broad datasets (i.e., for a total of eight datasets). In the first dataset, we required that "non-hybridizable" taxa were sympatric (>10% range overlap) with another species, thus possessing sufficient opportunity to hybridize. Consequently, species that were not sympatric with any other salamander taxa were excluded from the dataset. Therefore, species classified as non-hybridizing in this analysis may have limited opportunity to hybridize but have not been observed to do so. Although we cannot account for historical species distributions (e.g. species' ranges may have previously overlapped), our primary (narrow) definition of hybridizability necessitates contemporary hybridization. Thus, our designation of sympatry occurs at the same time scale as our designation of hybridizability.
In the second and third datasets, we tested for a family-specific effect of hybridization on diversification rates. As plethodontids are the most diverse family of salamanders (Fig. 1B) and have been the subject of extensive study 37 , the greatest number of instances of hybridization have been documented in this group (Table S1). Additionally, many plethodontids of the genus Plethodon in the Eastern United States have been described as species on the basis of a threshold genetic distance (e.g. [39][40][41][42]. Coincidentally, nearly half of these species have been observed to hybridize in nature. To test for an effect of this relative taxonomic over splitting, we excluded the Plethodon glutinosus group (Fig. 1B) and repeated analyses under both the narrow and broad definitions of www.nature.com/scientificreports www.nature.com/scientificreports/ hybridization (see Supplement for additional details). Thus, we produced datasets including either only members of family Plethodontidae, or the nine families that exclude plethodontid salamanders. In total, ten analyses were conducted, using five trees for each of the narrow and broad datasets, respectively. Information on the total number of species included in each tree, as well as the number of species in each state (i.e., H or NH) may be found in Supplementary Table S1. Assessment of trait dependent diversification. We applied the HiSSE (Hidden State Speciation & Extinction 61 ) trait-dependent diversification model to test for differences in speciation and extinction rates between hybridizing and non-hybridizing lineages. HiSSE infers speciation and extinction rates for a binary character while allowing for heterogeneity in diversification rate to exist within each character state. In the HiSSE model, hidden states are co-distributed with the trait of interest and account for unsampled traits that may simultaneously contribute to the diversification process. The hidden state need not be associated with any single trait but may instead be associated with a set of traits or suite of traits. Inclusion of the hidden state thus ameliorates the confounding effects of unsampled traits on diversification rate estimations by allowing for greater rate heterogeneity in the tree than in previous SSE models 61 . Thus, we are in essence measuring the impact of contemporary hybridization while controlling for other, correlated traits on diversification rate. Finally, by accounting for increased rate heterogeneity in character-dependent and character-independent (null) models through the inclusion of hidden states, model rejection properties are greatly improved 61 relative to previous SSE models.
We evaluated a total of 14 competing models using HiSSE, half of which represent a model of character-dependent diversification and the remaining half represented models of character independence (Supplementary Tables S2-S9). Models of trait dependence varied in the number of hidden states included and in Figure 1. Inferred speciation rate along the salamander phylogeny for both hybridizing and non-hybridizing species. Branch outlines depict speciation rates corresponding to the inset density plot. Branch interiors depict the probability that a lineage is hybridizable (white) or non-hybridizable (black). The arrow labeled A denotes family Plethodontidae, and the arrow labeled B denotes the Plethodon glutinosus group.
www.nature.com/scientificreports www.nature.com/scientificreports/ the number of free transitions among states. In character-dependent models, our four states where Hybridizing (H) and Non-hybridizing (NH), each being associated with one of the two hidden state (A & B) for four total character states/diversification regimes (H-A, H-B, NH-A, NH-B).
Additionally, to account for the fact that diversification rates may be biased by incomplete sampling of extant diversity within a phylogeny 62-64 , we assume that 20% of extant species of salamanders hybridize in nature. This value was chosen because it approximately equals the mean frequency of hybridization across our datasets (19.75%). To explore the effect this assumption had on our results, we repeated these analyses assuming (1) Fig. S4). Thus, while unable to designate all hybridizing taxa as such in our phylogeny due to their not being represented in the literature, we have explicitly addressed this uncertainty in our analysis. We elaborate upon our choice of sampling fractions in the Supplementary Materials.
To assess whether the inclusion of extremely young species impacted our diversification rate estimates, we ran a single analysis excluding species younger than 1MY, as per Beaulieu & O'Meara 61 using the narrow dataset and all species. This in turn led to the removal of 14 species. Results were qualitatively identical to those obtained including these young species, so all subsequent analyses were conducted including them. To improve the performance of the Maximum Likelihood optimization procedure implemented in HiSSE, we used simulated annealing to first traverse the likelihood surface to identify optimal starting values for subsequent ML-optimization. Rather than reporting the results of individual model fits, we instead take the approach of investigating model-averaged parameter estimates for each sampled character state 65 (see Supplementary Information for further justification). That is, parameter estimates obtained from each fitted model are averaged together such that their contribution to the average is proportional to their relative support (Akaike weights) among the set of candidate models (Supplementary Table S13). This leads the best supported models to have the greatest impact on the final model averaged parameter estimates. Diversification rates are returned for each sampled state respectively, as not all models include hidden states.
To test for significance among diversification rates inferred for each state, we calculated all possible ratios between non-hybridizing and hybridizing species' model-averaged parameter estimates and calculating the proportion of comparisons in which the value for the non-hybridizing lineage is greater than that of hybridizing lineage. Thus, we obtain empirical P-values (reported in Supplementary Table S14) in which a value of 0 means in every comparison, hybridizing lineages were inferred to have rates greater than those of non-hybridizing lineages, and vice versa. This test is extremely conservative and tests the null hypothesis that non-hybridizing species always experience diversification rates greater than hybridizing species.
In summary, four phylogenies (all species, sympatric species, plethodontids, & non-plethodontids) were analyzed using two datasets (narrow & broad definitions of hybridization). Each data/tree-set combination was tested assuming the three aforementioned differences in prevalence of hybridization in nature. In all, a total of 32 rounds of model testing were performed, comparing seven models of trait-independent diversification and seven models of trait-dependent diversification for each sampling fraction. Further information on how these data/tree-sets were produced and analyzed may be found in the Supplementary Materials. Lastly, a description of our test of sensitivity to phylogenetic uncertainty may also be found in the Supplement (Supplementary Figs. S6 & S7; Table S15) www.nature.com/scientificreports www.nature.com/scientificreports/ To complement our HiSSE analyses in a manner that is largely insensitive to the potentially confounding relationship between branch lengths and propensity to hybridize (see discussion), we conducted sister clade comparisons. Specifically, using all comparisons of three or more taxa (i.e. two sister species hybridize, and the sister lineage does not), we used the method of Barraclough, Harvey and Nee 66 to test the hypothesis that hybridizing clades had greater richness than non-hybridizing clades. To assess confidence, 1000 permutations of contrast signs were conducted. This test was repeated for each of the eight datasets described above (Supplementary  Table S16).
To test for possible circularity of causality between diversification rates, species richness and opportunity to hybridize, we quantified the relationships between the three. Specifically, we tested for a relationship between 1) mean diversification rates and the proportion family hybridizing, 2) mean diversification rates and species richness, and 3) proportion family hybridizing and species richness ( Supplementary Fig. S8). All analyses were conducted at the family level, and mean diversification rates were obtained using lineage-specific model-averaged diversification (speciation, extinction, net-diversification) rates. Simple linear regressions were conducted in R v3.6.1 67 .

Results
Across both datasets (narrow and broad) assuming 20% of species hybridize, three of four analyzed phylogenies consistently found that hybridizing species experience increased speciation rates, decreased extinction rates, and therefore increased net-diversification rates (all significant: Figs. 1-3; Tables 1, 2, Supplementary Tables S6 & S7). Net-diversification of hybridizing lineages in these trees were on average 4X greater than that of non-hybridizing lineages. Our results were insensitive to phylogenetic uncertainty (Supplementary Fig. S7; Table S15).
In all cases of inferred trait dependent diversification, a HiSSE model was the best supported (Supplementary  Tables S2-S9). Of these HiSSE models, two groups emerged: one in which hybridizing species did not harbor a second hidden state (All Species -Narrow, Plethodontids -Narrow & Broad), and one in which both hybridizing and non-hybridizing taxa had two hidden states ( Table 2). The former of these (without a second hidden state for hybridizing taxa) is interpretable as meaning that there is less diversification rate heterogeneity experienced by hybridizing taxa than by non-hybridizing taxa. Interestingly, non-hybridizing taxa were sometimes inferred to have slightly negative net-diversification rates ( Table 2). In the case of the best-fit model using all taxa assuming 20% of species hybridize, this leads to an expected waiting time of 128 million years before the next net-loss in diversity. These negative rates are not persistent however; examination of ancestral state reconstructions indicates that the hidden state responsible (NH-B) for these rates is distributed primarily along the tips (Supplementary Figs S9-S12).
Sister clade contrasts broadly supported results of the HiSSE analysis. That is, hybridizing clades were found to have significantly greater species richness than non-hybridizing clades for all datasets except for plethodontid salamanders. Details on significance of these tests may be found in Supplemental Table S16.

Discussion
Here, we show strong support that contemporary hybridization is correlated with elevated diversification rates in the order Caudata. Net-diversification of hybridizing species tends to be significantly greater than that of non-hybridizing species, driven primarily by a coincident increase in speciation and decrease in extinction rates (Figs. 2, 3; Tables 1, 2). The accelerated diversification of hybridizing salamanders appears ephemeral, however; rate differences rapidly become less pronounced deeper in the tree due to turnover of the hidden states (Fig. 1,  Supplementary Figs. S9-S12). Possible mechanisms leading to this result include frequent range expansions and contractions (i.e. 68 ) that have been documented in salamanders (e.g. [69][70][71] ) and the process of reinforcement which has long been recognized to contribute to the diversification process 14,72-74 . We outline the potential contribution of each below.
Salamanders often exhibit substantial genetic differentiation at small geographic distances (e.g., 200 m 75 ) owing to limited dispersal abilities and low rates of gene flow 76 , thus leading to an abundance of opportunities to evolve in allopatry. Additionally, terrestrial species such as the Plethodontid salamanders of the southeastern United States experience elevational range expansions and contractions associated with climatic change 71 . Perhaps this combination of the primarily sessile nature of many salamander species and frequent repeated secondary contact leads to hybridization occurring regularly across evolutionary timescales. Under these scenarios, hybridization may then play a creative evolutionary role in the diversification process similar to that observed in haplochromine cichlids 24,25,77 . Allopatric speciation of haplochromine cichlids has occurred in lakes that frequently have dried, split, and reformed, whereas sympatric speciation has occurred within lakes in which lineages exhibited extreme habitat specificity and have been reproductive isolated at fine spatial scales 78 . Under these circumstances, hybridization may have afforded genetic rescue from the consequences of small population size by providing increased standing genetic variation and thus expedited adaptation to novel stressors, as in Lake Victorian cichlids post-colonization 25 .
Such hypotheses of repeated contractions, expansions, and secondary contact of salamander populations have been well supported, across both North [79][80][81] and South America 50,82 , as well as in Europe [83][84][85] . Often associated with glaciation/deglaciation or orogeny of mountain ranges, geological events may act as species pumps for salamanders (i.e. 86 ). However, while periodic geographic range expansion and contraction may initiate speciation, the reproductive isolation that evolves may be incomplete, predisposing the young species/diverging lineages to hybridization (i.e. 23 ). As a result, hybridization may commonly occur in salamanders during periods of climatic fluctuations.  Table 2. Best-fit models assuming 20% of species hybridize. Included are models that received >5% Akaike weights for their respective analyses. For each dataset, the best fit model is bold. Maximum-likelihood parameter estimates for speciation, extinction, and net diversification are reported. Non-hybridizing is abbreviated as NH, Hybridizing as H; A and B indicate the two hidden states.

Narrow
www.nature.com/scientificreports www.nature.com/scientificreports/ There is now clear evidence for latitudinal and elevational range shifts mediated by climate change 87-91 and a consequent increase in frequency of hybridization among previously isolated taxa [92][93][94][95] . An informed understanding of the influence of hybridization on macroevolutionary diversification may thus provide invaluable context for contemporary processes. This possibility of climate-change mediated hybridization has already been demonstrated in plethodontid salamanders (P. shermani & P. teyahalee 71 ), as well as in ecologically divergent subspecies of salamandrid salamanders (S. salamandra 96 ). Thus, it seems likely that salamander species worldwide, particularly those found at high elevations due to their more limited potential geographic distributions, may experience a heightened frequency of hybridization as climate change advances. While generalizations regarding the outcome of hybridization should be made with caution 30 , our study indicates that perhaps speciation reversal [97][98][99] need not be the expectation. Rather, our study implies that hybridization may facilitate adaptation to novel conditions under climate change, leading to diversification of new salamander lineages.
Here, we show the novel result of a strong correlation of contemporary hybridization with elevated speciation and net diversification at a large taxonomic scale. However, reinforcement, defined as the strengthening of prezygotic reproductive isolation in sympatry 14 , is intrinsically intertwined with hybridization. Reinforcement has been documented both experimentally 72 and observationally 74 to accelerate the initiation and/or completion of the speciation process 73,100,101 . For instance, reinforcement is likely to play an important role in the speciation process due to strong interspecific sexual selection and mate choice in plethodontids 102 . Indeed, patterns of sexual isolation among populations of Plethodon jordani and P. teyahalee match expectations of reinforcement 14 , with sexual selection being stronger in sympatry than in allopatry 103 . Although we cannot currently quantify the contribution of reinforcement to diversification rate differences using our data, we urge further research measuring the degree of association between contemporaneous hybridization and reinforcement among taxa. Nonetheless, were reinforcement to play a role in the production of the patterns observed in this study, the very occurrence/ process of hybridization would be the ultimate driver (i.e., cannot have reinforcement without hybridization). Under such a scenario, our study design is well-suited to identify such a signal.
Although a generative role of hybridization is robustly supported across three of our four datasets, evidence for such a role outside of the Plethodontidae is more limited ( Table 1, Supplementary Tables S2-S14). We find two possible explanations for this finding. Firstly, the positive association between diversification rates and hybridization may be unique to Plethodontid salamanders. However, family Plethodontidae is the largest extant family of salamanders, comprising approximately 2/3rds of the present diversity (471 of 716 species: amphibiaweb.org). Thus, our observation of hybridization facilitating the diversification process applies to the majority of salamanders and implies that, at a minimum, contemporaneous hybridization does not impede the diversification process of extant salamanders.
Secondly, it is highly probable that our analysis of non-plethodontid salamanders is lacking in power. SSE models have long been known to lose much of their power when dealing with small number of OTUs (trees < 300 taxa 104 ). For example, for trees of 300 species, BiSSE attains a power of at most 50%, with power dropping below 15% of trees of 100 taxa 104 . Our phylogeny of non-plethodontids includes only 167 species; that we detected a positive relationship between hybridization and speciation rates using our narrow (most conservative) dataset despite such reduced power is a testament to the strength of the signal in our data. Whereas our larger datasets [complete (469 spp), sympatric (368 spp) and plethodontids (306 spp)] have greater power, our lack of detection of a relationship between hybridization and diversification rates in non-plethodontids using our broad definition of hybridization is perhaps unsurprising, given the low power of the analysis (also see Supplementary Materials for an elaboration of power). Although the power of HiSSE under such scenarios has not been specifically established, accuracy of parameter estimation does decay with decreasing tree size 61 . Consequently, we cautiously interpret the results of the analysis of non-plethodontids (163 spp). Interestingly, sister clade comparisons consistently supported a positive relationship between hybridization and species richness in non-plethodontids, despite not supporting such a relationship in plethodontids (Supplementary Table S16). These results are insensitive to branch-lengths, thereby ameliorating potential concerns related to the relationship between hybridization, branch-lengths, and diversification rates 23 .
Importantly, parameter estimates are largely reasonable. For instance, the greatest speciation rate inferred by any analysis (Plethodontids assuming 30% of taxa hybridize using the full tree: Supplementary Table S13), of 0.159 species/million years (MY) can be interpreted as a waiting time, such that on average, hybridizing species speciate every 6.29 MY. Extinction rates appear less reliably estimated however; some estimates functionally equal zero, leading to the large percent decrease in extinction rates observed for hybridizing relative to non-hybridizing species. In some cases, extinction rates exceed speciation rates in non-hybridizing taxa, leading to negative net-diversification rates. That being said, averaged extinction rates across all analyses for hybridizing and non-hybridizing taxa led to more reasonable waiting times of 127.9 and 54.2 MY respectively. Further, recent studies have documented even more negative net-diversification rates than inferred herein 105 . Taken together, it appears that extinction plays an important role in the diversification of salamanders, leading to a reduction in net-diversification rates towards the present relative to hybridizing species.
An important question regarding the interpretation of our results is the relationship between lineage diversification rates, species richness, and opportunity to hybridize. Because the relationship between lineage diversification rate and opportunity to hybridize are not necessarily independent, rapidly diversifying lineages may simply have greater opportunity to hybridize due to increased diversification rates. Although a legitimate concern, we did not find evidence that the increased diversification rates we observe are due to increased family-level species richness leading to increased opportunity to hybridize ( Supplementary Fig. S8). Further, it is unlikely that non-random taxonomic sampling has biased our results, as there is no relationship between clade specific sampling fraction and frequency of hybridization ( Supplementary Fig. S5).
Although the ability of methods to accurately infer extinction rates has been debated recently 63,106,107 , we emphasize that our results are robust to this concern. Our central result, that hybridizing lineages experience www.nature.com/scientificreports www.nature.com/scientificreports/ increased net diversification, is driven by both increased speciation rates and decreased extinction rates. Further, in nearly all cases, the magnitude of increase of speciation rate is greater than that of the decrease in extinction. Thus, our results are likely robust even to inaccuracies in the estimation of extinction rate.
An important distinction between our study and most previous studies investigating the influence that hybridization exerts on the diversification process is that of the time-scale at which hybridization is being assessed. Following Seehausen's 6 landmark paper "Hybridization and Adaptive Radiation," tests and discussion of his hypothesis, that ancient, widespread hybridization facilitates adaptive radiation became abundant in the literature (e.g. 18,22,23,108,109 ). Whereas much of the subsequent studies focused on ancient hybridization, our study instead focuses on the effects of contemporary hybridization.
In a pertinent study, Wiens et al. 23 , tested the hybrid swarm hypothesis in the Plethodon glutinosus group (indicated in Fig. 1B) using two nuclear and two mitochondrial genes. They did not recover strongly supported evidence of genealogical discordance at the base of this group; these results were interpreted as not being supportive of Seehausen's hypothesis. Further, they identified a positive relationship between age of species and reproductive isolation. They argue that the observed relationship between diversification rate and hybridization in this group was a consequence of this relationship. Although a legitimate concern, we argue that this hybridization is likely to still have biologically relevant consequences on diversification rates. Specifically, hybridization may either 1) facilitate the divergence of these young species i.e. through reinforcement/strengthening of prezygotic isolation, or 2) erode their divergence leading to species collapse. Whereas the former hypothesis predicts increased speciation rates, the latter predicts increased extinction rates. We find strong, consistent evidence in favor of the former.
We explicitly tested the hypothesis that contemporary hybridization plays a creative role in the diversification process in the broadest taxonomic and temporal scale study to date, and our observations strongly supported the predictions of this hypothesis. Specifically, hybridization was found to be correlated with both increased speciation rates and decreased extinction rates, resulting in increased net diversification rates relative to non-hybridizing lineages. Although other factors certainly contribute to the observed diversification dynamics, we have shown that hybridization plays a significant role, while accounting for hidden, correlated states in our analysis. Nearly all studies of hybridization have focused on individual case studies in which hybridization results in species collapse 98 or promotes diversification in a single species group 12,13,22,27 . Such studies are necessarily limited in the extent to which their results may be generalized 30 , particularly because results were equivocal across studies. Consequently, we advocate that our approach can be applied at broad taxonomic and evolutionary timescales to facilitate robust tests of the role of hybridization in the lineage diversification process. We anticipate our results are broadly generalizable to animal groups in which homoploid hybridization occurs because only 17 species of salamanders are known to be polyploid 110 , and our dataset includes only seven (Ambystoma mexicanum, A. barbouri, A. jeffersonium, A. laterale, A. texanum, A. tigrinum, and Lissotriton vulgaris) hybridizing polyploid taxa (none of which are plethodontids). Our study adds to the growing evidence that hybridization may fuel rapid diversification (e.g. 52 ) and is a compliment to speciation genomics studies characterizing the genomic basis of this process (e.g. 17,111 ). Herein we have shown that hybridization may act as a generative force across a phylogenetic order, and additional studies at such macroevolutionary scales are needed to determine if this pattern holds more generally across the tree of life.