# Phylogenetic measures reveal eco-evolutionary drivers of biodiversity along a depth gradient

David Eme, M. Anderson, E. Myers, C. Roberts, L. Liggins

<figure><img src="image_1.png" /></figure>

## To cite this version:

David Eme, M. Anderson, E. Myers, C. Roberts, L. Liggins. Phylogenetic measures reveal eco-evolutionary drivers of biodiversity along a depth gradient. Ecography, 2020, 43 (5), pp.689-702. 10.1111/ecog.04836. hal-03616036

HAL Id: hal-03616036

https://hal.inrae.fr/hal-03616036

Submitted on 22 Mar 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Research

# Phylogenetic measures reveal eco-evolutionary drivers of biodiversity along a depth gradient

D. Eme, M. J. Anderson, E. M. V. Myers, C. D. Roberts and L. Liggins

D. Eme (https://orcid.org/0000-0001-8790-0412) ✉ (david.eme3@gmail.com), M. J. Anderson, E. M. V. Myers (https://orcid.org/0000-0001-5343-349X), New Zealand Inst. for Advanced Study (NZIAS), Massey Univ., Albany Campus, Auckland, New Zealand. Present address of DE: IFREMER, unité Ecologie et Modèles pour l'Halieutique, Nantes, France. – C. D. Roberts, Museum of New Zealand Te Papa Tongarewa, Wellington, New Zealand. – L. Liggins (https://orcid.org/0000-0003-1143-2346), School of Natural and Computational Sciences (SNCS), Massey Univ., Auckland, New Zealand, and Auckland War Memorial Museum, Tāmaki Paenga Hira, Auckland, New Zealand.

## Ecography

43: 1–14, 2020

doi: 10.1111/ecog.04836

Subject Editor: Thiago Rangel

Editor-in-Chief:

Jens-Christian Svenning

Accepted 20 December 2019

<figure><img src="image_3.png" /></figure>

Energy and environmental stability are positively correlated with species richness along broad-scale spatial gradients in terrestrial ecosystems, so their relative importance in generating and preserving diversity cannot be readily disentangled. This study seeks to exploit the negative correlation between energy and stability along the oceanic depth gradient to better understand their relative contribution in shaping broadscale biodiversity patterns. We develop a conceptual framework by simulating speciation and extinction along energy and stability gradients to generate expected patterns of biodiversity for a suite of complementary phylogenetic diversity metrics. Using a time-calibrated molecular phylogeny for New Zealand marine ray-finned fishes and a replicated community ecological sampling design, we then modelled these metrics along large-scale depth and latitude gradients. Our results indicate that energy-rich shallow waters may be an engine of diversity for percomorpha, but also suggest that recent speciation occurs in ancient fish lineages in the deep sea, hence questioning the role of energy as a key driver of speciation. Despite potentially facing high extinction early in their evolution, ancient phylogenetic lineages specialized for the deep-sea were likely preserved by environmental stability during the Cenozoic. Furthermore, intermediate depths might be a 'museum' (or zone of overlap) for distinct lineages that occur predominantly in either shallow or deep-sea waters. These intermediate depths (500–900 m) may form a 'phylogenetic diversity bank', perhaps providing a refuge during ancient (Mesozoic) extreme anoxic events affecting the deep sea and more recent (Pliocene–Pleistocene) climatic events occurring in shallow ecosystems. Finally, the phylogenetic structures observed in fish communities at intermediate depths suggest other processes might restrict the co-occurrence of closely related species. Overall, by combining a conceptual framework with models of empirical phylogenetic diversity patterns, our study paves the way for understanding the determinants of biodiversity across the largest habitat on earth.

Keywords: community phylogenetics, deep sea, energy gradient, latitudinal gradient, museum/cradle, New Zealand marine Actinopterygii

<figure><img src="image_4.png" /></figure>
**Introduction**

Despite harbouring some of the strongest environmen-
tal gradients on earth, little is known regarding mecha-
nisms shaping biodiversity along the ocean's depth gradient
(Ramirez-Llodra et al. 2010). Along other broadscale gradi-
ents (e.g. latitude, elevation), biodiversity tends to increase
with increasing energy, temperature and environmental stabil-
ity (e.g. at low latitude) (Mittelbach et al. 2007, Jetz and Fine
2012, Fine 2015). Energy, temperature and environmental
stability are, however, positively correlated along latitudinal
and elevational gradients (Fine 2015, McClain and Schlacher
2015, Valentine and Jablonski 2015). In contrast, while energy
and temperature both decrease with increasing depth, envi-
ronmental stability increases (Fig. 1a) (Ramirez-Llodra et al.

2010, Fine 2015, Valentine and Jablonski 2015). Since the
global anoxic events (OAEs) that caused mass extinction
during the Mesozoic (Priede and Froese 2013) as recently
as 94–93 Ma ago (Jenkyns 2010), environmental conditions
in the deep sea have remained reasonably stable (Priede and
Froese 2013). Owing to the buffering effect of the water
column (seasonal variation rarely occurs deeper than 400 m
[Fiedler 2010]), deeper ocean environments would have been
less affected by the great seasonal variation of the Eocene/
Oligocene boundary (Ivany et al. 2000) and climatic oscil-
lations during the Miocene and Pleistocene (e.g. glacial/
inter-glacial periods) than shallow-water environments
(Pillans et al. 1998, McClymont et al. 2016). Thus, study-
ing patterns of biodiversity along the depth gradient offers an
opportunity to disentangle the often-confounded drivers of

<figure><img src="image_5.png" /></figure>

Figure 1. (a) Correlations between energy (black) and stability (gray) along latitudinal and elevation gradients and along the depth gradient. (b) Conceptual framework based on simulations (Supplementary material Appendix 1) proposed to disentangle the role of energy versus stability on evolutionary processes shaping biodiversity over depth. The central tree represents a hypothetical phylogenetic tree of a community (i), and each quadrant shows what phylogenetic structures and associated changes in phylogenetic diversity values that would be expected for communities at different positions along the energy and stability axes, given their expected effects on speciation and extinction rates. SR = species richness, PD = phylogenetic diversity, MPD = mean pairwise distance, VPD = variance of the pairwise distance, MNTD = mean nearest taxonomic distance. Right: the predicted change in phylogenetic diversity measures (PD, MPD, VPD, MNTD) versus depth when (ci, cii) energy alone, (d) stability alone, or (ei, eii) both stability and energy are driving patterns in biodiversity.
energy/temperature and stability on biodiversity (Fine 2015, Valentine and Jablonski 2015).

Examining the phylogenetic structures of communities, rather than just taxonomic species richness, has enhanced our understanding of broadscale biodiversity patterns and their drivers (Cavender-Bares et al. 2009, Mouquet et al. 2012). Several alpha (local) diversity measures have been proposed to capture different aspects of this phylogenetic structure (Vamosi et al. 2009, Pavoine and Bonsall 2011, Tucker et al. 2017). For example, phylogenetic diversity (PD) captures the overall amount of evolutionary history in a sampled community as a sum of the branch lengths linking species through a phylogeny (Faith 1992), whereas mean pairwise distance (MPD, Webb 2000) represents the average phylogenetic distance (or average time divergence) among species, and is influenced by ancient diversification events (Mazel et al. 2016). In contrast, mean nearest taxonomic distance (MNTD, Webb 2000) estimates the average phylogenetic distance among the closest relatives within a community, and is therefore influenced by recent speciation events (Mazel et al. 2016). Last, whereas the variance of the pairwise distance (VPD, Clarke and Warwick 2001) measures variability in phylogenetic distances among species, influenced by both recent and ancient divergences, the variance of the nearest taxonomic distance (VNTD, Tucker et al. 2017) focuses on closest relatives, so is not influenced by ancient divergence. As each measure draws on unique information from the phylogenetic relationships among species, documenting patterns for several of these measures simultaneously can provide novel insights into the mechanisms shaping biodiversity along broad-scale gradients, such as depth.

Here we present a simple conceptual framework (Fig. 1b) that predicts patterns of phylogenetic alpha diversity (hereafter phylo-α) under different scenarios of speciation and extinction along broad-scale geographical/environmental gradients (see Supplementary material Appendix 1 for simulation details). Our framework draws on macroecological hypotheses commonly proposed for both terrestrial and marine systems that suggest a positive relationship between energy/temperature and speciation, and a negative relationship between long-term environmental/geological stability and extinction (Rohde 1992, Allen et al. 2006, Davies et al. 2007, Mittelbach et al. 2007, Gillman et al. 2009, Condamine et al. 2013, 2019, Rolland et al. 2014, Siqueira et al. 2016). In our framework we explore a broad spectrum of extinction rates/scenarios for a given speciation level and vice-versa (Fig. 1b, Supplementary material Appendix 1), and acknowledge the potential for a negative relationship between energy and extinction (Evans et al. 2005), and that long-term environmental stability may also increase speciation rate through specialization (Jansson and Dynesius 2002).

We began with a hypothetical community phylogenetic tree generated with intermediate levels of speciation and extinction (Fig. 1bi) and portray how this structure, and associated phylogenetic diversity measures, would change

with increasing or decreasing speciation or extinction rates. For example, high speciation and low extinction will cause trees to have higher PD and MPD values, but shorter distances among close relatives (low MNTD) and lower VPD (Fig. 1bii). This pattern typifies a ‘cradle’ (generating new species) as well as a ‘museum’ (preserving old lineages), likely to occur in regions of high energy and low environmental/geological variability (e.g. tropical regions, Jablonski et al. 2006, Pyron and Wiens 2013, Rolland et al. 2014). An increase in the extinction probability (i.e. more frequent mass extinction events) in conjunction with a high speciation rate, as expected in energy-rich but unstable regions, would generate phylogenies with long basal branches creating distinct lineages with many short terminal branches. This phylogenetic structure would yield lower values for PD and MPD, increases in VPD, and MNTD would remain stable and low (Fig. 1biii). High extinction in association with low speciation, as expected in highly unstable and energy-poor regions, would generate depauperate phylogenies having a mixture of long and short branches, decreasing the PD and MPD, but increasing the VPD and the MNTD (Fig. 1biv). Finally, under a scenario of low speciation and extinction, as expected in stable but low-energy areas, species-poor phylogenies with an even branch-length distribution would be generated, yielding communities with relatively low PD and VPD, but higher MPD and MNTD (Fig. 1bv).

Considering this conceptual framework in the context of the depth gradient, we made specific predictions regarding expected changes in phylo-α according to underlying drivers. First, if energy/temperature boost speciation and decrease extinction (Allen et al. 2006, Pyron and Wiens 2013), PD and MPD should be highest in energy-rich surface waters and decrease with depth, while VPD and MNTD should be lowest in shallow waters, then increase with depth (Fig. 1bii, iv, ci). Alternatively, if energy is not the major determinant of speciation rates in marine systems (Rabosky et al. 2018, O’Hara et al. 2019), but rather only decreases extinction rates (Evans et al. 2005), MNTD should be constant along the depth gradient, PD and MPD should be highest in shallow waters and decrease with depth, while VPD should be lowest in shallow waters, then increase with depth (Fig. 1bii, biii, cii). Second, if long-term stability decreases extinction, and increases the opportunity for speciation through specialization (Jansson and Dynesius 2002, Kozak and Wiens 2012, Salisbury et al. 2012, Rolland and Salamin 2016, Gajdzik et al. 2019), we expect PD and MPD to increase with depth toward more stable environments, while VPD and MNTD should decrease with depth (Fig. 1bii, iv, d). Third, if both energy/temperature and long-term stability yield low extinction and high speciation, PD and MPD would be U-shaped versus depth, while MNTD and VPD would be unimodal (Fig. 1ei). Finally, if energy/temperature does not affect speciation (Rabosky et al. 2018, O’Hara et al. 2019) but only extinction, and long-term stability yields low extinction and high speciation, PD should be consistently low across shallow and intermediate depths before increasing
at deeper depths, and MPD would be U-shaped (extinction decreases in both shallow and deep waters), while MNTD and VPD would decrease with depth (Fig. 1b, bi, bii, bv, eii).

Previous studies of biodiversity versus depth describe either a decline in species richness with depth or a peak at intermediate depths (Ramirez-Llodra et al. 2010, Costello and Chaudhary 2017). Patterns in taxonomic diversity indices (e.g. AvTD) versus depth are also variable (Tolimieri and Anderson 2010, Campbell et al. 2011, Zintzen et al. 2011, 2012). Studies that have inferred evolutionary processes along large-scale gradients in the ocean using near-complete species phylogenies suggest that speciation rates may be higher in marine fishes (Rabosky et al. 2018) and ophiuroids (O'Hara et al. 2019) found in polar regions (possibly owing to instability since the development of icecaps around 41–34 Ma [Crame 2018]) and at depth, challenging the idea that energy is the primary driver of speciation rates in the ocean.

Our study provides the first analysis of complementary phylo-α measures using a time-calibrated phylogeny to examine phylogenetic structures for communities of marine ray-finned fishes sampled along broad-scale ocean gradients. Specifically, we modeled variation in each of five phylo-α indices versus depth (50–1200 m; replicated across latitude) within New Zealand waters to test our predictions regarding the relative roles of speciation and extinction (Fig. 1c-e), potentially driven by energy/temperature and long-term environmental stability. We used multivariate cluster analysis of our phylo-α measures to identify ‘phyloregions’ within which fish communities demonstrate similar phylo-α signatures, hence are likely to have been influenced by similar underlying drivers.

# Material and methods

## Data acquisition

Phylo-α measures of fish communities were calculated based on presence/absence data for 149 teleost species identified from baited remote underwater video (BRUV) apparatuses deployed in a structured replicated sampling design (Zintzen et al. 2017). Fish species were recorded from 2h of video footage obtained by each BRUV deployment. The BRUVs were deployed at seven locations spanning 21° latitude within New Zealand waters: Rangitahua (the Kermadec Islands, KER), the Three Kings Islands (TKI), Great Barrier Island (GBI), White Island (Wis), Kaikōura (KKA), the Otago Peninsula (OTA) and the Auckland Islands (AUC). For each location, 5–6 replicate videos were taken at each of seven targeted depths: 50, 100, 300, 500, 700, 900 and 1200 m.

## Phylogenetic tree building

We built a Bayesian posterior distribution of time-calibrated molecular phylogenies for New Zealand marine ray-finned fishes (class Actinopterygii) following the procedure of Eme et al. (2019). We used their same 15 gene supermatrix,

adding DNA sequences for 14 species never sequenced before, and also added 6 taxa that lacked sequence data. Taxa lacking sequence information were added to the tree reconstruction using conservative topological constraints based on their taxonomic status. The overall supermatrix included 803 species (797 species with DNA and 6 outgroup species not present in New Zealand) and 3679 DNA sequences. The tree reconstruction included 165 soft topological constraints to avoid spurious phylogenetic relationships due to data deficiency (incomplete taxon/gene sampling). The tree was calibrated in absolute time using 44 fossil calibration points. Four independent MCMCs were constructed using BEASTv2.4.7 (Bouckaert et al. 2014) to obtain a posterior distribution of time-calibrated trees and a maximum clade credibility tree (MCCT). We randomly sampled 100 trees from the posterior distribution to quantify phylogenetic uncertainty in downstream analyses (Rangel et al. 2015) (Supplementary material Appendix 2 for laboratory and tree-building methods, and Data availability statement).

## Phylogenetic diversity metrics

We calculated five phylo-α indices to capture different aspects of the phylogenetic structure of communities (Tucker et al. 2017). Although phylogenetic diversity (PD, Faith 1992) is strongly correlated with species richness (SR), it provides insights regarding rates of speciation. For instance, low PD associated with high SR indicates a high rate of speciation (Davies et al. 2007). Mean phylogenetic distance among species (MPD, Clarke and Warwick 1998, Webb 2000) is independent of SR, but is sensitive to large phylogenetic distances contributed by co-occurring species belonging to old and distinct lineages (Mazel et al. 2016), whereas mean nearest taxonomic distance (MNTD) is influenced by recent speciation only (Webb 2000, Mazel et al. 2016). Variation in phylogenetic distances among species within a community (VPD) is independent of SR and MPD (Clarke and Warwick 2001). This regularity index captures more complex phylogenetic structure, such as the presence of distinct old lineages with recent clusters of closely related species (Zintzen et al. 2011), a signature of high extinction in the past (Sanmartín and Meseguer 2016), and/or strong environmental filtering (Cavender-Bares et al. 2009, Zintzen et al. 2011). Finally, the variability of phylogenetic distances among close relatives (VNTD, Tucker et al. 2017) was used to capture the temporal breadth of recent speciation. We also calculated correlations among these metrics for our dataset (Supplementary material Appendix 2).

## Simulation of community phylogenetic structures

To build our conceptual framework (Fig. 1b), we calculated SR and the phylo-α measures from hypothetical trees generated under varying speciation and extinction rates using four sets of simulations. Simulation sets differed by relaxing the extinction rate through time (including mass extinction events) and sampling different fractions of extant
species. Because MNTD and VNTD showed similar trends in our simulation framework, only MNTD is presented in Fig. 1, but both metrics were included in the empirical study. For details regarding simulation procedures and results see Supplementary material Appendix 1 Fig. A1–A8, Table A1; note that Fig. 1b is based on simulations shown in Supplementary material Appendix 1 Fig. A5.

## Statistical analysis

There were 47 depth-by-latitude cells in our full sampling design (there were no BRUV deployments at 1200 m for either White Island or the Auckland Islands, due to poor weather conditions, see Zintzen et al. 2017). From the full BRUV dataset, we pooled the species sampled from BRUVs within each cell, hereafter called the ‘pooled 47-cell dataset’.

We used the MCCT tree to calculate phylo-α values for each metric in each cell and performed model selection with depth and latitude as explanatory variables using Akaike's information criterion with small-sample-size correction (AICc). Model selection performed on average values obtained over 100 trees did not alter the results (r = 0.98–0.99 for all metrics). For each metric, we considered four sets of models (Supplementary material Appendix 2 Table A2): 1) a two-factor ANOVA design with depth and latitude as factors; 2) an ANCOVA design with latitude as a factor and depth as a covariate; 3) an ANCOVA design with depth as a factor and latitude as a covariate; and 4) a regression model with depth and latitude as continuous predictors including quadratic terms for both and their interactions. With each set of models, we also considered simpler sub-set models and retained the model with the lowest AICc. Model selection was done in R (R Core Team) using ordinary least squares models. To quantify phylogenetic uncertainty, we fitted the best model for each phylo-α metric for all 100 trees and computed the mean and variance of parameter estimates and their contributions towards the explained variance. To estimate the proportion of variation explained by phylogenetic uncertainty, we also performed a linear mixed-effects model with trees as a random effect (Supplementary material Appendix 2).

We used clustering methods to identify bioregions of similar phylogenetic characteristics using all five phylo-α measures based on the MCCT. First, we used the ‘find.cluster’ function of the adegenet R package (Jombart 2008) to perform sequential k-means analyses for k = 2–15 clusters, and chose the optimal number of clusters using the recommended ‘diffNgroups’ option (i.e. the biggest difference in Bayesian information criterion values between successive values of k, Jombart et al. 2010). To check the general robustness of our clusters, we also performed hierarchical agglomerative group-average clustering (UPGMA). Next, we repeated k-means clustering for the optimal number of clusters for all 100 trees to accommodate phylogenetic uncertainty. Analyses were repeated excluding PD to decrease the influence of SR. Finally, to identify biodiversity patterns at a finer spatial scale, all analyses were also done using values obtained for each

phylo-α metric at the individual BRUV level, then averaged for each of the 47 cells.

## Results

The MCCT showed general agreement with the tree published by Eme et al. (2019) (Supplementary material Appendix 2 Fig. A9).

### Modeling phylogenetic diversity versus depth and latitude

The phylo-α metrics varied more strongly with depth than with latitude (Fig. 2). All metrics, except VNTD, showed a significant relationship with depth. SR and PD (which were highly correlated, Supplementary material Appendix 2 Fig. A10, A11) decreased sharply with depth between 50 and 300 m, then showed a secondary peak at 900 m. SR and PD had significant unimodal relationships with latitude, peaking around 37°S. MPD showed an unexpected unimodal relationship with depth, peaking at 700–900 m (Fig. 2c). MNTD also had a unimodal relationship with depth, showing a maximum between 300 m and 700 m and a minimum at 1200 m (Fig. 2d), matching the expected relationships in Fig. 1ei. VPD increased sharply with depth, especially beyond 500 m, while VNTD remained fairly constant throughout the depth and latitudinal ranges examined here (Fig. 2f, 1). These patterns echoed phylogenetic structure in species’ occurrences versus depth, with a dominance of Percomorpha in shallow waters being replaced by species belonging to more ancient lineages (Anguilliformes, Notacanthiformes, Alepocephaliformes, Stomiiformes, Gadiformes) at depth (Supplementary material Appendix 2 Fig. A12). In contrast, there was no obvious phylogenetic structure in species’ occurrences versus latitude (Supplementary material Appendix 2 Fig. A13).

Model selection always favoured models with depth and latitude as continuous variables. The parameter estimates, their significance and the explained variance of individual terms retained in the best models remained very consistent across the distribution of 100 trees for all phylo-α metrics (Fig. 3, Table 1), indicating robustness of observed patterns to phylogenetic uncertainty. The best model for all phylogenetic indices, except VNTD, included depth (Table 1) as either a linear or quadratic term (for MPD, MNTD and VPD), and its explained variance was always higher than for latitude, except for PD (Fig. 3, Table 1). A significant interaction between depth and latitude was identified for PD and MNTD (Table 1, Supplementary material Appendix 2 Fig. A14). PD values decreased from shallow to deep communities, particularly at low and intermediate latitudes (Supplementary material Appendix 2 Fig. A14a). For MNTD, the significant interaction between depth and latitude (Table 1) revealed a peak at intermediate depths, which was less pronounced at intermediate latitudes (Supplementary material Appendix 2 Fig. A14c).
<figure><img src="image_6.png" /><figcaption>Figure 2. Relationships of diversity indices with depth (top) and latitude (bottom), when considering the average value estimated over 100 trees on the pooled 47 cells dataset for: species richness (a, g), phylogenetic diversity, PD (b, h); mean pairwise distance, MPD (c, i); mean nearest taxon distance, MNTD (d, j); variance of the pairwise distance, VPD (e, k); and variance of the nearest taxon distance, VNTD (f, l). Black horizontal bars, black dots and boxes show the median, outliers and interquartile range respectively. Whiskers are up to 1.5 times the interquartile range. Blue dots represent the average value per depth/latitude. For visualization purposes only, the blue curves represent the fit of the best model (dashed lines represent non-significant trends) selected based on AICc among three alternative models: model with a single term, a model with single and quadratic terms, and a model with single, quadratic and cubic terms for depth (a–f) or latitude (g–l). The model selection was performed separately for depth and latitude. For species richness, a generalized linear model with a negative binomial error model with a log link was used, for all other metrics we used an ordinary least square models with a Gaussian error, and identity link.</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>Figure 3. Variability of the proportion of explained variance for the different parameters retained in the best models for the phylogenetic diversity measures estimated over the 100 trees, for the pooled 47 cells dataset. Notch, black dots, color dots and boxes show the median, the average, the outliers and the interquartile range respectively. Whiskers are up to 1.5 times the interquartile range.</figcaption></figure>

For PD, the models explained on average 50% of the variance and the interaction between depth and latitude was the most important term in the model (17.9% of explained variance) followed by the quadratic effect of latitude (Table 1, Fig. 3). For MPD, the models explained 51% of the variance and depth was the most important (i.e. 44%, summing depth and depth²). Latitude and its quadratic term remained significant but their contribution was minor (i.e. below 7% when summing latitudinal effects, Table 1). For MNTD, the model explained 38% of the variance and the quadratic term for depth contributed the most (23%), followed by the depth-by-latitude interaction (10%). For VPD, depth was the most important predictor, explaining 70% of model variation (summing depth and depth²). There were sharp increases in VPD along the depth gradient across all latitudes (Fig. 2e, Supplementary material Appendix 2 Fig. A14d). For VNTD, none of the model parameters were statistically significant.

The linear mixed-effects models (with trees as a random factor) yielded identical results regarding the most important model terms (Supplementary material Appendix 2 Table A3). Phylogenetic uncertainty (i.e. the random factor) contributed little to overall variation (0–2%, Supplementary material Appendix 2 Table A3).

Similar results were obtained when analyses were done based on smaller sampling units (individual BRUV units) for all metrics, except for VNTD (Supplementary material Appendix 2 Fig. A15–A17, Table A4, A5). The latter showed a significant positive relationship with depth (Supplementary

material Appendix 2 Fig. A15), explaining 36% of the variance in the best model (when summing depth and depth², Supplementary material Appendix 2 Table A4).

## Phylogenetic bioregions across depth and latitude

K-means clustering of the phylo-α measures for the pooled 47-cell dataset identified four clusters that were robust to phylogenetic uncertainty (Fig. 4a–b, Supplementary material Appendix 2 Fig. A18) (VNTD was removed because it was uncorrelated with latitude or depth). The first phylogenetic bioregion (hereafter ‘phyloregion’) corresponded to low-latitude shallow waters (in red, Fig. 4) characterized by high PD and low MPD, MNTD and VPD, suggesting a scenario of high speciation, and intermediate-to-high extinction by reference to the conceptual framework (Fig. 1, Supplementary material Appendix 1 Fig. A5). This suggested that energy might be a driver of recent and sustained speciation in shallow and subtropical waters (see Discussion). A second phyloregion corresponded to intermediate depths (300–700 m) at low latitudes plus shallow waters (50–500 m) at high latitudes, which had high MNTD and low MPD, PD and VPD (in orange, Fig. 4). A third bioregion (in blue) corresponded to depths around 700–900 m and some shallower depths (500–300 m) at higher latitudes; it was characterized by high MPD, MNTD and intermediate VPD and PD values. These two phyloregions (orange and blue) characterized by intermediate stability and energy (i.e. mesotrophic), suggests a
<table><caption>Table 1. Parameter values and proportion of explained variance for the best models for the five phylogenetic diversity measures with depth and latitude, estimated across 100 trees extracted from the posterior distribution, for the pooled 47 cells dataset. PD: phylogenetic diversity, MPD: mean pairwise distance, MNTD: mean nearest taxon distance, VPD: variance of the pairwise distance, VNTD: variance of the nearest taxon distance. Significant p-values (lower than 0.05) are indicated in bold font (except intercept) as well as parameters with high explained variance.</caption><thead><tr><th>Metrics</th><th>Parameter</th><th>Estimate (SD)</th><th>Estimate range</th><th>p-value (SD)</th><th>p-value range</th><th>Explained variance % (SD)</th><th>Explained variance % range</th></tr></thead><tbody><tr><td rowspan="5">PD</td><td>Intercept</td><td>-1720.3 (118.6)</td><td>-2024.6 to -1352.4</td><td>0.312 (0.033)</td><td>0.24-0.412</td><td></td><td></td></tr><tr><td>Depth</td><td>-3.0 (0.1)</td><td>-3.3 to -2.8</td><td><strong>0.000</strong> (2.8E-05)</td><td>6.4E-05-2.2E-04</td><td><strong>6.8</strong> (0.9)</td><td>4.6-9.0</td></tr><tr><td>Latitude</td><td>205.8 (7.2)</td><td>184.8-222.5</td><td><strong>0.018</strong> (0.003)</td><td>0.0120.028</td><td><strong>9.7</strong> (0.8)</td><td>7.8-11.4</td></tr><tr><td>Latitude<sup>2</sup></td><td>-3.2 (0.1)</td><td>-3.5 to -3.0</td><td><strong>0.003</strong> (5.0E-04)</td><td>0.002-0.005</td><td><strong>15.7</strong> (0.6)</td><td>14.2-17.1</td></tr><tr><td>Depth × Latitude</td><td>0.07 (2.0E-03)</td><td>0.06-0.08</td><td><strong>4.0E-04</strong> (8.0E-05)</td><td>2.2E-0.4-6.3E-04</td><td><strong>17.9</strong> (0.7)</td><td>16.4-19.6</td></tr><tr><td rowspan="5">MPD</td><td>Intercept</td><td>-126.0 (21.7)</td><td>-176.1 to -46.7</td><td>0.41 (0.066)</td><td>0.313-0.712</td><td></td><td></td></tr><tr><td>Depth</td><td>0.1 (9.0E-03)</td><td>0.1-0.14</td><td><strong>0.004</strong> (0.003)</td><td>0-0.013</td><td><strong>41.7</strong> (3.7)</td><td>31.8-50.7</td></tr><tr><td>Latitude</td><td>17.8 (1.2)</td><td>13.8-20.6</td><td><strong>0.025</strong> (0.005)</td><td>0.016-0.043</td><td>0.2 (0.2)</td><td>0-0.9</td></tr><tr><td>Depth<sup>2</sup></td><td>0 (0)</td><td>0-0</td><td>0.152 (0.065)</td><td>0.027-0.299</td><td><strong>2.5</strong> (0.9)</td><td>1.1-5.6</td></tr><tr><td>Latitude<sup>2</sup></td><td>-0.23 (0.02)</td><td>-0.26 to -0.2</td><td><strong>0.023</strong> (0.009)</td><td>0.015-0.039</td><td><strong>6.6</strong> (0.7)</td><td>4.8-8.2</td></tr><tr><td rowspan="5">MNTD</td><td>Intercept</td><td>38.6 (10.7)</td><td>11.8-62.1</td><td>0.493 (0.114)</td><td>0.26-0.82</td><td></td><td></td></tr><tr><td>Depth</td><td>0.4 (0.02)</td><td>0.4-0.5</td><td><strong>0</strong> (0)</td><td>0-0.002</td><td><strong>4.4</strong> (1.6)</td><td>0.9-8.1</td></tr><tr><td>Latitude</td><td>2.3 (0.2)</td><td>1.7-2.8</td><td>0.098 (0.04)</td><td>0.029-0.219</td><td>0.5 (0.5)</td><td>0-2.3</td></tr><tr><td>Depth<sup>2</sup></td><td>0 (0)</td><td>0-0</td><td><strong>0</strong> (0)</td><td>0-0.001</td><td><strong>23.5</strong> (3.5)</td><td>16.1-32.7</td></tr><tr><td>Depth × Latitude</td><td>-6.0E-03 (0)</td><td>-7.0E-03 to -5.0E-03</td><td><strong>0.014</strong> (0.007)</td><td>0.004-0.037</td><td><strong>10.0</strong> (1.4)</td><td>6.8-14.0</td></tr><tr><td rowspan="3">VPD</td><td>Intercept</td><td>4289.4 (379.2)</td><td>3474.267-5441.93</td><td>0.001 (0.001)</td><td>0-0.004</td><td></td><td></td></tr><tr><td>Depth</td><td>-3.8 (1.3)</td><td>-6.829 to -1.006</td><td>0.457 (0.162)</td><td>0.145-0.842</td><td><strong>61.1</strong> (1.7)</td><td>56.7-65.0</td></tr><tr><td>Depth<sup>2</sup></td><td>0.015 (2.0E-03)</td><td>0.011-0.019</td><td><strong>0.001</strong> (0.001)</td><td>0-0.009</td><td><strong>8.8</strong> (1.363)</td><td>5-12.3</td></tr><tr><td rowspan="2">VNTD</td><td>Intercept</td><td>5047.7 (312.2)</td><td>4204.9-5889.5</td><td>0 (0)</td><td>0-0</td><td></td><td></td></tr><tr><td>Depth</td><td>1.0 (0.4)</td><td>-0.2 to 2.2</td><td>0.532 (0.153)</td><td>0.17-0.997</td><td>1.0 (0.7)</td><td>0-4.0</td></tr></tbody></table>

scenario of low speciation and low to intermediate extinction (see the Discussion regarding other potential processes that may explain this unexpected high MPD and low VPD). The fourth phyloregion included the deepest depth strata across all latitudes. It was characterized by low MNTD and PD, but with high MPD and VPD (in black, Fig. 4). This group aligns with a scenario of intermediate to high recent speciation and high past extinction (Fig. 1, Supplementary material Appendix 1 Fig. A5), suggesting that the lack of energy might not preclude recent speciation, and that deep water might not always have been stable. These phyloregions were robust to 1) the removal of PD (Supplementary material Appendix 2 Fig. A19) and 2) the inclusion of VNTD when considering smaller-scale sampling units (Supplementary material Appendix 2 Fig. A20, A21). Results obtained using hierarchical agglomerative clustering generally mirrored the results found using the k-means approach, although there were some minor discrepancies for cells occurring at the interfaces between groups (Supplementary material Appendix 2 Fig. A18–A21).

## Discussion

We have revealed new broadscale phylogenetic diversity patterns for marine ray-finned fishes. Our study indicates that both energy and environmental stability are important drivers generating phylogenetic structures in fish communities along the depth gradient, where energy and stability

are inversely related. Although there are some discrepancies between the observed and expected relationships (Fig. 1) for some phylogenetic diversity metrics, overall our results suggest that energy-rich shallow waters are an engine (i.e. ‘cradle’) of diversity for percomorphs, and that recent speciation for ancient lineages occurring in the deep sea is not driven by energy. Despite potentially facing high extinction early in their evolution, ancient phylogenetic lineages specialized for the deep sea were likely preserved during the Cenozoic by environmental stability. Furthermore, our findings of high MPD and MNTD at intermediate depths suggest that these depths might be a ‘museum’ (or zone of overlap) for the distinct phylogenetic lineages that occur either predominantly in shallow or deep-sea ecosystems. Contrary to previous findings (Brown and Thaje 2014), our study proposes that these intermediate depths might have served as a ‘refuge’ during extreme historical environmental changes such as the Mesozoic anoxic events, and the more frequent and recent climatic events in the shallow ocean scattered through the Cenozoic, rather than being a cradle of speciation. Finally, the phylogenetic community structures observed in these mesotrophic waters suggest some other processes might restrict the co-occurrence of closely related species.

### New phylogenetic diversity patterns versus depth for vertebrates

For the first time, our study documents phylogenetic diversity patterns for the most speciose group of vertebrates
<figure><img src="image_8.png" /></figure>

Figure 4. Phylogenetic bioregions across depth and latitude: KER, Rangitahua (Kermadec Islands); TKI, Three Kings Islands; GBI, Great Barrier Island; Wis, White Island; KKA, Kaikōura; OTA, Otago Peninsula; and AUC, Auckland Islands. (a) The best groupings (k=4 clusters) across depth and latitude of the multivariate phylogenetic diversity measures including phylogenetic diversity (PD), mean pairwise distance (MPD), variance of the pairwise distance (VPD) and mean nearest taxonomic distance (MNTD) based on k-means approach realized over the 100 trees on the pooled 47 cells dataset. The pies represent the proportion of trees supporting the samples in the four clusters. (b) Two dimensional PCA ordination plot showing the four retained clusters using phylogenetic diversity measures extracted from the MCC tree; the convex hull polygons show the envelope of the groups (see the two-dimensional PCA ordination plot showing the four clusters based on the 100 trees, in Supplementary material Appendix 2 Fig. A18).

(Actinopterygii) across broad-scale depth and latitude gradients using a molecular phylogeny. In this study, phylogenetic diversity (PD) matched species richness, as expected (Davies and Buckley 2011, Fritz and Rahbek 2012). Both generally decreased with depth (this decrease was more pronounced at low latitudes), but peaked at intermediate depths for some high-latitude locations, consistent with previous work (Leathwick et al. 2006, Ramirez-Llodra et al. 2010, Brown and Thatje 2014, Costello and Chaudhary 2017). Such an interaction between depth and latitude is likely the reason a cubic relationship (including a secondary peak at intermediate depths) was found to be the best model when depth alone was considered as an explanatory variable (Fig. 2a). The increase of VPD with depth also confirmed previous findings for fishes, supporting the co-occurrence of highly divergent fish lineages likely specialized for the deep sea that have subsequently undergone recent speciation, creating clusters of related species (Zintzen et al. 2011, 2012, Rabosky et al. 2018). In contrast, MPD and MNTD presented new and unexpected patterns with depth for vertebrates; MPD was maximal at 900 m while MNTD showed a broader unimodal pattern suggesting low levels of close relatives at intermediate depths (300–700 m), a pattern also observed for brittlestars (O’Hara et al. 2019). Finally, an increase in VNTD with depth observed only at the smaller spatial scale may be attributed to consistent co-occurrence of certain lineages (Anguilliformes, Notacanthiformes, Aulopiformes), but high

species turnover of close relatives within these lineages. Such fine spatial partitioning among close relatives may be driven by small-scale niche differentiation or biotic interactions such as competition (Cavender-Bares et al. 2009).

Clustering of sampled cells according to phylo-α measures revealed four main regions having similar phylogenetic structures. Phyloregions were primarily stratified by depth and secondarily by latitude (Fig. 4a). Shallow depth strata (50, 100 m) had the highest turnover of phyloregions across latitude (red, yellow, blue and black in Fig. 4a), while at increasing depths phyloregions become more homogenous across latitude (blue and black in Fig. 4a).

## Eco-evolutionary hypotheses behind phylogenetic patterns

Our results generally support our third hypothesis (Fig. 1eii): both energy/temperature and long-term stability act in concert as drivers of biodiversity, and energy/temperature is most likely not only influencing extinction but also speciation, at least for the percomorphs. However, the complex relationships of the phylogenetic diversity patterns along depth and latitude that deviate from the conceptual models suggest that this hypothesis is too naïve to fully explain the observed patterns and other factors undoubtedly play a role. Discrepancies between the observed and expected patterns may reflect uncertainty (species composition, phylogenetic

<figure><img src="image_9.png" /></figure>
relationships) inherent in observational studies and/or oversimplification of our theoretical models. For instance, all species are assumed to have the same physiological limits in our conceptual framework, but due to phylogenetic inertia some lineages may be unable to adapt to the environmental conditions of the deep sea and may be filtered along spatial gradients, leading to highly complex and unanticipated patterns of phylogenetic diversity. Furthermore, our conceptual framework considered only variation along the depth gradient, whereas our results (Fig. 2) indicate that patterns of phylogenetic diversity also vary with latitude. Deviations from our simple conceptual models point the way towards novel hypotheses regarding the distribution of phylogenetic diversity at large spatial scales.

**_The shallow ocean is a cradle of diversity for percomorph fishes_**

Our results (namely, that percomorphs were dominant in shallow waters with high SR and PD, and low MPD, MNTD and VPD) are consistent with the known prevalence of percomorph fishes in coastal and warm tropical shallow waters (e.g. central Indo-Pacific, Price et al. 2015, Miller et al. 2018) and the general increase in coastal marine diversity in warm marine regions (Tittensor et al. 2010). Shallow and low-latitude regions are thought to have high speciation rates caused by: 1) warm, energy-rich waters that may speed up the evolutionary process and maintain viable population sizes (reducing extinction); 2) increased habitat diversity providing opportunities for diversification (Alfaro et al. 2007, Price et al. 2015, Siqueira et al. 2016); and 3) sea-level variations during glacial and inter-glacial periods (e.g. eustatism during the Pleistocene) that can restrict gene flow among populations to create incipient species (Cowman 2014, Ludt and Rocha 2015, Leprieur et al. 2016a). Recent studies on fishes tend to reject the hypothesis that elevated speciation rates are the main driver of species richness in shallow and tropical waters (Price et al. 2015, Miller et al. 2018, Rabosky et al. 2018). Instead, it is suggested that the older age of tropical ecosystems (e.g. central Indo-Pacific) and their long-term climatic stability (34 My) has allowed more time for speciation and colonization, along with decreased extinction risk (Pellissier et al. 2014, Leprieur et al. 2016b, Miller et al. 2018). However, prior studies have not considered depth gradients (but see Rabosky et al. 2018, O'Hara et al. 2019), where there is a sharper decrease in energy with depth than is found with increasing latitude. In our study, we cannot rule out the simultaneous contributions of habitat diversity, instability and fragmentation as potential drivers of speciation in shallow seas. Nonetheless, the rapid decline of SR, PD and the increasing time lag between speciation events with increasing depth (except at the deepest depth we sampled), suggest that energy might also be an important driver of diversity along the depth axis of marine systems (Woolley et al. 2016).

Disentangling the relative contribution of speciation and extinction in the diversification process is difficult without an abundant fossil record (Rabosky 2010). Nevertheless, our findings support the idea that the enhanced diversity found in

shallow waters (50 m) could be the result of increased recent speciation, predominantly in percomorphs (Supplementary material Appendix 2 Fig. A12). Superficially, such a conclusion may seem to contradict recent studies that found the highest rate of speciation in deep and polar waters (Rabosky et al. 2018, O'Hara et al. 2019). However, we observed the lowest values of MNTD at 1200 m indicating that recent speciation is highest in the deep sea, followed by shallow waters at 50 m (Fig. 2d). Recent studies (Rabosky et al. 2018, O'Hara et al. 2019) used a coarser resolution of depth strata (i.e. above and below 200 m) that may have failed to detect raised speciation rates at shallow depths. Nonetheless, the mean DR statistic and relative PD (corrected for species richness) for the shallowest strata in tropical and temperate waters reported by O'Hara et al. (2019) indicated high recent speciation rates that were comparable to deep polar waters, suggesting speciation rates at shallow depths might only be enhanced in conjunction with warm photic surface water.

**_Intermediate ocean depths are a museum of diverse phylogenetic lineages_**

Our findings suggest that intermediate depths of the lower mesophotic zone (500–900 m) are a museum rather than a cradle of diversity, as proposed elsewhere (Brown and Thatje 2014, Gaither et al. 2016, Costello and Chaudhary 2017). In our study, communities were species-poor but represented a wide diversity of phylogenetic lineages characteristic of both shallow (Percomorpha) and deep (Gadiformes, Aulopiformes, Notacanthiformes) regions leading to high MPD values. Such high MPD values can result from either a high or low speciation rate along with a low to moderate extinction rate (Supplementary material Appendix 1 Fig. A5, A7) making it difficult to define the explicit contribution of speciation. Nonetheless, the high MNTD values found in parallel with high MPD values at intermediate depths (Fig. 2d) suggest that low speciation and extinction typify these depths, characteristic of a museum. Intermediate depths might provide a transitional zone where the lower and upper physiological boundaries for shallow and deep-sea specialists, respectively, overlap (Carney 2005, Brown and Thatje 2014). Moreover, these fish communities lacked closely related species – a pattern previously described for brittlestars (O'Hara et al. 2019) – indicating that other factors may limit species diversity at these depths. Intermediate depths are characterized by rapid declines in temperature and particulate organic matter (Levin 2003), but also the presence of oxygen minimum zones (Levin 2003; although unlikely in New Zealand, Chiswell et al. 2015). This may present a challenging environment for fishes, with limited resources and increased competition among close relatives (Priede and Froese 2013). However, over a long period of time, the mesophotic zone might act as an important refuge for diverse lineages; both shallow and deep-sea specialists may survive at intermediate depths when faced with extreme events, such as climatic variation in the shallows (e.g. Pleistocene climatic oscillations; Moffitt et al. 2015) or rarer and ancient (Mesozoic) oeanic anoxic events (OAE) in the deep sea (Priede and Froese
2013, Song et al. 2014, Rogers 2015), enabling subsequent post-event recolonization (Song et al. 2014). As a potential 'museum', intermediate depths form a 'phylogenetic diversity bank' that could be crucial in the future to mitigate anthropogenic threats affecting shallow and/or deep-sea marine ecosystems (Rogers 2015, Danovaro et al. 2017).

## The deep sea as a filter for ancient fish lineages

The deepest fish communities sampled here resembled a selective museum of old phylogenetic lineages that have undergone recent speciation (also suggested by Zintzen et al. 2011, 2012). Given that marine diversity most likely originated in shallow coastal waters (Sallan et al. 2018), several clades of fishes must have invaded the deep sea early in their evolutionary history (Haedrich 1996, Priede and Froese 2013, Sallan et al. 2018) and subsequently radiated (Barros-García et al. 2018). Potentially strong phylogenetic inertia for the adaptive traits required to cope with food scarcity, extreme pressure and low temperatures at depth may explain the low prevalence of recently derived, shallow-origin fish lineages in the deep sea (Priede 2017). Differential filtering of phylogenetic lineages may explain the increase of VPD and deviations from the relationships expected in our conceptual framework (Fig. 1), where all species were considered equal. The phylogenetic structure of deep-sea communities furthermore suggests that ancient actinopterygian lineages might have suffered extinction in the past (i.e. long internal branches). It is unclear whether old deep-sea lineages survived the Mesozoic mass extinction events within in situ deep-sea refugia (Guinot et al. 2013, Thuy et al. 2014) or if they recolonized the deep sea from intermediate depths (Priede and Froese 2013, Priede 2017). Our results showed that, on average, closely related species at 1200 m diverged around 45 Ma ago, which post-dates the last global OAE (OAE2, 96 Ma ago, Takashima et al. 2006) and one of the most recent regional OAE (56 Ma ago, Dickson et al. 2014), after which the deep-sea environment has remained relatively stable (Priede and Froese 2013).

Phylogenetic structures of fishes at depth suggested recent, and potentially high, rates of speciation. This pattern supports other recent findings of high speciation rates in deep-sea environments or cold polar regions for fishes, ophiuroids and octocorals (Chrysogorgiidae) (Pante et al. 2012, Bribiesca-Contreras et al. 2017, Rabosky et al. 2018, O'Hara et al. 2019), and in energy-poor environments of terrestrial systems (Weir and Schlüter 2007, Schlüter and Pennell 2017, Quintero and Jetz 2018). The increase of recent speciation rates in depauperate, cold and polar regions strongly suggests that energy might not be the key driver of speciation in marine systems (Rabosky et al. 2018, O'Hara et al. 2019). The instability of polar seas during the Cenozoic, including several cooling phases that formed ice-caps (41–34 Ma and 14 Ma, Crame 2018), may have driven extinctions, leaving vacant niches that fostered recent speciation (Crame 2018, O'Hara et al. 2019). In temperate deep-sea environments, perhaps long-term stability, coupled with localised physical seafloor disturbances which create spatial heterogeneity (Levin et al. 2001) increased the potential for specialization

and speciation. These recently formed hypotheses demand further examination; the depth gradient may yet prove to be highly useful for decoupling the relative roles of environmental stability versus energy.

## Conclusion

Our development of a simple conceptual framework, and associated description of multiple phylogenetic diversity patterns across broad depth and latitudinal gradients for marine ray-finned fishes (Actinopterygii) offers new insights regarding potential links between energy and long-term stability with the evolutionary mechanisms driving biodiversity. Our findings support the hypothesis that high energy is not a strict requirement for high speciation rates (Rabosky et al. 2018). Furthermore, long-term environmental stability might be crucial to preserve phylogenetically distinct lineages of fishes and thus buffer their extinction (Mittelbach et al. 2007). Last, we found that mean phylogenetic distances among species peaked at intermediate depths, despite their being species poor. This was primarily due to the co-occurrence of shallow and deep specialists, suggesting that intermediate depths might have provided a refuge against frequent climatic instability in shallow waters and rarer deep-sea extinction events (Song et al. 2014).

From the complexity of the phylogenetic diversity patterns we observed, it is apparent that several other abiotic and biotic factors play a role in shaping phylogenetic structures of fish communities, not only through their effects on speciation and extinction processes, but also on dispersal. Scenarios and hypotheses presented here neglect the crucial role of dispersal (Miller et al. 2018), which might be investigated through in-depth documentation of patterns in phylogenetic beta-diversity (Graham and Fine 2008). A better understanding of eco-evolutionary mechanisms shaping distributions of fishes will require integrated studies of taxonomic, phylogenetic and functional diversity, as well as emerging mechanistic models such as the spatial diversification model of lineages through time (SPLIT, Leprieur et al. 2016b, Descombes et al. 2018, Pontarp et al. 2019). Such methods model speciation, extinction and dispersal in spatially-explicit and temporally dynamic environments (e.g. tectonic movements, Leprieur et al. 2016b, Descombes et al. 2018). However, disentangling the contribution of energy and environmental stability along continuous gradients will also require consideration of the evolvability and heritability of species' niche breadths through time (Brayard et al. 2005), allowing for phylogenetic niche conservatism, patrhopatric/parapatric speciation and selective extinction. In addition, to test for the potential indirect contributions of energy, such mechanistic models could be expanded to include diversity-dependent speciation models driven by ecological opportunities (Etienne et al. 2019), indexed on the amount of available energy. We consider that our conceptual framework and quantitative tests provide a first step toward the development of these more elaborate approaches. Integrating empirical and theoretical models in this way, we can gain a
greater understanding of biodiversity determinants for the largest habitat on earth, which will be crucial to achieve better management of ocean ecosystems under increasing future anthropogenic threats (Rogers 2015, Danovaro et al. 2017).

## Data availability statement

All the new sequences have been deposited in GenBank (accession numbers MN915253-MN915266 and MN918523-MN918528) and corresponding metadata is deposited in the Genomics Observatory Metadatabase (GEOME; Expedition code: ActinopterygiiNZ2_CO1_lligins, GUID: ; Expedition code: ActinopterygiiNZ_RAG1_lligins, GUID:). All DNA sequence metadata, alignments and phylogenetic trees are available from the Dryad Digital Repository: <https://doi.org/10.5061/dryad.0rxwdbrw1> (Eme et al. 2020).

**Acknowledgements** – We thank the contribution of the New Zealand eScience Infrastructure (NeSI; <www.nesi.org.nz>) and the Simurg cluster (<www.massey.ac.nz/-theochem/simurg/doku.php>) owned by the New Zealand Inst. of Advanced Study (NZIAS) and the Centre for Theoretical Chemistry and Physics (CTCP), for their high performance computing facilities. We thank Carl Struthers, Jeremy Barker, Salme Kortet and Andrew Stewart (Te Papa) for access to tissue samples, and Georgia Breckell and Vanessa Arranz-Martinez (SNCS) for extraction and PCR of additional taxa. We thank Te Papa's Acquisition Program for assisting the development of an extensive fish tissue collection over the last 15 years, which enabled this research. We thank the editor and two anonymous referees for their valuable comments on an earlier draft of this manuscript.

**Funding** – This work was supported by the Royal Society of New Zealand, Marsden Grant 15-MAU-132. DE and EMVM were supported by the Royal Society of New Zealand, Marsden Grant 15-MAU-132. MJA was supported by a James Cook Fellowship from the Royal Society of New Zealand (Grant JCF-15-MAU-006). LL was supported by a Rutherford Foundation Postdoctoral Fellowship (Grant RFT-15-MAU-003-PD). CDR was supported (in part) by the NZ National Inst. of Water and Atmospheric Research (NIWA) Core Funded Coasts and Oceans Programme 2: biological resources subcontract for fundamental knowledge of marine fish biodiversity.

**Conflicts of interest** – None declare.

**Author contributions** – DE, MJA, LL designed the study, LL obtained the DNA sequences, DE built the phylogenetic trees, performed all the statistical analyses and lead the writing. LL took a lead editing role. All authors contributed to interpretation of the results and writing.

## References

Alfaro, M. E. et al. 2007. Do reefs drive diversification in marine Teleosts? Evidence from the pufferfish and their allies (order Tetraodontiformes). – Evolution 61: 2104–2126.

Allen, A. P. et al. 2006. Kinetic effects of temperature on rates of genetic divergence and speciation. – Proc. Natl Acad. Sci. USA 103: 9130–9135.

Barros-García, D. et al. 2018. Phylogenetic analysis shows the general diversification pattern of deep-sea notacanthiforms (Teleostei: Elopomorpha). – Mol. Phylogenet. Evol. 124: 192–198.

Bouckaert, R. et al. 2014. BEAST 2: a software platform for Bayesian evolutionary analysis. – PLoS Comput. Biol. 10: 1–6.

Brayard, A. et al. 2005. Latitudinal gradient of taxonomic richness: combined outcome of temperature and geographic middomains effects? – J. Zool. Syst. Evol. Res. 43: 178–188.

Bribiesca-Contreras, G. et al. 2017. The importance of offshore origination revealed through ophiuroid phylogenomics. – Proc. R. Soc. B 284: 20170160.

Brown, A. and Thatje, S. 2014. Explaining bathymetric diversity patterns in marine benthic invertebrates and demersal fishes: physiological contributions to adaptation of life at depth. – Biol. Rev. 89: 406–426.

Campbell, N. et al. 2011. Species richness, taxonomic diversity and taxonomic distinctness of the deep-water demersal fish community on the Northeast Atlantic continental slope (ICES Subdivision VIa). – ICES J. Mar. Sci. 68: 365–376.

Carney, R. S. 2005. Zonation of the deep biota on continental margins. – Oceanogr. Mar. Biol. Annu. Rev. 43: 211–278.

Cavender-Bares, J. et al. 2009. The merging of community ecology and phylogenetic biology. – Ecol. Lett. 12: 693–715.

Chiswell, S. M. et al. 2015. Physical oceanography of the deep seas around New Zealand: a review. – N. Z. J. Mar. Freshwater Res. 49: 286–317.

Clarke, K. R. and Warwick, R. M. 1998. A taxonomic distinctness index and its statistical properties. – J. Appl. Ecol. 35: 523–531.

Clarke, K. R. and Warwick, R. M. 2001. A further biodiversity index applicable to species lists: variation in taxonomic distinctness. – Mar. Ecol. Prog. Ser. 216: 265–278.

Condamine, F. L. et al. 2013. Macroevolutionary perspectives to environmental change. – Ecol. Lett. 16: 72–85.

Condamine, F. L. et al. 2019. Assessing the causes of diversification slowdowns: temperature-dependent and diversity-dependent models receive equivalent support. – Ecol. Lett. 22: 1900–1912.

Costello, M. J. and Chaudhary, C. 2017. Marine biodiversity, biogeography, deep-sea gradients and conservation. – Curr. Biol. 27: R511–R527.

Cowman, P. F. 2014. Historical factors that have shaped the evolution of tropical reef fishes: a review of phylogenies, biogeography and remaining questions. – Front. Genet. 5: 1–15.

Crame, J. A. 2018. Key stages in the evolution of the Antarctic marine fauna. – J. Biogeogr. 45: 986–994.

Danovaro, R. et al. 2017. An ecosystem-based deep-ocean strategy. – Science 355: 452–454.

Davies, R. G. et al. 2007. Environmental predictors of global parrot (Aves: Psittaciformes) species richness and phylogenetic diversity. – Global Ecol. Biogeogr. 16: 220–233.

Davies, T. J. and Buckley, L. B. 2011. Phylogenetic diversity as a window into the evolutionary and biogeographic histories of present-day richness gradients for mammals. – Phil. Trans. R. Soc. B 366: 2414–2425.

Descombes, P. et al. 2018. Linking species diversification to palaeoenvironmental changes: a process-based modelling approach. – Global Ecol. Biogeogr. 27: 233–244.

Dickson, A. J. et al. 2014. The spread of marine anoxia on the northern tethys margin during the Paleocene–Eocene thermal maximum. – Paleogeography 29: 471–488.

Eme, D. et al. 2019. An integrated pathway for building regional phylogenies for ecological studies. – Global Ecol. Biogeogr. 28: 1899–1911.

Eme, D. et al. 2020. Data from: Phylogenetic measures reveal evolutionary drivers of biodiversity along a depth gradient. – Dryad Digital Repository, <https://doi.org/10.5061/dryad.0rxwdbrw1>.
Etienne, R. S. et al. 2019. A minimal model for the latitudinal diversity gradient suggests a dominant role for ecological limits. – Am. Nat. 194: E122–E133.

Evans, K. L. et al. 2005. Species-energy relationships at the macroecological scale: a review of the mechanisms. – Biol. Rev. 80: 1–25.

Faith, D. P. 1992. Conservation evaluation and phylogenetic diversity. – Biol. Conserv. 61: 1–10.

Fiedler, P. C. 2010. Comparison of objective descriptions of the thermocline. – Limnol. Oceanogr. Methods 8: 313–325.

Fine, P. V. A. 2015. Ecological and evolutionary drivers of geographic variation in species diversity. – Annu. Rev. Ecol. Evol. Syst. 46: 369–392.

Fritz, S. A. and Rahbek, C. 2012. Global patterns of amphibian phylogenetic diversity. – J. Biogeogr. 39: 1373–1382.

Gaither, M. R. et al. 2016. Depth as a driver of evolution in the deep sea: insights from grenadiers (Gadiformes: Macrouridae) of the genus Coryphaenoides. – Mol. Phylogenet. Evol. 104: 73–82.

Gajdzik, L. et al. 2019. Diversification and functional evolution of reef fish feeding guilds. – Ecol. Lett. 22: 572–582.

Gillman, L. N. et al. 2009. Latitude, elevation and the tempo of molecular evolution in mammals. – Proc. R. Soc. B 276: 3353–3359.

Graham, C. H. and Fine, P. V. A 2008. Phylogenetic beta diversity: linking ecological and evolutionary processes across space in time. – Ecol. Lett. 11: 1265–1277.

Guinot, G. et al. 2013. Cretaceous stem chondrichthyans survived the end-Permian mass extinction. – Nat. Commun. 4: 1–8.

Haedrich, R. L. 1996. Deep-water fishes: evolution and adaptation in the earth's largest living spaces. – J. Fish Biol. 49: 40–53.

Ivany, L. C. et al. 2000. Cooler winters as a possible cause of mass extinctions at the Eocene/Oligocene boundary. – Nature 407: 887–890.

Jablonski, D. et al. 2006. Out of the tropics: evolutionary dynamics of the latitudinal diversity gradient. – Science 314: 102–106.

Jansson, R. and Dynesius, M. 2002. The fate of clades in a world of recurrent climatic change: Milankovitch oscillations and evolution. – Annu. Rev. Ecol. Syst. 33: 741–777.

Jenkyns, H. C. 2010. Geochemistry of oceanic anoxic events. – Geochem. Geophys. Geosyst. 11: 1–30.

Jetz, W. and Fine, P. V. A. 2012. Global gradients in vertebrate diversity predicted by historical area-productivity dynamics and contemporary environment. – PLoS Biol. 10: e1001292.

Jombart, T. 2008. adegenet: a R package for the multivariate analysis of genetic markers. – Bioinformatics 24: 1403–1405.

Jombart, T. et al. 2010. Discriminant analysis of principal components: a new method for the analysis of genetically structured populations. – BMC Genet. 11: 94.

Kozak, K. H. and Wiens, J. J. 2012. Phylogeny, ecology and the origins of climate-richness relationships. – Ecology 93: 167–181.

Leathwick, J. R. et al. 2006. Variation in demersal fish species richness in the oceans surrounding New Zealand: an analysis using boosted regression trees. – Mar. Ecol. Prog. Ser. 321: 267–281.

Leprieur, F. et al. 2016a. Historical and contemporary determinants of global phylogenetic structure in tropical reef fish faunas. – Ecography 39: 825–835.

Leprieur, F. et al. 2016b. Plate tectonics drive tropical reef biodiversity dynamics. – Nat. Commun. 7: 11461.

Levin, L. A. 2003. Oxygen minimum zone benthos: adaptation and community response to hypoxia. – Oceanogr. Mar. Biol. Annu. Rev. 41: 1–45.

Levin, L. A. et al. 2001. Environmental influences on regional deep-sea species diversity. – Annu. Rev. Ecol. Evol. Syst. 32: 51–93.

Ludt, W. B. and Rocha, L. A. 2015. Shifting seas: the impacts of Pleistocene sea-level fluctuations on the evolution of tropical marine taxa. – J. Biogeogr. 42: 25–38.

Mazel, F. et al. 2016. Influence of tree shape and evolutionary timescale on phylogenetic diversity metrics. – Ecography 39: 913–920.

McClain, C. R. and Schlacher, T. A. 2015. On some hypotheses of diversity of animal life at great depths on the sea floor. – Mar. Ecol. 36: 849–872.

McClymont, E. L. et al. 2016. Pliocene-Pleistocene evolution of sea surface and intermediate water temperatures from the southwest Pacific. – Paleooceanography 31: 895–913.

Miller, E. C. et al. 2018. Explaining the ocean's richest biodiversity hotspot and global patterns of fish diversity. – Proc. R. Soc. B 285: 20181314.

Mittelbach, G. G. et al. 2007. Evolution and the latitudinal diversity gradient: speciation, extinction and biogeography. – Ecol. Lett. 10: 315–331.

Moffitt, S. E. et al. 2015. Response of seafloor ecosystems to abrupt global climate change. – Proc. Natl Acad. Sci. USA 112: 4684–4689.

Mouquet, N. et al. 2012. Ecophylogenetics: advances and perspectives. – Biol. Rev. 87: 769–785.

O'Hara, T. D. O. et al. 2019. Contrasting processes drive ophiuroid phylodiversity across shallow and deep seafloors. – Nature 565: 636–639.

Pante, E. et al. 2012. Deep-sea origin and in-situ diversification of Chrysogorgiid octocorals. – PLoS One 7: e38357.

Pavoine, S. and Bonsall, M. B. 2011. Measuring biodiversity to explain community assembly: a unified approach. – Biol. Rev. 86: 792–812.

Pellissier, L. et al. 2014. Quaternary coral reef refugia preserved fish diversity. – Science 344: 1016–1019.

Pillans, B. et al. 1998. A review of the Milankovitch climatic beat: template for Plio-Pleistocene sea-level changes and sequence stratigraphy. – Sediment. Geol. 122: 5–21.

Pontarp, M. et al. 2019. The latitudinal diversity gradient: novel understanding through mechanistic eco-evolutionary models. – Trends Ecol. Evol. 34: 211–223.

Price, S. A. et al. 2015. Phylogenetic insights into the history and diversification of fishes on reefs. – Coral Reefs 34: 997–1009.

Priede, I. G. 2017. Deep-sea fishes: biology, diversity, ecology and fisheries. – Cambridge Univ. Press.

Priede, I. G. and Froese, R. 2013. Colonization of the deep sea by fishes. – J. Fish Biol. 83: 1528–1550.

Pyron, R. A. and Wiens, J. J. 2013. Large-scale phylogenetic analyses reveal the causes of high tropical amphibian diversity. – Proc. R. Soc. B 280: 20131622.

Quintero, I. and Jetz, W. 2018. Global elevational diversity and diversification of birds. – Nature 555: 246–250.

Rabosky, D. L. 2010. Extinction rates should not be estimated from molecular phylogenies. – Evolution 64: 1816–1824.

Rabosky, D. L. et al. 2018. An inverse latitudinal gradient in speciation rate for marine fishes. – Nature 559: 392–395.

Ramirez-Llodra, E. et al. 2010. Deep, diverse and definitely different: unique attributes of the world's largest ecosystem. – Bio-geosciences 7: 2851–2899.

Rangel, T. F. et al. 2015. Phylogenetic uncertainty revisited: implications for ecological analyses. – Evolution 69: 1301–1312.

Rogers, A. D. 2015. Environmental change in the deep Ocean. – Annu. Rev. Environ. Resour. 40: 1–38.
Rohde, K. 1992. Latitudinal gradients in species diversity: the search for the primary cause. – Oikos 65: 514–527.

Rolland, J. and Salamin, N. 2016. Niche width impacts vertebrate diversification. – *Global Ecol. Biogeogr.* 25: 1252–1263.

Rolland, J. et al. 2014. Faster speciation and reduced extinction in the tropics contribute to the mammalian latitudinal diversity gradient. – *PLoS Biol.* 12: e1001775.

Salisbury, C. L. et al. 2012. The latitudinal gradient in dispersal constraints: ecological specialisation drives diversification in tropical birds. – *Ecol. Lett.* 15: 847–855.

Sallan, L. et al. 2018. The nearshore cradle of early vertebrate diversification. – *Science* 464: 460–464.

Sanmartín, I. and Meseguer, A. S. 2016. Extinction in phylogenetics and biogeography: from timetrees to patterns of biotic assemblage. – *Front. Genet.* 7: 1–17.

Schlüter, D. and Pennell, M. W. 2017. Speciation gradients and the distribution of biodiversity. – *Nature* 546: 48–55.

Siqueira, A. C. et al. 2016. Evolutionary processes underlying latitudinal differences in reef fish biodiversity. – *Global Ecol. Biogeogr.* 25: 1466–1476.

Song, H. et al. 2014. Anoxia/high temperature double whammy during the Permian-Triassic marine crisis and its aftermath. – *Sci. Rep.* 4: 4132.

Takashima, R. et al. 2006. Greenhouse world and the Mesozoic Ocean. – *Oceanography* 19: 64–74.

Thuy, B. et al. 2014. First glimpse into Lower Jurassic deep-sea biodiversity: in situ diversification and resilience against extinction. – *Proc. R. Soc. B* 281: 20132624.

Tittensor, D. P. et al. 2010. Global patterns and predictors of marine biodiversity across taxa. – *Nature* 466: 4–9.

Supplementary material (available online as Appendix ecog-04836 at <www.ecography.org/appendix/ecog-04836>). Appendix 1–2.

Tolimieri, N. and Anderson, M. J. 2010. Taxonomic distinctness of demersal fishes of the California Current: moving beyond simple measures of diversity for marine ecosystem-based management. – *PLoS One* 5: 22–24.

Tucker, C. M. et al. 2017. A guide to phylogenetic metrics for conservation, community ecology and macroecology. – *Biol. Rev.* 92: 698–715.

Valentine, J. W. and Jablonski, D. 2015. A twofold role for global energy gradients in marine biodiversity trends. – *J. Biogeogr.* 42: 997–1005.

Vamosi, S. M. et al. 2009. Emerging patterns in the comparative analysis of phylogenetic community structure. – *Mol. Ecol.* 18: 572–592.

Webb, C. 2000. Exploring the phylogenetic structure of ecological communities: an example for rain forest trees. – *Am. Nat.* 156: 145–155.

Weir, J. T. and Schlüter, D. 2007. The latitudinal gradient in recent speciation and extinction rates of birds and mammals. – *Science* 315: 1574–1576.

Woolley, S. N. C. et al. 2016. Deep-sea diversity patterns are shaped by energy availability. – *Nature* 533: 7603–7620.

Zintzen, V. et al. 2011. Increasing variation in taxonomic distinctness reveals clusters of specialists in the deep sea. – *Ecography* 34: 306–317.

Zintzen, V. et al. 2012. Diversity and composition of demersal fishes along a depth gradient assessed by baited remote underwater. – *PLoS One* 7: e48522.

Zintzen, V. et al. 2017. Effects of latitude and depth on the beta diversity of New Zealand fish communities. – *Sci. Rep.* 7: 1–10.