Patterns and Variation in Benthic Biodiversity in a Large Marine Ecosystem

While there is a persistent inverse relationship between latitude and species diversity across many taxa and ecosystems, deviations from this norm offer an opportunity to understand the conditions that contribute to large-scale diversity patterns. Marine systems, in particular, provide such an opportunity, as marine diversity does not always follow a strict latitudinal gradient, perhaps because several hypothesized drivers of the latitudinal diversity gradient are uncorrelated in marine systems. We used a large scale public monitoring dataset collected over an eight year period to examine benthic marine faunal biodiversity patterns for the continental shelf (55–183 m depth) and slope habitats (184–1280 m depth) off the US West Coast (47°20′N—32°40′N). We specifically asked whether marine biodiversity followed a strict latitudinal gradient, and if these latitudinal patterns varied across depth, in different benthic substrates, and over ecological time scales. Further, we subdivided our study area into three smaller regions to test whether coast-wide patterns of biodiversity held at regional scales, where local oceanographic processes tend to influence community structure and function. Overall, we found complex patterns of biodiversity on both the coast-wide and regional scales that differed by taxonomic group. Importantly, marine biodiversity was not always highest at low latitudes. We found that latitude, depth, substrate, and year were all important descriptors of fish and invertebrate diversity. Invertebrate richness and taxonomic diversity were highest at high latitudes and in deeper waters. Fish richness also increased with latitude, but exhibited a hump-shaped relationship with depth, increasing with depth up to the continental shelf break, ~200 m depth, and then decreasing in deeper waters. We found relationships between fish taxonomic and functional diversity and latitude, depth, substrate, and time at the regional scale, but not at the coast-wide scale, suggesting that coast-wide patterns can obscure important correlates at smaller scales. Our study provides insight into complex diversity patterns of the deep water soft substrate benthic ecosystems off the US West Coast.


Introduction
The existence of consistent, global patterns in biodiversity has long intrigued ecologists, as general patterns hint at common underlying mechanisms in the search for universal laws in ecology. A robust pattern observed across taxa and over large spatial scales is the inverse relationship between species biodiversity and latitude, with increasing diversity towards the equator (latitudinal diversity gradient; [1][2][3]). This gradient has largely been observed in terrestrial systems, with important contributions from the study of extinct marine taxa [4][5][6] and increasing interest in patterns of extant marine diversity [7][8][9][10]. The latitudinal diversity gradient is hypothesized to be caused by many different ecological and evolutionary processes that interact over space and time (see [11][12][13] for recent summaries). Given the challenge of inferring process from pattern, disentangling these hypothesized drivers remains one of the great challenges in macroecology [14]. Deviations from the latitudinal diversity gradient [15] provide an opportunity to examine what conditions "break the rule". Marine systems, in particular, provide such an opportunity, given that marine diversity does not always follow a strict latitudinal gradient and two of the main hypothesized drivers of the latitudinal diversity gradient (i.e., high diversity in tropics results from high temperature and high productivity) are uncorrelated in some marine environments [16].
Notably, marine diversity often shows inconsistent patterns across latitude. Although a landmark analysis of~200 studies found a strong latitudinal gradient in marine biodiversity, the pattern differed in strength among habitat types, organisms, and spatial scales [17]. Other studies have found inconsistent relationships with latitude when limited to a few taxa [18,19] or spatial scales smaller than an entire continent [20,21]. Such studies have found no relationship between marine diversity and latitude, or even a reverse relationship, i.e. higher diversity at higher latitudes, often the patterns are attributed to depth or other physical gradients (e.g., [22,23]). By explicitly examining and accounting for factors that contribute to deviations in marine latitudinal diversity gradients, we can gain insight into the processes that drive biodiversity patterns in marine systems.
Patterns of marine diversity might not follow a strict latitudinal gradient because the drivers of marine diversity themselves are not always correlated with latitude. In terrestrial systems, latitude and elevation are proxies for gradients in temperature, divergence time and evolutionary stability that drive the gradient in diversity [6,[24][25][26]. In marine systems, depth acts as a similar proxy to latitude, as increasing depth co-varies with decreasing light availability, decreasing temperature, increasing pressure, decreasing productivity and greater stability in salinity and temperature [27]. Further, depth may be a more influential driver at smaller spatial scales, just as elevational gradients strongly influence terrestrial diversity [8,28]. This pattern has been observed in groundfish assemblages in the North Pacific, where depth is a stronger predictor of diversity than latitude [29,30]. At greater depths, benthic organisms are often tightly associated with particular substrates, such as sand, mud, and rock, which can influence their broader distributional patterns. In particular, benthic marine faunal diversity is often highest around rocky outcroppings, where habitat complexity is high, leading to more physical space for organisms [31,32]. Thus, in marine systems, if substrate type varies latitudinally, this may act as another important driver of marine biodiversity.
Deviations from the latitudinal gradient in diversity may also result from the implicit timeaveraged nature of most biodiversity studies. Macroecological patterns of biodiversity have largely been evaluated using long-term time series or data from a snapshot in time [33], but diversity is unlikely to be spatially consistent from year-to-year, particularly in highly dynamic systems that include species with broad-scale dispersal [34]. In marine systems, the environmental and climatological factors that play an important role in structuring communities vary inter-annually (e.g. El Niño Southern Oscillation and the Pacific Decadal Oscillation; [35][36][37]). Thus, explicitly evaluating inter-annual variation in diversity on ecological timescales could explain deviations from traditional latitudinal gradients and uncover the underlying processes that shape marine biodiversity.
Given the difficulty in disentangling these underlying drivers of marine biodiversity, ecologists have turned to comparisons of multiple metrics of biodiversity to help tease apart a priori hypotheses about the different drivers [38][39][40][41]. While traditional metrics of biodiversity focus on taxonomic distinctions among species, analysis of functional diversity is a way to understand the mechanisms underlying the distribution of organisms by relating functional traits to environmental drivers [41,42]. In particular, by describing the variation in species functional traits, functional diversity captures patterns of biodiversity that species richness and abundance do not. Examination of broad-scale latitudinal patterns in functional diversity and the relationships among different metrics of biodiversity offer a way forward towards quantifying biodiversity gradients and understanding their drivers [15,41,[43][44][45][46][47].
We tested for the presence of a latitudinal diversity gradient in temperate marine communities with a large-scale, spatio-temporal analysis of patterns in marine benthic faunal diversity. Because patterns of marine biodiversity across latitude are sensitive to a myriad of environmental drivers, temporal variability, and the biodiversity metric used, these factors must be accounted for in any analysis of latitudinal gradients. Thus, using a public scientific monitoring dataset spanning the California Current Large Marine Ecosystem off the US West Coast, we asked, 1) Is there a latitudinal gradient in marine faunal diversity? 2) What is the role of depth and substrate in structuring such a latitudinal gradient? and 3) Are these patterns dependent on space and time? We hypothesized that marine benthic diversity will be higher at lower latitudes [17], in deeper waters [16], and around rocky substrates [31,32]. Additionally, we hypothesized that these patterns will be sensitive to inter-annual variability in the coastal environment, though a priori we did not have an expectation about the direction of the impact on marine benthic diversity patterns. Following Hillebrand [3,17], we expected that the strength of the latitudinal diversity gradient will be stronger at coast-wide scales than at smaller, regional scales, where oceanographic processes tend to influence community structure and function [48][49][50]. To test these hypotheses, we examined spatio-temporal patterns of diversity for benthic fish and invertebrate species across eight years using three diversity metrics (taxonomic richness, taxonomic diversity, and functional diversity). Given the importance of biodiversity to maintaining ecosystem function and the gap in knowledge of biodiversity patterns in the California Current marine system, understanding the variability of species and functional trait distributions over large spatial scales is essential to conserving the structure and function of this and other systems in the face of environmental change.

Study Area
The California Current Large Marine Ecosystem (CCLME) comprises the coastal and offshore region south of British Columbia to the southern tip of Baja California, encompassing 2.2 million km 2 surface area (Fig 1; [51,52]). The CCLME spans temperate to subtropical climatic regions and experiences strong seasonal upwelling, which brings cold, nutrient rich water to the surface during the summer and increases productivity. The benthic habitat includes hard substrate (submarine canyons, ridges, and rocky reefs) within an expansive matrix of unconsolidated sand, mud, and gravel [53]. Loose sediment particles range from silt to boulders. The continental shelf is fairly narrow, so the continental slope (approximately >200m depth) is relatively close to the coast, with coastal proximity to deep slope greatest in the south [52]. High inter-annual variability in temperature and productivity characterizes the CCLME. During the period of data collection, summer bottom temperatures in the study area ranged from 2.91°C-14.7°C.

Species Abundance Data
We used fish and invertebrate data from the West Coast Groundfish Bottom Trawl Survey (WCGBTS) conducted by the Northwest Fisheries Science Center (NWFSC), National Oceanic and Atmospheric Administration (NOAA), US [54]. WCGBTS is a fishery-independent survey led by NOAA scientists and conducted annually in conjunction with chartered commercial fishing vessels to collect and identify species of fish and invertebrate taxa. We examined data for the years 2003-2010, as 2003 was the first year of data collection by the NWFSC and 2010 was the most recent year of data available. Scientific survey trawls are conducted annually from mid-May until late October. The survey follows a stratified-random sampling design based on geographic location (north or south of Point Conception, California) and depth (55-183 m, 184-549 m, or 550-1,280 m). The survey samples both the continental shelf (depth 55-183 m) and the continental slope (depth 184-1,280 m) from Cape Flattery, WA (47°20 0 N) to the US-Mexico border (32°40 0 N). Tows employed a Aberdeen-type net, with a small mesh (3.81cm 2 ) codend liner to retain smaller fish and invertebrates [54]. Trawling occurred over soft to rocky benthic substrate within randomly selected cells (1.5 nautical miles (nm) longitude by 2.0 nm latitude) in the three depth strata. Additional information on the survey methodology can be found in Bradburn et al. [54]. The WCGBTS was designed to provide information on biomass, abundance, length, and age of commercially harvested fish for population assessments, but all individuals captured (regardless of commercial interest) are identified to the lowest taxonomic level possible (S1 Table). Species that are unidentified onboard were labeled, frozen, or preserved in formalin and saved for later identification [54]. Due to time constraints between tows, staff did not enumerate individual organisms; instead, the biomass of each identified taxon was determined on deck. We standardized the biomass of each species to the area swept of the trawl net (kg/ha), a measure of catch per unit effort (CPUE), to control for differences in the sampling area across trawls. We minimized the risk of double counting species by combining weights of species rarely identified to species level into a single category by genus (S1 File). We also eliminated pelagic species to avoid inflating diversity estimates with species incidentally caught during net deployment or retrieval.

Environmental Data
We included three other predictors known to influence diversity that co-vary with latitude: depth, year, and substrate. Depth was recorded at the mid-point of each trawl. Bottom temperature was also recorded for each trawl, but it was strongly collinear with depth. We chose to use depth as a predictor instead of temperature, because depth captures additional dependent factors beyond temperature. Year was included as a predictor not only as a proxy for unmeasured environmental conditions, but also to account for the temporal structure of the data collection across eight years.
We extracted seafloor substrate type from the Surficial Geologic Habitat Map (SGHM) GIS layer for the US West Coast [55]. The SGHM combines data collected from bathymetry, fieldsampled point surveys, acoustic imagery, and sub-bottom analysis to create a mixed-resolution map of the seafloor, interpreted using the highest level of detail that could be justified through each method. Thus, habitat identification represented on the map is accurate, but varies in precision (C. Romsos, pers. communication). We used the SGHM to create a binary variable classifying the sea-floor bottom as either soft (silt, mud, sand, and gravel) or rocky (cobbles to small boulders). We then overlaid trawl survey points onto the habitat layer and sampled the substrate type at the point of the trawl. Each point was assigned the primary habitat type over which the trawl was towed, thus mixed-habitat trawls are not represented in this dataset.

Biodiversity Measures
We examined two taxonomic biodiversity indices, species richness (S), and Shannon diversity (H 0 ) for fish and invertebrates. We calculated H 0 as the sum of proportional biomass/area swept (kg/ha) of all species in a sample: where p i is the proportion of biomass (CPUE) of species i in a trawl [56]. We chose this measure for two reasons: our data were measured as biomass, rather than counts of individuals/species, which precluded the use of many biodiversity indices; and H 0 , in which biomass units have been applied as a measure of abundance [57][58][59][60][61], is a commonly used and readily comparable index of biodiversity [61,62]. In addition, we used a third measure, Rao's quadratic entropy index (Q) to quantify functional diversity [63], or the diversity of traits in a system, using life history and trophic characteristics of the species observed. We only calculated Q measures for fish, because limited data on invertebrate traits were available. For fishes, we acquired information on a total of 18 traits using FishBase (S1 File and S2 Table; [64,65]). We used these traits to calculate Q for each trawl. Q is a measure of the abundance-weighted trait diversity of a sample: where S is the total species richness, d ij is the dissimilarity in traits between species i and species j, and p i and p j are the relative abundances of species i and j. We calculated Q because trait diversity may show patterns that are not otherwise seen using measures of taxonomic richness and diversity [42,66].

Statistical Analysis
We used model-averaging to fit a series of generalized linear models (GLMs) to test our hypotheses regarding fish and invertebrate biodiversity gradients in relation to latitude and environmental co-variates for the five response variables: S and H 0 for fish and invertebrates, and Q for fish species only. This approach is common in the ecological literature for its ability to flexibly account for many environmental co-variates and differences in the distribution of response variables while preserving the interpretability of a linear framework [67]. All GLMs included the same four potential explanatory variables (latitude, year of the trawl, depth, and substrate type) based on previous work done in this system, as well as studies using fishery survey data to describe abundance patterns [29,30,68]. We restricted the candidate set to linear models 1) to retain interpretability and facilitate comparison with previous tests of a latitudinal diversity gradient, 2) because processes hypothesized to cause a latitudinal diversity gradient provide no mechanistic framework for determining a priori the appropriate non-linear functional form, and 3) because we sought to keep our methods consistent across predictors. We then used a model averaging information theoretic approach in conjunction with the GLMs [69,70]. This approach tests biological hypotheses in the form of candidate models and provides a quantitative measure of relative support for competing hypotheses based on a robust statistical framework [71]. Model selection is preferred to traditional null hypothesis testing for observational data, where causative factors cannot be isolated by sampling design or when multiple models have similar levels of support [71][72][73]. Traditional, stepwise model selection routines have been shown to lead to spurious conclusions, as inference is conditional on the single "best" model and variables included [74]. In contrast, Akaike's Information Criterion (AIC)based model averaging incorporates the inherent uncertainty in model selection and thus improves parameter estimation when competing models have similar support from the data, and protects against spurious conclusions regarding important variables in the observed system [75].
We treated each combination of maximum likelihood-estimated parameters in the regression models as competing hypotheses describing the interaction between the explanatory variables (latitude, depth, year and substrate) and S, H 0 , and Q, respectively. We considered the model (out of a possible 2 13 = 8,192 models given by the number of environmental variables and their interactions) with the lowest AICc value the top model and used to calculate ΔAICc for each subsequent model ranked by AICc, where the best model's ΔAICc = 0. We then calculated Akaike weights, w i , for each model [69,76]. We considered all models with a ΔAICc < 2 to be supported by the data and therefore included them, creating a confidence set of top models [69,76]. From the confidence set of models, we calculated model averaged parameter estimates (β i ) as weighted sums of parameter coefficients given by the product of the parameter estimate (β i ) in model i times its Akaike weight, w i [76]. Similarly, we calculated 95% confidence intervals (95% CI) around the model averaged parameters from the confidence set of models [69,76]. We inferred that variables with confidence intervals which did not include zero were related to the biodiversity metric being modeled. We then used the sign of the model averaged parameters to infer the direction of the relationship of the variables correlated with biodiversity in the CCLME for both individual variables as well as interactions among them.
The use of generalized linear models in this framework allowed us to account for differences in the functional form of our response variables. Fish and invertebrate H 0 and Fish Q fit the assumptions of normality, and thus were modeled with a normal distribution. We modeled invertebrate S with a negative binomial distribution, rather than the more commonly used Poisson distribution, because diagnostic plots of the average square residuals to the mean revealed overdispersion [77]. We modeled fish S with Poisson because although the fish S data were estimated to be slightly overdispersed (C = 1.099), a negative binomial model did not provide a better fit to the data than a Poisson model (Likelihood ratio test, χ2 = 1.714, p = 0.1905). We examined histograms of residuals for all models and used variance inflation factors to determine that all predictors were non-collinear. We checked all responses (S, H 0 , and Q) for significant outliers using leverage scores and studentized residuals and found no points with both high leverage and high residuals. We also analyzed the residuals of our best models for spatial autocorrelation using Moran's I. Values fell near zero, indicating that models adequately accounted for spatial autocorrelation without further correlation structures (see S1 File). We calculated R 2 (for H 0 and Q) and McFadden's pseudo-R 2 (for S) for all models. Calculating R 2 wasn't possible for our final averaged models, so we calculated R 2 using a model that included all the (un-averaged) parameters selected for in the final model. Thus, these R 2 and pseudo-R 2 estimates should be interpreted with caution. All analyses were conducted in R, version 2.12.0, with spatial patterns of biodiversity displayed using the Fields and Maps packages, and model selection from the MuMIn package [78,79].

Comparisons to Models Accounting for Latitude
To test whether latitude is an important driver of marine biodiversity, we used the evidence ratio of the AICc-selected best model with and without latitude (Eq 3). The ratio, ρ, provides an estimate of the relative support in the data for two competing models based on the AICc weights [76]: where w i and w j are the AICc weights of models i and j. We compared the top-ranking model that included latitude to the top-ranking model that did not include latitude for all five responses. This metric allowed us to determine, given the data, the factor by which the complex best model is potentially more likely compared to a model without latitude.

Comparisons Across Biogeographic Regions
Previously published research suggests that there may be discrete regions in the California Current system in regard to biodiversity and that the factors affecting biodiversity patterns within regions may vary [80]. To elucidate patterns at the regional scale that may be obscured by our coastwide analysis, we divided the coastline into three biogeographic regions, north of Cape Mendocino, CA (North), between Cape Mendocino and Point Conception, CA (Central), and south of Point Conception (South; Fig 1). We performed the same GLM/model selection routines for these three regions as we did for the entire CCLME dataset. With the addition of these three regional models for each of our three diversity responses to those of the entire system, we carried out 20 total model selection procedures (four total spatial datasets five diversity responses) in this study.

Results
Overall, the trawl survey encountered 233 fish and 310 invertebrate taxa in 171 families across all years, from Washington to California (S1 Table). Maps of richness (S), Shannon diversity (H 0 ), and functional diversity (Q) suggested that spatial patterns across both depth and latitude differ between fish and invertebrate taxa (Figs 2 and 3). The relationships of fish and invertebrate S with latitude and depth over time revealed some large scale patterns and temporal variability (Fig 4, S1 and S2 Figs). We observed highly consistent patterns across years in the relationship between depth and fish S, but more variability among years in latitude and fish S and invertebrate S (Fig 4).

Patterns Across the CCLME
Latitude represented an important driver of all metrics of diversity for fish and invertebrates across the entire spatial region of the CCLME, as all AICc-selected models with latitude were orders of magnitude more likely than those without latitude (all ρ ! 131.35, Table 1). While latitude appeared in all top-ranking models, its relationship with the different metrics of diversity varied and was dependent on depth, year, and substrate. In other words, diversity was not always highest at low latitudes, but instead showed complex patterns in space and time.
Overall, fish S increased with increasing latitude, but this relationship also depended on the effects of depth, year, and substrate ( Table 2, Fig 3 and S1 Fig). For example, the effect of latitude on fish S was stronger in rocky areas than in soft-sediment areas (Table 2). Furthermore, fish S appeared to be more influenced by depth than latitude, with a peak in richness at approximately 200 m and then declined sharply with increasing depth across most latitudes (Figs 3  and 4). Although there were obvious peaks in fish H 0 and Q at certain depth and latitude combinations (Figs 2 and 3), there were no significant predictors for fish H 0 or Q, indicating that these responses did not vary linearly based on our explanatory variables ( Table 2). Fish H 0 and Q did not change across the latitudinal gradient, suggesting that our explanatory variables did not AFFECt these abundance-weighted metrics ( Table 2).
Consistent with the results of fish S, both invertebrate S and H 0 were influenced by a positive interaction among latitude, depth, and year (Table 3), thus invertebrate S and H 0 increased at higher latitudes, but also in deeper waters and through time. The clear positive relationship with depth is evident when invertebrate S and H 0 is plotted against both depth and latitude (Fig 3). For both invertebrate S and H 0 , the effect of latitude was stronger in rocky substrates than soft substrates (Table 3). We did not examine Q for invertebrates due to lack of biological trait data.

Patterns Across Biogeographic Regions
When we divided the spatial extent of the survey into three biogeographic regions, we did not see a consistent weakening of the effect of latitudinal gradients, but instead the smaller-scale models revealed patterns that were obscured at the coast-wide scale (Table 4 and see S3 Table  for full results). For example, fish H 0 , which was not associated with any co-variates at the coast-wide scale, was positively associated with the interaction of latitude, depth, and year in the South region (Table 4 and S3 Table). However, no model co-variates strongly explained the patterns of fish H 0 in the North and Central regions. Likewise, several important predictors for fish Q emerged at the regional scale, in contrast to the coast-wide analysis, which did not find   This was determined using the evidence ratio (ρ), which is calculated by first determining the best model using Akaike's Information Criterion correction (AICc), and then dividing the AICc weight of the model when the term is in the model (w i ) by the AICc weight of the model when the term is removed (w j ). The greater the ρ, the more important latitude is as a predictor in the model. North region, the interaction between latitude, depth, and year was positive (Table 4 and S3  Table). In the Central region, when the substrate was rocky, fish Q decreased with latitude (Table 4 and S3 Table). In the South region, there was a greater effect of latitude at shallower depths, and a greater effect of depth in rocky substrates (Table 4 and S3 Table). Within each biogeographic region, all metrics of diversity were related to latitude, but, as on the whole-coast scale, the metrics did not follow a consistent latitudinal gradient. For example, the effect of the three-way interaction between depth, latitude, and substrate for invertebrate S was positive in the North region, but negative in the Central region (Table 4 and S3 Table). Our analysis revealed greater complexity in the patterns of both fish and invertebrate diversity at the regional scale compared to the coast-wide scale.

Discussion
Our study documents diversity patterns and evaluates important abiotic correlates of benthic fish and invertebrate diversity at the scale of a large marine ecosystem. We show that benthic biodiversity in this system does not conform to a simple latitude-diversity gradient, as we expected. Rather, the effect of latitude depends upon the interplay of depth and substrate, and varies depending on year, taxa, and biogeographic region. We also found that fish and invertebrate diversity displayed contrasting responses to depth, which likely reflects key physiological and behavioral differences. Further, by analyzing the influence of latitude, depth, substrate, and The Akaike Information Criterion correction (AICc) was used to rank models and any model that ranked <2 ΔAICc was averaged to obtain final estimates presented. Bolded model terms have a 95% confidence interval (CI) that did not include zero. Relative importance (RI) refers to the proportion of output models that contained the term before the model estimates were averaged. doi:10.1371/journal.pone.0135135.t003 year within and across each of the three biogeographical regions, we identified smaller-scale relationships that were obscured at the coast-wide scale. These findings contribute to a growing body of work that suggests sublittoral to bathyal marine ecosystems do not always closely adhere to a strict latitudinal diversity gradient [81,82]. Although it is possible that our maximum latitudinal range of 15°is insufficient to document a strong and unidirectional relationship with diversity in this system, our range is similar to those used in many other studies of latitudinal gradients of marine biodiversity [17]. Furthermore, latitude was an important predictor of diversity even at the smaller regional scales, suggesting that even modest differences in latitude can explain some of the variability in benthic diversity. Overall, our analysis supports a more nuanced view of the relationships between macroecological diversity patterns, latitude, and other environmental drivers of diversity.
Our findings are in line with research conducted in the North Atlantic, which has documented that patterns of diversity can be better attributed to depth and substrate than latitude Table 4. Summary of model-averaged variables included in the confidence model set for fish richness (S), Shannon diversity (H 0 ) and Rao's quadratic entropy (Q) and invertebrate S and H 0 for coast-wide and the three biogeographic regions, North of Cape Mendocino, the Central Region, and South of Point Conception. The Akaike Information Criterion correction (AICc) was used to rank models and any model that ranked <2 ΔAICc was averaged to obtain final estimates. Only model terms which have a 95% confidence interval (CI) that did not include zero are included. Complete model averaged results for the biogeographic regions are presented in S3 Table. doi:10.1371/journal.pone.0135135.t004 [7,28,83]. In addition, depth may act as a proxy for increasing heterogeneity of sediment grain size, a known correlate of infaunal diversity [7]. Etter and Grassle [84], in a study of the western North Atlantic, found that the bathymetric patterns of species diversity were largely attributable to changes in sediment particle diversity, both of which peaked at~1,500 m and declined thereafter. We have no trawl-specific data on substrate grain size or variability, but our observation of an increase in invertebrate diversity with depth is consistent with this previous work. Comparing diversity patterns across broad taxonomic groups, we observed divergent responses to depth by fish and invertebrates. Consistent with findings in the Atlantic Ocean, invertebrate diversity displayed a linear increase with depth (Fig 3; [7 ,28,83]). Deep water bathyal habitats are cold and dark as light only penetrates to a maximum of 1000 m [85]. In deep-water habitats, temperature and productivity are less variable than in nearshore ecosystems [86][87][88]. Climatic stability has been hypothesized to drive high diversity in terrestrial tropical ecosystems, because it promotes speciation and an increased diversification rate [2,6,[89][90][91][92]. Thus, deep-water habitats may be analogous to tropical terrestrial ecosystems in that environmental conditions remain stable over time.
In contrast, fish diversity (H 0 ) across the study range did not show a linear relationship with depth but was highest near the continental shelf break (~200 m, Fig 3). Several non-mutually exclusive processes could explain this mid-depth peak. First, the continental shelf zone is highly dynamic, and influenced by upwelling, currents, tides, and internal waves [93]. Along with seasonal changes in water temperature, upwelling brings cold, nutrient-rich water to the surface and enhances primary productivity. High productivity along the shelf break has been shown to fuel aggregations of zooplankton, micronekton, and fish [94] and support distinct assemblages of species that can attain high biomass [95,96], an explanation that is consistent with the longstanding hypothesis that high primary productivity drives the terrestrial latitudinal diversity gradient [2,97]. Second, as a transition zone between shallow shelf lithosphere and the continental slope, the shelf break is associated with high topographical relief and habitat features such as submarine canyons. Thus, mean trawl depth may mask more important predictors of fish diversity, including depth gradient within a single trawl and habitat complexity, both of which would be expected to affect fish diversity [32,98]. Finally, many benthic fish species display ontogenetic shifts in habitat, whereby larval stages reside in estuaries and nearshore habitats and migrate to deeper waters of the continental shelf as adults [99,100]. Thus, the shelf break may be associated with high fish diversity due to its role as an oceanographic, topographic, and biotic transition zone.
The mid-depth peak in fish richness points to one of the limitations of the use of GLMs. Our model-selection and-averaging procedure restricted the candidate set to linear models, as we want to test our hypotheses regarding latitudinal gradients in diversity, but could have potentially obscured non-linear relationships with predictors. For example, invertebrate S and H 0 display higher values at the latitudinal extremes of the study area, with lower values in the mid-latitude trawls (see Figs 2 and 3). While some of this variability is accounted for by our covariates (i.e. depth, substrate, year), the pattern illustrates that a linear latitudinal gradient may not fully describe the complexity of diversity in the CCLME. While outside the scope of this effort, future work should also compare linear to nonlinear models to better elucidate the complex patterns of marine biodiversity.
Our analysis also revealed a measure of scale-dependency in the influence of predictors on biodiversity. Within each biogeographic region, as we found at the coast-wide scale, diversity was influenced by latitude, but did not follow a simple latitudinal gradient. Rather, our results were more consistent with high diversity areas being driven by regional-scale features and biogeographic transition zones. In the CCLME, biogeographic breaks at Cape Mendocino and Point Conception can restrict gene flow within populations [101,102] and may limit dispersal, altering species compositions and ultimately biodiversity [83,103]. Furthermore, the northern and southern areas of our trawl survey, which represent transition zones with neighboring large marine ecosystems, were particularly rich in fish and invertebrates [80]. The area of high biodiversity off the coast of Washington and Northern Oregon is an ecotone between the Gulf of Alaska LME and the CCLME. Here, the Subarctic Current branches, bringing water from the Northern Gulf of Alaska towards the west coast and the CCLME [85]. This delivery of subarctic water is one mechanism for larval dispersal of subarctic species to habitats in the northern CCLME and may explain our observed peak in biodiversity in this region. The California Bight, at the southern-most range of the survey, comprises another area of relatively high biodiversity. This geographic feature is at the southern range limit of many temperate species and the northern range of subtropical species [51]. Thus, the high local species richness (α-diversity) observed at the transition zones in our study may be explained by high regional richness (γ-diversity) in these areas. Similarly, in a variety of other systems, trends in local richness can be at least partially explained by variation in the regional species pool [3,104,105]. Using only a coast-wide scale approach, these regional nuances could be obscured or misinterpreted as a simple large-scale diversity gradient.
In contrast to traditional metrics of species diversity, fish functional trait diversity was related to several environmental gradients at the regional scale, but not at the coast-wide scale. In each biogeographic region, functional diversity was related to latitude, but this relationship was dependent on different environmental factors in each region. This suggests that functional diversity is sensitive to environmental variation in the CCLME, but the correlates of functional diversity are regional.
In some cases, these regional patterns of functional diversity reflect similarities to patterns in traditional metrics of species diversity (e.g. fish S and Q in the Central region were both positively related to latitude and substrate predictors). However, incorporating functional diversity revealed additional relationships between the marine environment and patterns of community structure that were not found using traditional metrics, suggesting that species functional roles do not necessarily follow patterns of taxonomic distinctions among species. This study contributes to a growing body of literature suggesting that fish functional diversity may not vary consistently along continuous gradients, but instead high diversity may be found in regional "hotspots" [44][45][46].
Inter-annual variation was also an important factor for describing overall biodiversity patterns along the CCLME (Fig 4, S1 and S2 Figs). Year was likely consistently included in our top-ranking models because of strong inter-annual variability of the marine environment along the US west coast. The CCLME is a dynamic ecosystem with divergent climatic conditions occurring on annual timescales. Our study period included a shift from more mild El Niño conditions in the early years of the survey to more turbulent La Niña conditions in the later years [106], which could account for some temporal variation in the data. Larval mobility and settlement in fish and invertebrates may mirror changes in ocean climate conditions [107]; inter-annual variability in these processes may be important in determining biodiversity in any given year. Understanding the drivers of these temporal patterns in biodiversity will be crucial if we are to predict how these areas of high and low biodiversity change annually, over longer time scales, and with a warming climate [33,108].
Our study provides insight into testing macroecological trends in biodiversity using a nontraditional data source. The WCGBTS, a long term, fishery-independent monitoring survey, is meant to collect data to support stock assessments for fisheries management, and basic ecological research is not a primary objective. Despite this, we were able to extract information to answer our fundamental ecological questions on a broad scale. Given the high costs of conducting marine trawl surveys, utilizing previously conducted data to address basic scientific questions is worthwhile and can produce important findings at a low cost, especially when the surveys provide consistent time-series data [109,110]. We propose that these types of datasets represent under-exploited resources for ecologists to conduct basic scientific research. Largescale monitoring projects, frequently conducted by government institutions, can be used to address broad-scale patterns in species distribution and abundance, and investigate the link between these factors to environmental variability. However, care in interpretation of results should be exercised when using these public datasets, which may give an incomplete picture of biodiversity because the sampling methods are designed to meet other objectives.
In marine ecosystems, complex spatial and temporal patterns of biodiversity have important implications for management, particularly spatial management. Currently, just 3% of US waters are protected by no-take MPAs [111]. Most of these exist in small nearshore pockets that may overlook important diversity elsewhere. For example, of marine reserves in California, offshore soft-bottom habitats are still somewhat underrepresented when compared to rocky habitats [112]. These spatial discrepancies in protection may influence fish biodiversity through removals of commercially important species and invertebrate diversity through habitat impacts [68,113]. Commercial fishing has made a transition to deeper water as shallower, nearshore fisheries have become depleted or closed [86]. We have shown here that benthic fish and invertebrates display different relationships between depth and biodiversity. Deepwater commercial fishing could especially impact invertebrate biodiversity, which we found had highest diversity in deeper waters. Given that fishing practices can influence fish and invertebrate biodiversity in various ways (e.g. [114][115][116][117][118]), future research should explicitly test how fishing pressure may interact with the factors that we found to affect biodiversity in the CCLME.
We used a multi-faceted approach to assess macroecological patterns in benthic biodiversity that combined multiple biodiversity metrics, scales and taxa. As such, we were able to document a high degree of spatial and temporal patterning of biodiversity in the California Current. Our study demonstrates that biodiversity in soft bottom deep water habitats in this system is influenced by latitude, but the effect of latitude is mediated by depth, substrate, and interannual variation. Ultimately, the variability in marine benthic biodiversity should be considered in future conservation efforts and for spatial planning.  Table. Biological traits used in functional diversity metrics for fish. (PDF) S3 Table. Generalized linear model-averaged results for fish richness (S), Shannon diversity (H 0 ) and Rao's quadratic entropy (Q) and invertebrate S and H 0 for the three biogeographic regions, North of Cape Mendocino, the Central Region, and South of Point Conception. The Akaike Information Criterion correction (AICc) was used to rank models and any model that ranked <2 ΔAICc was averaged to obtain final estimates presented. Bolded model terms have a 95% confidence interval (CI) that did not include zero. Relative importance (RI) refers to the proportion of output models that contained the term before the model estimates were averaged. (PDF)