Changing Ecological Opportunities Facilitated the Explosive Diversification of New Caledonian Oxera (Lamiaceae)

Abstract Phylogenies recurrently demonstrate that oceanic island systems have been home to rapid clade diversification and adaptive radiations. The existence of adaptive radiations posits a central role of natural selection causing ecological divergence and speciation, and some plant radiations have been highlighted as paradigmatic examples of such radiations. However, neutral processes may also drive speciation during clade radiations, with ecological divergence occurring following speciation. Here, we document an exceptionally rapid and unique radiation of Lamiaceae within the New Caledonian biodiversity hotspot. Specifically, we investigated various biological, ecological, and geographical drivers of species diversification within the genus Oxera. We found that Oxera underwent an initial process of rapid cladogenesis likely triggered by a dramatic period of aridity during the early Pliocene. This early diversification of Oxera was associated with an important phase of ecological diversification triggered by significant shifts of pollination syndromes, dispersal modes, and life forms. Finally, recent diversification of Oxera appears to have been further driven by the interplay of allopatry and habitat shifts likely related to climatic oscillations. This suggests that Oxera could be regarded as an adaptive radiation at an early evolutionary stage that has been obscured by more recent joint habitat diversification and neutral geographical processes. Diversification within Oxera has perhaps been triggered by varied ecological and biological drivers acting in a leapfrog pattern, but geographic processes may have been an equally important driver. We suspect that strictly adaptive radiations may be rare in plants and that most events of rapid clade diversification may have involved a mixture of geographical and ecological divergence.

Oceanic islands are widely regarded as laboratories of evolution mainly owing to their isolated past biogeographical history (Losos and Ricklefs 2009). Their geographical isolation limits dispersal events from outside, resulting in ecological niches being filled by species diversification rather than colonization. Consequently, islands are expected to have sheltered many adaptive radiations, where single lineages have diversified rapidly into several distinct niches, resulting in disproportionately high phenotypic diversity (Schluter 2000) when compared with their continental relatives (Silvertown et al. 2005). This phenomenon can be particularly apparent on large and high-elevation islands that exhibit diverse and sharp environmental gradients (Whittaker et al. 2008). Although adaptive radiations may not necessarily result from an acceleration of speciation rates (see Givnish 2015 for a discussion), all adaptive radiations are spurred by ecological opportunities; that is, the release of vacant niches due to unused resources, species extinctions or the acquisition of new traits (Losos 2010). The existence of adaptive radiations thus posits a central role of natural selection causing ecological divergence and speciation, and some plant radiations have been highlighted repeatedly as paradigmatic examples of this scenario, such as the silversword alliance and lobeliads in Hawaii (Carlquist et al. 2003;Givnish et al. 2009), the Macaronesian Aeonium and Echium (Jorgensen and Olesen 2001), and Veronica (Hebe) in New Zealand (Wagstaff and Garnock-Jones 1998).
However, adaptive radiations sensu stricto appear to be quite uncommon in the plant kingdom, with the aforementioned well-documented cases probably representing rare exceptions. In fact, many radiations likely occur, at least in part, in a non-adaptive way (Givnish 1997, Rundell andPrice 2009). Under such a scenario, speciation is not primarily driven by ecological divergence, but mainly by neutral divergence in allopatry. Ecological differentiation may occur following nonecological speciation, through random trait divergence, adaptation to different climatic regimes or habitats, or due to trait divergence favoring local coexistence if secondary sympatry occurs (e.g., Tobias et al. 2014). While molecular phylogenies will often provide little power to discriminate between different scenarios, it remains possible to assess the relative importance of ecological differentiation and geographic divergence in rapid species diversifications. Considering that theoretical evidence suggests that most rapid evolutionary radiations may have occurred through a mixture of ecological speciation and neutral divergence in allopatry (Aguilée et al. 2012(Aguilée et al. , 2011, we can legitimately wonder The New Caledonian biodiversity hotspot is a remote archipelago in the southwest Pacific, composed of an old-, large-, and high-elevation main island (= Grande Terre, ca. 37 Ma, ca. 16,600 km 2 , 1600 m max. elevation) exhibiting exceptional levels of species richness and endemism (ca. 3400 species, 75% endemic; Morat et al. 2012;Munzinger et al. 2016). The singular geological and climatic New Caledonian history led to the implementation of complex and abrupt environmental gradients, resulting in a mosaic of highly distinct habitats. The archipelago topography is especially complex on Grande Terre, with many relatively steep valleys and high summits (Bonvallot 2012) presenting both physical barriers to dispersal and strong selection gradients. The archipelago also has a suite of diverse bedrock types, in particular volcano-sedimentary, metamorphic and ultramafic (metal-rich soils with chemical and physical properties that constrain plant growth). These latter bedrock types, which act as strong environmental filters and exhibit a patchy distribution, have probably played a key role in plant speciation and flora evolution (Pillon et al. 2010). In addition, Pliocene and Pleistocene climatic fluctuations considerably affected the dynamic of New Caledonian biotas, leading to the origination of new habitats such as the unique shrubby sclerophyllous vegetation (i.e., maquis; Jaffré 1980), during periods of forest contraction (Hope and Pask 1998;Stevenson 2004;Stevenson and Hope 2005). This likely triggered speciation in palms (Pintaud et al. 2001), and has contributed to the persistence of old angiosperm lineages in forest refugia (Pouteau et al. 2015). All the aforementioned characteristics have contributed to the tremendous plant diversity of New Caledonia (Jaffré 1993), and have also driven numerous plant radiations (e.g., Psychotria; Pycnandra; Munzinger et al. 2016). Though plant radiations are increasingly highlighted in the archipelago through molecular investigations, few studies have focused on the relative effect of drivers implied in their in situ diversification (e.g., Barrabé et al. 2014;Pillon et al. 2014;Paun et al. 2016).
The old age, geographic isolation, and geologic and topographic complexity of New Caledonia suggest that adaptive radiations may be common among many genera endemic to New Caledonia. However, few of these genera exhibit the joint signatures of morphological and ecological divergence expected from adaptive radiations, and rather seem to be relictual lineages (Pillon et al. 2017). The woody genus Diospyros (Ebenaceae) was the first purported clear case of an adaptive radiation of a plant group within New Caledonia due to its wide ecological diversity (Paun et al. 2016), although the radiation was not related to obvious morphological and/or physiological differences. The tree genus Geissois (Cunoniaceae) was described as a cryptic adaptive radiation, in which species exhibit notable differences in leaf element composition likely linked to the occupation of varied soils ), but does not fully satisfy the criteria of an adaptive radiation as outlined by Givnish (2015), as no physiological adaptation has yet been highlighted. It thus seems that New Caledonia has harbored few, if any, adaptive radiations in a strict sense (sensu Givnish 2015 andSchluter 2000). This has been explained as a consequence of the relatively old age of the island, the prevalence of woody species, the paucity of potential pollinator species, and/or a reduced number of ecological opportunities as compared with other island systems (Pillon et al. 2017).
In this work, we report on a molecular systematic investigation of the genus Oxera (Lamiaceae), which has been hypothesized to have diversified in New Caledonia through adaptive radiation ). The genus is composed of 33 species endemic to New Caledonia (Gâteblé unpublished data, available at http://endemia.nc/flore/fiche588), and was recently enlarged to include five other Australasian and Pacific species (Barrabé et al. 2015). It was demonstrated that the entire New Caledonian clade originated from a single colonization and constitutes the tenth largest plant radiation in the archipelago, which is puzzling given that Lamiaceae are an under-represented family within the flora of New Caledonia (Pillon et al. 2010). Distinct subclades were partially circumscribed within Oxera and can arguably be considered as several independent micro-radiations within the genus (Barrabé et al. 2015). Oxera exhibits strongly divergent morphology in terms of life form, flowers and fruits, and occupies a vast diversity of distinct habitats. In relation to this remarkable morphological diversity, we hypothesize that shifts of functional, reproductive (pollination and dispersal), and environmental niches may have played a key role during the diversification of Oxera. This role can be two-fold: 1) niche shifts may have triggered species differentiation and speciation, and 2) certain niche types may have accelerated the rate of genetic divergence and speciation.
The dazzling morphological and ecological diversity of New Caledonian Oxera could reflect adaptations to different ecological factors (see Givnish 2010 for a discussion of life form adaptations). First, the great array of life forms observed in Oxera could be associated with different environment preferences, as observed in the silversword alliance (Carpenter et al. 2003), and this may have created conditions that facilitated parapatric isolation and speciation. Robust lianas mainly grow in closed rainforests with leaves usually reaching the canopy, while slender lianas are mostly encountered in open sclerophyllous vegetation and forest edges; monocaulous trees (i.e., single-stemmed) are mainly confined to rainforest understories; and finally, shrub species are often able to establish in both closed and more open vegetation (Gâteblé, pers. obs.; see Fig. 1a and b). Such contrasting life forms clearly reflect functional strategies regarding growth under different light conditions (Givnish 1995;Santiago and Wright 2007;Selaya and Anten 2008). Second, most species of Oxera are confined to narrow environmental conditions within New Caledonia. They occur either in rainforests, dry forests or maquis, or occur at different elevations and/or on different bedrock types (limestone, ultramafic, and volcano-sedimentary rocks). Although some species of Oxera are clearly specialized toward particular habitats or have very limited geographic ranges in narrow valleys or on remote summits, others are encountered throughout the archipelago due to wider environmental tolerances. These contrasting habitat and geographic patterns may be hypothesized as the signature of parapatric or allopatric speciation processes, and provide an interesting basis for testing whether habitat or geography (or both) have been major drivers of clade diversification.
Third, the floral diversity of Oxera was previously linked to distinct flower pollinators (de Kok 1998) based upon field observations of animal visits (Table 1).
These pollinator preferences could potentially be drivers of evolution of reproductive isolation within Oxera. Honey-eaters (Meliphagidae) and white-eyes (Zosteropidae) have been observed visiting Oxera species displaying showy and curved flowers, profusely producing nectar, and bearing arched anthers typically extending out of the corolla (Table 1). Long-tubular white and/or green bilabiate flowers, sometimes releasing a strong fragrance, have been observed to be visited by butterflies and moths (Pieridae and Sphingidae; Table 1 ;de Kok 1997;Haxaire and Salesne 2016). Despite the lack of recorded visitor observation, the flower attributes of the Australasian O. splendida, which bears large, flared, white flowers that open nocturnally and release a strong fragrance, are very similar to that of Fagraea, renowned for its bat-pollination (Momose et al. 1998).
Fourth, the large variety of fruit morphology exhibited by Oxera suggests different dispersal agents, although these mechanisms have not been fully identified in the field (Table 1). We assume that the different behaviors and flight abilities of animal dispersers would have affected the spatial scale of gene flow, and thus produced different rates of genetic divergence between Oxera lineages, due to varied dispersal abilities and habitat selection of dispersers (Givnish 2010, Theim et al. 2014). The largest white fruits encountered in Australasia are consumed by bats and cassowaries (Richards 1990;Hutson et al. 1992;Cooper and Cooper 2004). Other large colored fruits (ca. >2 cm) could be exclusively eaten by imperial pigeons (Ducula, Columbidae), which is the only bird able to swallow and effectively disperse such large fruits (Meehan et al. 2002;Carpenter et al. 2003;Cibois et al. 2017;P. Bachy pers. comm. for O. palmatinervia). It is noteworthy that the sedentary behavior of pigeons confines them mostly to forest understories (Colombo 2008;Wotton and Kelly 2012). Smaller fruits (ca. <2 cm) are consumed by more diminutive birds such as Sturnidae (Table 1; Tassin et al. 2010), which probably move over short distances, but are able to shift between different vegetation types (Barré et al. 2006). Finally, large grey or brown fruits having armed appendages are considered to be consumed by large geckos (Rhacodactylus, Mniarogekko, and Correlophus genera) because they are not discouraged by the fruits' repelling structures (Whitaker and Bauer, pers. comm.). Large geckos also have limited spatial movements (Clark et al. 2009), which may have constrained gene flow in Oxera.
Given the large diversity of morphological traits, ecological attributes, and geographical ranges exhibited by Oxera, we suspect that species diversification within subclades of this genus has been driven by several ecological and biological factors acting at different time periods, which is in a leapfrog pattern (as originally defined by Chase and Palmer 1997). Available information on species coexistence within similar regions or habitats may also help test whether sympatric speciation has occurred in relation to different pollination and dispersal niches. Indeed, such a scenario suggests that sister species with distinct biotic niches should co-occur either in the same region or habitat. Sister species divergence based only on ecological preference or geographic range would otherwise suggest that speciation mainly occurred through parapatric or allopatric speciation, respectively. We thus posit that phylogenetic patterns of species ranges and niches bear the signature of both selective and neutral forces as engines of clade diversification. Here, we report the first detailed phylogenetic reconstruction of this exceptional radiation of Lamiaceae, the genus Oxera, within the New Caledonian biodiversity hotspot, and attempt to identify the major factors that have shaped its current species diversity. We used a combination of phylogenetics, molecular dating, and comparative analyses of species diversification and niche evolution in order to 1) reconstruct the tempo and rate of species diversification in Oxera, 2) extricate and identify the different drivers of its diversification, and 3) appraise their respective effects. This study specifically asks a pivotal question: Is Oxera the single clear adaptive radiation in New Caledonia?

Species Sampling and Molecular Data Sets
In order to both establish species relationships and to estimate divergence times within the genus Oxera, the focal ingroup was composed of the 38 recognized species according to the most recent systematic investigations in New Caledonia (Gâteblé unpublished data, available at http://endemia.nc/flore/fiche588), and taxonomical revisions by de Mabberley (1999a, 1999b) for Australasia and Pacific. We added to the data set the sister lineage Hosea, as well as 10 species belonging to the closely related Clerodendrum clade identified in Yuan et al. (2010). Three more divergent species of Lamiaceae were used as external outgroups (i.e., Ajuga chamaepitys, Teucrium pyrenaicum, Rotheca sp. Wen 9487; Supplementary Appendix S1a; available on Dryad at https://doi.org/10.5061/dryad.mm1ft40). A total of 51 species were used in this sampling, which we refer to as the Oxera data set, including 22 species with previously published sequences (downloaded from GenBank), and 29 species with new sequences generated for this study. This molecular data set consisted of the 12 DNA loci sequenced in Barrabé et al. (2015). This allowed us to improve 1) the phylogenetic resolution and 2) branch length estimates relative to Barrabé et al. (2015). All new DNA sequences were generated using the DNA extraction, amplification, and sequencing  Barrabé et al. (2015). All accessions used in this study are listed in the Supplementary Appendix S1a available on Dryad.
For divergence time estimation, we also used two additional alignments composed of four DNA plastid regions (matK, rps16, trnL-F, and trnL-trnF; Supplementary Appendix S2a available on Dryad) from two data sets, with taxon sampling that spanned the Lamiales and Lamiaceae. The Lamiales data set included: 1) a subset of the Oxera data set (i.e., five species; Supplementary Appendix S1c available on Dryad); and 2) 177 representatives of various lineages within Lamiales, whose sequences were downloaded from GenBank (Supplementary Appendix S1c available on Dryad). The Lamiaceae data set was thus composed of 1) the aforementioned Oxera subset and 2) 59 taxa representing all major Lamiaceae lineages (downloaded likewise from GenBank; Supplementary Appendix S1b available on Dryad). This large-scale sampling allowed us to incorporate Lamiales and Lamiaceae fossils as calibration points (see below).
All DNA alignment matrices used in this study are available on Dryad.

Phylogenetic Inference
For the three molecular data sets, phylogenetic inferences were conducted using Bayesian Markov Chain Monte Carlo (MCMC) as implemented in MrBayes v3.2.1 (Ronquist et al. 2012) based on single-locus and concatenated data sets. Best-fit nucleotide sequence evolution models were identified using jModelTest (Posada 2008) based on the Akaike criterion (more details on settings are provided in Supplementary Appendix S2a available on Dryad). The combined data sets were partitioned to allow each locus to possess specific model parameters (Nylander et al. 2004). The methodological approach of the Bayesian MCMC analyses followed that described in Barrabé et al. (2015) and the details of parameter settings are provided in Supplementary Appendix S2b available on Dryad. All analyses were run on the Cyber Infrastructure for Phylogenetic Research cluster (CIPRES; http://www.phylo.org/; last accessed 18 October 2017). The post-burn-in trees resulting from the MCMC stationary phase were used to construct a majority-rule consensus tree and calculate Bayesian posterior probabilities (PPs); clades were considered well-supported when PP values were ≥0.95). The three resulting consensus trees are available in TreeBASE (ID 23353) and Dryad.
In order to check for phylogenetic incongruence between nuclear and plastid loci, we also ran two phylogenetic reconstructions separately for all six concatenated nuclear regions and all six concatenated plastid regions, using the same MrBayes parameter settings and models of nucleotide evolution as above. These two analyses were only applied to the Oxera data set. We then assessed whether there was any supported (PP ≥0.95) incongruence between two trees.

Fossil Dating and Divergence Time Estimation
We performed a two-step approach to estimate the temporal evolution of Oxera. This allowed us to test the impact of different dating calibrations and taxonomic sampling on divergence time estimates. In the first step, we performed a fossil calibration on the Lamiales and Lamiaceae data sets using fossils described in a recent Lamiaceae molecular dating study (Yao et al. 2016). We incorporated three Lamiaceae (Melissa, Ocimum, and Stachys) and three Lamiales fossils (Bignoniaceae, Acanthaceae, and Oleaceae) for dating the phylogeny inferred with the Lamiales data set, and two Lamiaceae fossils (Melissa and Ocimum) for dating the phylogeny inferred with the Lamiaceae data set (for further details on the design and parameters of these calibrations, see Supplementary Appendix S2c, available on Dryad and explanations provided inYao et al. 2016). For root calibrations we used divergence times estimated in Magallón et al. (2015): 1) the divergence between Lamiales and Gentianales-Solanales for the Lamiales data set and 2) the divergence between Plocospermataceae and other Lamiales for the Lamiaceae data set. In the second step, the divergence time estimate between Rotheca (sp. Wen 9487) and other Lamiaceae recovered during the first step (i.e., the dating of the Lamiaceae data set) was used as a secondary root calibration point (see Results section) for the Oxera data set.
For the three molecular data sets, molecular divergence times were estimated using the Bayesian MCMC approach implemented in BEAST v.1.8.0 (Drummond and Rambaut 2007). DNA regions were combined and partitions were set as in the above Bayesian MCMC analyses (Supplementary Appendix S2c available on Dryad). An uncorrelated relaxed molecular clock model was selected following a lognormal distribution, and the Birth-Death process implemented for the tree prior. Other details on analyses, calibration and parameter settings are provided in Supplementary Appendix S2c available on Dryad. The post burn-in trees were summarized, and Bayesian PPs, median height (= age estimate), and the 95% highest posterior density heights interval of each node (95% HPD) assessed, using a Maximum Clade Credibility target tree (named as MCCT tree for the Oxera data set) in Treeannotator v.1.8.0 (Drummond and Rambaut 2007). For the Oxera data set a subset of 100 dated trees (RD trees hereafter) were also randomly sampled from the stationary phase to integrate phylogenetic uncertainty in some of the following analyses (see below). All divergence time analyses were conducted using CIPRES. The three resulting Maximum Clade Credibility trees are available in TreeBASE (ID 23353) and Dryad.

Comparative Analyses of Species Diversification
To obtain an overall view of the Oxera diversification, we first generated a lineage through time diagram, based on a pruned version of the MCCT tree including only Oxera species and a single sample per species (pMCCT tree), and also on the RD trees (pruned as in the pMCCT tree). We also calculated net diversification rates, following the conservative approaches implemented in Pillon (2012) and Barrabé et al. (2014), to allow comparisons with other oceanic insular plant radiations, using median crown ages and their 95% HPD estimated from the molecular dating, and two extreme values for relative extinction (null and equal to 0.9). Net diversification rates were likewise estimated per unit of area and log(area).
To assess more precisely the rate and tempo of species diversification in Oxera through time, we carried out diversification analyses with maximum likelihood (ML) and Bayesian modeling approaches. The ML analyses were conducted in the TreePar R package (Stadler 2011), using an optimization algorithm to infer past rates of speciation and extinction and their temporal variation. We fitted various diversification models: a pure birth model, a constant birth-death model, several more complex models allowing one to several rate shifts through time (up to 10 shifts), and finally several diversity-dependent models. All models were fitted on the pMCCT tree (see Supplementary Appendix S2d, available on Dryad for further details on parameter settings) and on the RD trees to account for phylogenetic uncertainty. The best-fit diversification model was identified using the corrected Akaike criterion, and a likelihood-ratio test was used for comparing nested pure birth, birth-death, and rate shift models.
The Bayesian diversification analyses were computed using BAMM (Rabosky et al. 2014). This allows modelling complex dynamics of speciation and extinction on phylogenetic trees using a reversible jump MCMC algorithm. We pruned the MCCT tree in the same manner as the pMCCT tree, but retained species belonging to the Clerodendrum clade to assess early variation in Oxera diversification rates. The sampling fraction was completed with an estimate of species richness, especially for the Clerodendrum clade (retrieved from the Kew Checklist website: http://wcsp.science.kew.org; last accessed 18 October 2017). As we expected two rate shifts to have taken place during the evolution of Oxera and its relatives (considering the large sizes of the genera Clerodendrum and Oxera as compared with the other small genera), the prior on the number of diversification shifts was set to 2. All other priors were estimated in the BAMMtools R package (Rabosky et al. 2014), with the default settings used for all other parameters. The MCMC was run for 50×10 6 generations and sampled every 10,000 generations. The analysis was conducted with three independent Markov chains. The first 10% generations of each run were discarded after checking for chain convergence and adequate MCMC sampling. The remaining generations were summarized to generate the 95% credibility interval for rate shift configurations, the marginal shift probability tree (where branch lengths are proportional to the probability that a shift occurred on a given branch), the best rate shift configuration with the highest maximum a posteriori probability, the phylorate plot (showing variation of mean diversification rates with a colored gradient) and finally the rates-through time plots (only built for the Oxera focal group).

Biological, Ecological, and Geographic Range Data sets
To identify and disentangle suspected drivers of species diversification in Oxera, we built a data set by compiling data on five major life traits: 1) flower morphology, 2) fruit morphology, 3) life forms, 4) geographical occurrences, and 5) environmental preferences. The flower and fruit morphology (i.e., organ dimensions, colors, shapes, and textures) were compiled from herbarium and fieldwork observations on organ dimensions, colors, shapes, and textures. Those traits were coded as 13 discrete and three continuous characters for flower morphology, and three discrete and three continuous traits for fruits (Supplementary Appendix S3a, d, and e available on Dryad). The life form data set was composed of a single discrete trait corresponding to the four life forms usually ascribed to Oxera, namely monocaulous, robust liana, slender liana, and shrub/small tree (Supplementary Appendix S3a and b available on Dryad). Note that no functional or physiological traits related to light requirements or nutrient use strategy were available for Oxera.
The geographical distributions of New Caledonian species were retrieved from the databases of the herbaria of Nouméa (NOU, VIROT), the Museum national d'Histoire naturelle of Paris (P, SONNERAT), the University of Zurich (Z), and additional field observations. The geographical coordinates of herbarium specimen records were then databased, error-corrected for distribution, and finally incorporated into a Geographical Information System (GIS) under QGIS 2.8.1. For Australasian and Pacific species distributions were assessed from de Kok and Mabberley (1999b) and manually added into the same GIS layer.
The environmental data set was composed of three variables that best explain plant distributions across New Caledonia (Jaffré 1993;Morat 1993; Supplementary Appendix S3a and c available on Dryad): geological substrates, vegetation types, and elevation. Species' rock types were coded as occurring on ultramafic rocks, on volcano-sedimentary and metamorphic rocks, or ubiquists. These geological attributes were determined by coupling geographical occurrences with a layer of geological substrate provided by the Direction de l'Industrie, des Mines et de l'Energie de la Nouvelle-Calédonie (New Caledonia). Species' vegetation types were coded as occurring in wet forests, in sclerophyllous vegetation types (i.e., dry forests and/or maquis vegetation), and ubiquists, based on herbarium and field observations (Gâteblé, pers. obs.). Species' elevation ranges were extracted from the GIS layer generated by a digital elevation model (resolution of 10 m) provided by the Direction des Infrastructures de la Topographie et des Transports Terrestres (New Caledonia), and 468 SYSTEMATIC BIOLOGY VOL. 68 then averaged across geographical occurrences of each species.
To identify ecological syndromes among Oxera species, we performed multivariate analyses on the environmental, flower, and fruit data sets. As we used both continuous and discrete traits, we conducted the PCA analyses by using the Hill and Smith (1976) principal component method under the ade4 R package (Chessel et al. 2004). This analysis allowed clustering species among pollinator syndromes, dispersal syndromes, and habitat preferences, which we subsequently treated as new discrete traits according to the species clustering (i.e., the niche data sets; Supplementary Appendix S3f available on Dryad). We posit that the ecological syndrome of a species reflects the ecological niche it occupies (Johnson 2010, Pauw 2013. No data on species growth and architecture were available, we therefore, could not evaluate life forms in a similar manner; this trait was treated as a discrete variable without any step of ordination analysis.

Comparative Analyses of Niche Evolution
To assess the mode and tempo of evolution of each niche during the Oxera radiation (i.e., of pollination and dispersal syndromes or habitat preferences), we measured their respective amount of phylogenetic signal by estimating and Pagel statistics (Pagel 1999). Values of and close to 1 indicate a strong phylogenetic signal and gradual trait evolution (i.e., according to a model of Brownian motion), and deviation from this expectation provides insight into temporal patterns of trait evolution. Values of close to 0 depicts punctual evolution, where trait divergence occurs independently from branch lengths in the phylogeny. Values of close to 0 reveal a low phylogenetic signal, where most trait change occurred late in evolutionary history and closely related species thus share very little trait similarity. For the habitat, pollination, and dispersal niches and statistics were estimated using the "fitContinuous" and "phylosig" functions of the geiger and phytools R packages, respectively (Harmon et al. 2008;Revell 2012). Estimations were performed on the RD trees, and by averaging estimated values for each of the first n axes of each PCA explaining 90% of the cumulative variance. For life form, each statistic was simply averaged through the RD trees using the "fitdiscrete" function (geiger) and the single discrete trait.
To reconstruct the evolutionary history of habitat preference, pollination, dispersal niches, and life form, we conducted ancestral state reconstruction analyses. We first assessed ancestral niche states using the "ace" function in the ape R package (Paradis et al. 2004) applied to the pMCCT tree, and the niche and life form data sets. We fitted three ML evolution models: the "ER" (with equal state transition rates), "SYM" (symmetrical), and "ARD" (with all rates unequal) models, identifying the best-fitting model using comparisons performed on the corrected Akaike criterion, and likelihood ratio test. We then estimated their respective absolute shift number and timing by performing stochastic character mappings (Huelsenbeck et al. 2003), using the "make.simmap" function in phytools. We launched 500 simulations using an estimated prior distribution on the root node, the best-fit ML model recovered in the previous analyses and applied first to the pMCCT tree, and then to the RD trees to obtain most likely envelopes of past shift number for habitat preference, pollination, dispersal niches and life-forms. For 11 of the 100 RD trees, we encountered intractable optimization issues when using the "make.simmap" function. We subsequently removed the 11 trees from these analyses and used only 89 RD trees. Absolute shift timings were extracted using a custom R-script from all simulations and summarized by plotting the median number of state changes in different time bins of a 0.5 Myr time-grid.
To investigate whether a particular niche state (pollination, dispersal, or habitat) or life-form would have affected clade diversification rates of Oxera, we conducted multi-state trait-based analyses that estimated simultaneously trait-dependent speciation, extinction, and transition rates on a phylogeny and character distribution (Fitzjohn et al. 2009). These were performed using the "make.musse" function under the diversitree R package (Fitzjohn et al. 2009), applied to the pMCCT tree, and the niche and life form data sets. We first fitted 17 models from a null model (rates independent of trait states) to a full model (all rate components dependent of trait states; Supplementary Appendix S4 available on Dryad), conserving the best-fitting using the corrected Akaike criterion. We performed subsequent Bayesian analyses through a MCMC, run for 10,000 generations, seeding them with rate estimates recovered from the previous ML analyses as priors, and discarding the first 1000 generations as burn-in. PP distributions of all parameters were summarized using the diversitree Rpackage. We also fitted the best model and compared it to the second best model across all RD trees to check that model selection was robust to phylogenetic uncertainty.

Geographical Evolution versus Habitat Evolution
To investigate the relative roles of geographical and environmental divergence during the diversification of Oxera, we computed age-range correlations (Fitzpatrick and Turelli 2006;Warren et al. 2008) where metrics of range overlaps, range asymmetries, and habitat distances computed between species pairs are compared with their estimated time of divergence (from the phylogeny). These metrics have been recurrently used to diagnose different speciation and post-speciation scenarios (Anacker and Strauss 2014;Grossenbacher et al. 2014).
Species geographical ranges were retrieved from the geographical data set by creating for each species a single convex polygon enclosing its locations using the "convex hull" tool under QGIS. For narrow endemic species, with one or two occurrences, we applied a 469 buffer of 0.5 km around each location, and used a range of 3.141592 × 0.5 2 km 2 for the former, and of 6.283184 × 0.5 2 km 2 for the latter. From these distribution data, range overlaps, range asymmetries, and habitat distances metrics were computed between all possible pairs of species. Range overlap was defined as the area occupied by two species divided by the area of the more narrowly ranging species, ranging from 0 (no co-occurrence) to 1 (complete co-occurrence). Range asymmetry was calculated as the area of the widerranging species divided by that of the smaller (ranging from 1 to infinity). The habitat similarity between two species was assessed by calculating their Euclidian distance (habitat distance) based on the eigenvector values of the first n axes explaining 90% of the cumulate variance retrieved from the environmental PCA (see above). These geographical variables (areas of all Oxera species, range asymmetries and range overlaps of all Oxera species pairs) are available on Dryad.
Age-range correlation analyses were performed using an extended version of the "age.range.correlation" function of the phyloclim R package (Heibl and Calenge 2013; https://github.com/danlwarren/arc-extensions/ blob/master/age.range.correlation.2.R; last accessed 18 October 2017), applied to our pMCCT and the RD trees. We launched 1000 iterations for the Monte Carlo resampling procedure to create the null hypothesis of no relationship between phylogenetic relatedness and range overlap, range asymmetry, or habitat distance, and tested it. We calculated linear regressions, whose slopes and intercepts indicate speciation mode (allopatric vs. sympatric) or habitat evolution (conservatism vs. divergence) through time. For each metric a "super-p-value" was calculated, corresponding to the tree proportion across the RD trees for which the linear regression was significant. We also assessed the frequency of range overlaps across all random species pairs and across all sister species pairs using a kernel density plot.
Finally, we estimated the respective importance of geographical versus habitat divergence and whether local co-existence is possible between closely related Oxera species. To do so, we superimposed species geographical overlaps to their habitat preferences (summarized into diagrams combining information on their geological, elevation, and vegetation attributes) within each Oxera subclade (as delimited in Barrabé et al. 2015) to determine whether closely related species could possibly occur in sympatry.

Phylogenetic Inference
The Bayesian MCMC analyses of the three species sampling (Lamiales, Lamiaceae, and Oxera data sets) recovered the same following robust phylogenetic relationships (Supplementary Appendices S5-S7 available on Dryad). The sister genera Hosea and Oxera formed a well-supported clade that was sister to the Clerodendrum clade (PP = 1). Analyses based on the Oxera data set provided more resolution and detail on the internal Oxera relationships, which are congruent with the ones established in Barrabé et al. (2015) (Supplementary Appendix S5 available on Dryad). The Australasian O. splendida was placed as sister to the New Caledonian radiation (PP = 0.72), which included all New Caledonian and nested Pacific subclades. We retrieved seven well-supported New Caledonian subclades (PP = 1; Fig. 2a), namely the oblongifolia, neriifolia, pulchella, baladica, subverticillata, robusta, and sulfurea subclades.
A few supported incongruences (PP ≥ 0.95) were observed between the phylogenetic trees reconstructed separately from the six nuclear loci and from the six plastid loci (Supplementary Appendix S8 available on Dryad), but these incongruencies were all located near the tips (i.e., between close species). This discordance did not affect the robustness of the seven New Caledonian subclades previously recovered. These minor incongruences underscore the importance of accounting for phylogenetic uncertainty by using the 100 RD trees in most of our evolution and geographical analyses (see Material and Methods section).

Oxera Diversification
The lineage through time plot highlighted a gradual and exponential increase of the number of Oxera lineages through time (Fig. 2a). The net conservative diversification rates are summarized in Tables 2 and 3. For the New Caledonian radiation, these rates were 0.2 species/Myr for an extinction of 0.9 and 1.04 species/Myr for a null extinction. With respect to per-unit-area and perunit-log(area) they were estimated between 1.1×10 −5 and 5.6×10 −5 species/Myr/km 2 , and between 0.02 and 470 SYSTEMATIC BIOLOGY VOL. 68 FIGURE 2. Evolutionary history of Oxera in New Caledonia, and of its habitat preferences, pollination and dispersal niches, and life forms. a) BEAST chronogram of the Oxera data set with superimposed lineage through time plot. Median age estimates and blue node bars (corresponding to the 95% highest posterior densities) are indicated for each node of interest. Bayesian PP from the BEAST analyses (on the right) and MCMC Bayesian analyses (on the left) are indicated with an asterisk for each node of interest when >0.95, and with a hyphen when <0.95. Vertical light grey rectangles indicate periods of intense weathering as inferred from Chevillotte et al. (2006). b) Absolute shift timings of the three niches and life forms inferred from stochastic character mappings, with global delta 18 O variation through time (in grey; from Zachos et al. 2001) and net diversification rate variation through time of Oxera, as inferred from BAMM analyses (in dark blue), superimposed. c) Phylogenetic signal for habitat preferences, pollination and dispersal niches, and life forms inferred from estimations of the and Pagel statistics. 0.106 species/Myr/log(km 2 ), respectively. The BAMM analyses inferred two early instances of changes in the diversification rate. The phylorate plot showed an increase in mean diversification rates at the crown node of Oxera and the second deepest node of the Clerodendrum clade, whereas Hosea demonstrated a decrease in these rates (Supplementary Appendix S12a available on Dryad). Five of the six most probable 95% credible   showed significant rate shifts mainly located on the two consecutive branches leading to the New Caledonian Oxera radiation, and those leading to most of the Clerodendrum clade (Supplementary Appendix S13 available on Dryad). Both shifts were depicted in the marginal shift probability tree, where these latter branches were significantly longer than most others, indicating a high probability that these shifts occurred on their respective branches (Supplementary Appendix S12b available on Dryad). The Oxera rates through time plots highlighted that speciation rates remained constant (with a mean value of ca. 0.71 species/Myr). Extinction rates decreased slightly at the radiation onset but remained constant after (mean value of ca. 0.31 species/Myr). The resulting net diversification rates increased early in the radiation onset and later remained quite constant (mean value of ca. 0.4 species/Myr; Fig. 2b, Supplementary Appendix S12c available on Dryad). The ML diversification analyses selected the pure birth model as best-fitting our data for the genus Oxera, with a net diversification rate of 0.575 species/Myr (Supplementary Appendix S14 available on Dryad). Accordingly, the pure birth model was the best selected model on 93% of the RD trees.

Niche Ordination
The multivariate analyses allowed the separation of four, three, and five clearly distinct ecological syndromes from the species' environmental, flower, and fruit data sets, respectively (Fig. 1c). The first four, nine, and six axes of the PCA analyses (explaining 90 % of the cumulate variance), respectively, were retained for subsequent comparative analyses. The three major habitats occupied by Oxera were forests on non-ultramafic rocks, forests on ultramafic rocks, sclerophyllous vegetation, and a fourth class was defined for ubiquist species. The three pollination syndromes characterized species with actinomorphic, gullet, and bilabiate floral morphs each adapted to bat, bird, moth/butterfly pollination, respectively. For dispersal syndromes the PCA axis markedly separated species with small fruits ≤ 2 cm (on the left) from those with large fruits >2 cm (on the right; Fig. 1c). The five dispersal syndromes discriminated species with potato-like (consumed by cassowaries and bats), paddle-like (by Columbidae), small rounded (by Columbidae and/or small birds such as Sturnidae), large rounded (by Columbidae), and thorny fruit morphs (presumably by large geckos).

Niche Evolution
Our estimates of both Pagel statistics showed distinct evolutionary patterns ( Fig. 2c; Supplementary Appendix S15 available on Dryad). We found low values of , indicating low phylogenetic signal (i.e., a recent accelerated evolution), for the habitat preference (Fig. 2c), and high phylogenetic signal for the dispersal niche and life form, suggesting that these niche characteristics have diversified early. For the pollination niche, estimates of were more spread out, with a lower 5% quantile of 0.2, a higher 95% quantile of 1.01, and a median value of 0.99; this also indicated a relatively high phylogenetic signal (Fig. 2c). Estimates of the statistic were close to 0 for the habitat preference (indicating a punctual evolution), and close to 1 for the dispersal niche and life form, respectively (indicating more gradual evolution). For the pollination niche, scaled between 0.35 (5% quantile) and 0.99 (95% quantile), with a median value of 0.69.
The best-fit evolutionary model identified in our ancestral state reconstruction analyses was the "ER" model for each niche/trait data set (Supplementary Appendix S16 available on Dryad). The reconstructed ancestral habitat preference of Oxera was ambiguous: forests on both ultramafic and non-ultramafic rocks were together reconstructed with a high likelihood (forested vegetation was recovered as ancestral for this node with cumulated likelihoods of both forest types of ca. 0.57). However, the ancestral pollination, dispersal, and life form were inferred with much less ambiguity; the Oxera ancestor was very likely a robust liana with gullet flowers producing small rounded fruits (Supplementary Appendix S17 available on Dryad). Analyses of stochastic character mapping allowed us to estimate that very few pollination, dispersal, and life form shifts occurred during the evolution of Oxera (up to two shifts per time bin of 0.5 Myr; Fig. 2b), and that most of these trait shifts occurred between 2 Ma and 0.5 Ma for pollination, and following 3 Ma and 3.5 Ma for dispersal and life form, respectively. Habitat shifts occurred between 9 Ma and present, but from 4 Ma onward their number increased markedly, finally reaching 17 shifts between 0.5 Ma and the present.
The best-fitting model of trait-dependent species diversification identified for habitats had speciation and transition rates all equal, a null extinction rate (Supplementary Appendix S4 available on Dryad), and no discernible effect of any habitat type. For the pollination syndromes, dispersal types and life forms, the best-fitting model was that with a null extinction, transition rates all equal, and a particular character state in which speciation rates were significantly different from other states (Supplementary Appendix S4 available on Dryad). These results were very robust to phylogenetic uncertainty, as 100% of the RD trees consistently yielded the same best models. Slender lianas had on average significantly higher speciation rates comparing to other life forms, actinomorphic flowers, and potato-like fruits had on average significantly lower speciation rates comparing to other floral and fruit types (Supplementary Appendix S18 available on Dryad).

Geographical Evolution
The range overlap between all species pairs was generally low (mean value of 0.187, and a standard deviation of 0.357), suggesting a high level of allopatry between all species. The kernel density plot highlighted two range overlap frequency peaks, the highest located for a range overlap of 0 and the shortest to 1 (Fig. 3). Age range correlation analyses showed that the slope of range asymmetry through time was significantly positive for linear regressions (all P-values <0.05 across RD trees, Supplementary Appendix S19 available on Dryad). For range overlap and habitat distance, linear regression slopes with divergence times were negative in both cases but not significant.
Different Oxera subclades highlighted distinct patterns of geographical overlap ( Fig. 4a and b). The baladica subclade was entirely and strictly allopatric since none of its species overlapped geographically. The oblongifolia subclade was completely sympatric as all its species overlapped. The other subclades (i.e., "neriifolia" "pulchella" "robusta" "subverticillata" and "sulfurea" subclades) exhibited more complex patterns (Fig. 4a). Overall, of the 96 possible pairs of species encountered within the seven New Caledonian subclades, we recorded 23 cases of partial or entire geographical overlap ( Fig. 4a  and b), but most of them involved species occurring in markedly distinct habitats. In fact, among all 96 species pairs, only five pairs involved geographically overlapping pairs of species growing in the same habitat, suggesting possible local co-existence (i.e., O. oblongifolia and O. glandulosa, O. coriacea, and O. palmatinervia, O. merytifolia, and O. subverticillata, O. gmelinoides, and O. pancheri). However, none of these species pairs have been observed co-occurring in the field (Gâteblé, pers. obs.).

An Early Burst of Diversification
Our study shows that the genus Oxera is a quite recent and rapid plant radiation that initiated during the late Miocene, and that most of the New Caledonian subclades within Oxera began to diversify in the early Pliocene (Table 2, Fig. 2). It is noteworthy that this New Caledonian radiation occurred long after the emergence of Grande Terre about 37 Ma (Cluzel et al. 2012). Thus, the radiation of Oxera originated in the archipelago via long-distance or stepping-stone dispersal, but was not triggered by Gondwanan fragmentation. This continental dispersal to the New Caledonian archipelago was probably aided by the ancestral dispersal type (i.e., small rounded fleshy fruits; Supplementary Appendix S17 available on Dryad), which may be the most efficient dispersal mode in the group. Following this dispersal to New Caledonia, the successful establishment of ancestral Oxera may have been facilitated by the climbing habit of its ancestor (robust woody liana identified as the ancestral state; Supplementary Appendix S17 available on Dryad), which is believed to allow the occupation of a broad range of local environments (Gianoli 2004). The finding that the ancestral colonizer of extant New Caledonian Oxera was reconstructed as having a liana life-form with fleshy fruits makes it a quite unique case since the liana form is generally associated with wind dispersal (Givnish 2010), implying dispersal limitations over water, and that the liana niche in islands remains often empty.
Diversification analyses identified an early burst of diversification (BAMM analyses; Fig. 2a and b, Supplementary Appendix S12 available on Dryad), likely located between the divergence of the Oxera ingroup and its sister lineage Hosea (ca. 12 Ma) and the crown node of Oxera within New Caledonia (ca. 4.5 Ma). This late Miocene/early Pliocene burst of species diversification indicates that Oxera underwent an initial period of rapid cladogenesis, especially in New Caledonia, concordant with an initial decreasing of extinction rates; Figure 2b, Supplementary Appendix S12 available on Dryad. This coincides with the origination and early diversification of other major New Caledonian radiations (e.g., Psychotria clade NC2), and especially lineages typical of sclerophyllous vegetation types (e.g., Thiollierea; Table 2).
This period of intense species diversification is concordant with the global onset of Pliocene glaciations (Zachos et al. 2001) as well as the temperature declines and severe aridification recorded in the South Pacific during the Pliocene (Karas et al. 2011). It also coincides with the timing of intense soil morphogenesis and weathering which is estimated to have occurred around 5.5 Ma in New Caledonia ( Fig. 2a and b; Chevillotte et al. 2006). Such climatic events have been suggested to have impacted the evolutionary history of many plant lineages by leading to the extinction of certain species (thus alleviating competitive effects), and the emergence of new biotas (e.g., more sclerophyllous and open vegetation), especially since the Miocene (Donoghue and Edwards 2014). It is important to note that some of these novel habitats exhibit a patchy spatial distribution, which has likely increased genetic divergence. This ostensibly has generated vacant niches, created ecological opportunities, increased spatial divergence and thus triggered speciation in many Australasian lineages, as highlighted in Australia for Banksia, Eucalyptus, and Allocasuarina (Crisp et al. 2004;Gallagher et al. 2003). The early rapid diversification of Oxera suggests that similar climatic and physical events fostered diversification within Oxera shortly after its arrival into New Caledonia.  Numbers in brackets refer to tree node numbers in Figure 4. Clades with similar combinations are highlighted in grey.

An Explosive Species Diversification, with Early Divergence in Biotic Niches
The rate of species diversification estimated for Oxera is exceptionally high, and represents the fastest known New Caledonian plant radiation, rivalled only by Thiollierea (Manns et al. 2012; Table 2). Estimated speciation rates within Oxera are comparable to the fastest island plant radiations documented in previous studies (e.g., Bidens in Hawaii or Echium in Macaronesia; Table 3). The detected early burst of diversification has been followed by a rather constant diversification, as inferred from BAMM and ML Treepar analyses even under weak divergence origination scenarios (Supplementary Appendices 12 and 14 available on Dryad); this may reflect the ongoing challenge within the evolutionary biology community to realistically model diversification processes, especially extinction (Stadler 2011;Morlon 2014). Nevertheless, it seems that since its origination Oxera has continued to diversify under constant speciation and rather null extinction, as exemplified by the exponential increase of species lineages through time (Fig. 2a).
Trait-dependent diversification models suggest that no specific adaptation has disproportionately favored speciation in Oxera (all biotic traits had equal speciation rates; Supplementary Appendices S4 and S18 available on Dryad), excepting perhaps the slender liana life form (acquired with the emergence of the neriifolia, pulchella, and oblongifolia subclades; see Table 4, Fig. 2a, Supplementary Appendix S18 available on Dryad), which is found in open sclerophyllous vegetation and exhibits high speciation rates probably because open habitats have increased in concert with Pliocene and Pleistocene cooling while remaining spatially fragmented (see below). However, the proliferation of species within Oxera might be related to the capacity of the genus to adopt different growth forms and floral and fruit morphologies, as hypothesized for angiosperms in general (Ricklefs and Renner 1994), probably in response to varied ecological opportunities available after its colonization of New Caledonia.
The early diversification of Oxera seems to have been mainly associated with niche shifts driven by pollination syndromes, dispersal strategies, and lifeform variety (Fig. 2b). This suggests that the early diversification burst within the clade coincides with ecological diversification, particularly in regards to pollen vectors, dispersal agents, and light requirements. Although it is not possible to determine whether these niche shifts were directly involved in the process of lineage divergence (ecological speciation), or have rather been caused by post-speciation adaptation to diverse ecological opportunities, these results suggest that diverse ecological opportunities have likely facilitated the diversification of Oxera within New Caledonia. The decelerated diversification of pollination syndromes, dispersal niches, and life-forms (revealed by a rather high phylogenetic signal; Fig. 2c) likely indicates that all or most of these niches previously vacant were ultimately filled during the early history of the clade. Our comparative analyses indicate that from the late Pliocene to early Pleistocene, the major Oxera lineages (i.e., the seven New Caledonian subclades plus O. morierei) had already originated (Fig. 2a), most of them occupying a unique niche combination. It is important to note that all recovered Oxera subclades were robust to contrasting phylogenetic signal between nuclear and plastid genes (Supplementary Appendix S8 available on Dryad). Hence, each Oxera subclade robustly holds a unique trait combination (Table 4), excepting perhaps the neriifolia and oblongifolia subclades, which similarly exhibit slender lianas, bilabiate flowers and small rounded fruits. This evolutionary divergence in life forms seems not to be correlated with habitat divergence, since these two characteristics have evolved largely independently across Oxera (Fig. 2b and c), suggesting that different growth forms occupy distinct functional niches related to light requirements within habitats.
Our study thus detected strong ecological segregation among Oxera subclades and species (Fig. 1c), but the precision of our morphological data may not be sufficient to discriminate finer ecological syndromes encountered in Oxera. Compared with our results, the genus may exploit an even larger diversity in pollinators and dispersers than highlighted in this study. For example, the actinomorphic flowers of Pacific species are sometimes visited by honeyeaters (Belcher and Sibson 1982;de Kok 1997). Gullet flowers with small and narrow tubes, such as the ones of O. sessilifolia, O. garoensis sp. nov. ined., and O. doubetiae sp. nov. ined., may be more efficiently pollinated by Diptera or Coleoptera rather than birds (Jourdan, pers. comm.). It would be illuminating to conduct detailed investigations on the pollination biology of different Oxera species using more precise measurements of floral morphology characters. This could reveal that pollination syndromes have in fact spurred diversification during the evolution of Oxera.

Recent Diversification Driven by the Interplay of Habitat
Shifts and Allopatry Our comparative analyses further suggest that since the end of Pliocene and the onset of the Pleistocene, the diversification of Oxera has been impacted by several punctual shifts in habitat preference (Fig. 2b  and c). These events coincided with another intense period of soil morphogenesis and weathering events in New Caledonia which occurred between 3 and 2.5 Ma ( Fig. 2a and b; Chevillotte et al. 2006), as well as several climatic oscillations documented in the south-western Pacific (Dodson and Macphail 2004). These events have likely caused myriad spatial rearrangements of New Caledonian biotas (Hope and Pask 1998, Stevenson et al. 2001, Stevenson 2004Stevenson and Hope 2005), thus generating novel ecological opportunities, and hence the recent diversification of Oxera.
Phylogenetic patterns of geographic range overlap and asymmetry between pairs of species show that geographical isolation seems to have played a prominent role during the diversification of Oxera. This suggests that neutral divergence may have been at least as important as ecological speciation. Most species pairs within Oxera have clear allopatric distributions (null range overlap more frequent through all species pairs; Fig. 3). This is exemplified by the baladica subclade, which appears to have diversified almost entirely in allopatry (Fig. 4). Few species pairs clearly originated from sympatric speciation, as suggested by complete range overlaps and trait divergence (Figs. 3 and 4). Except for the oblongifolia subclade, the geographical range overlap between species remains remarkably low (<0.4; Supplementary Appendix S19 available on Dryad) indicating that given the recent radiations of Oxera, most sister species have remained in allopatry, and that secondary sympatry may not have yet occurred in most sister species pairs.
The examination of geographic and ecological overlap between species pairs suggests that local species coexistence may in fact be very rare in Oxera. Among the 96 putative pairs of species (within the seven New Caledonian subclades), only five pairs potentially coexist locally, although this has never been noted after years of field observations (Gâteblé, pers. obs.). For example in the oblongifolia subclade (which exhibits largely sympatric species), each species grows in a specific set of environmental conditions, excepting the ubiquist O. oblongifolia (Fig. 4). Contemporary species coexistence is probably prevented by the recent interplay of allopatry and parapatry due to habitat shifts and environmental oscillations that have caused cycles of forest expansion and retraction. With Oxera being mainly a forest-adapted group (25 of 38 species confined to forests, and forest habitat identified as ancestral; Supplementary Appendix S17 available on Dryad), it is likely that vegetation cycles along the complex topographical gradients of New Caledonia have triggered spatial divergence and probably speciation. This process has been further amplified by the limited dispersal ability of Oxera since its dispersers are relatively sedentary. Such processes of vegetation oscillations along steep physical gradients, sometimes creating barriers to dispersal, is now considered to be a major driver of clade diversification, as demonstrated in recent theoretical simulations (Aguilée et al. 2012).
CONCLUSION: IS OXERA A CLEAR ADAPTIVE RADIATION FROM NEW CALEDONIA? The genus Oxera constitutes a remarkable example of a rapid diversifying lineage that originated from a single ancestor. We observed that the group experienced an early burst of species diversification, followed by an exceptionally high and constant species diversification rate (Fig. 2b). This signature is considered typical of adaptive radiations for some authors Slowinski and Guyer 1993;Soulebeau et al. 2015) but not all (Givnish 1997(Givnish , 2015, mainly because it is usually thought that early diversification bursts are driven by the exploitation of new ecological opportunities, which allows lineages to enter new adaptive zones. The early history of Oxera seems to have been paved by niche shifts toward novel pollination, dispersal and life-form traits. These niche shifts later decelerated, perhaps due to a saturation of available niche space. However, our analyses do not highlight any diversity dependent diversification processes, which are often considered as another hallmark of adaptive radiations (Rabosky 2009). This is likely because other ecological opportunities have helped further diversification of the group after its initial phase of biotic niche exploration. More recent opportunities likely emerged through landscape formation and environmental oscillations which may have triggered allopatric and parapatric speciation processes. In fact, the explosive diversification of Oxera within New Caledonia appears to have been triggered by varied ecological opportunities during the course of evolution. Adaptive radiation should not be narrowly defined through quantitative tests or by quantifying diversification tempos, since some well-known adaptive radiations have little or no effect on speciation, or even a negative effect (Givnish 2015), and are therefore more identified by high rates of morphological and ecological diversification (Givnish 1997;Sanderson 1998); this remains to be investigated for the Oxera diversification.
The explosive radiation experienced by Oxera thus constitutes a remarkable model for studying the drivers involved in the evolution of the New Caledonian flora, especially for recent epochs. The genus is an unclassifiable, multifaceted radiation seemingly triggered by both adaptive and neutral processes (adaptive, geographical, climatic, etc.), and driven by several classes of biotic 478 SYSTEMATIC BIOLOGY VOL. 68 and abiotic factors. Oxera could be therefore regarded as a partial adaptive radiation since changing ecological opportunities have contributed to its radiation, with an early diversification of biotic traits (life-form, pollination, dispersal), but more recently driven by the joint action of habitat changes and allopatric speciation leading to scant local coexistence between species. Although different diversification drivers may have been overlapping in time, the early stage of Oxera diversification could be then viewed as a leapfrog radiation (sensu Chase and Palmer 1997), where diversification drivers have relayed to one another at different time periods, but geographic processes may have later been an equally important driver of speciation. Indeed, past climatic history of New Caledonia has also acted as a major evolutionary force, with two pivotal events, a dramatic aridity at the early Pliocene leading to the release of newly vacant niches, and the accentuation of glacialinterglacial cycles from the early Pleistocene, causing repeated biota rearrangements. From our study, it is legitimate to ask whether any other case of adaptive radiation can be clearly demonstrated in old islands system such as New Caledonia, since the distinctive hallmark of past adaptive divergence is progressively blurred by more recent environmental fluctuations creating new ecological opportunities.

SUPPLEMENTARY MATERIAL
Data available from the Dryad Digital Repository: http://dx.doi.org/10.5061/dryad.mm1ft40 for DNA alignment matrices and geographical variables. The resulting trees are available from TreeBASE (ID 23353).