The rise of a marine generalist predator and the fall of beta diversity

Determining the importance of physical and biological drivers in shaping biodiversity in diverse ecosystems remains a global challenge. Advancements have been made towards this end in large marine ecosystems with several studies suggesting environmental forcing as the primary driver. However, both empirical and theoretical studies point to additional drivers of changes in diversity involving trophic interactions and, in particular, predation. Moreover, a more integrated but less common approach to the assessment of biodiversity changes involves analyses of spatial β diversity, whereas most studies to date assess only changes in species richness (α diversity). Recent research has established that when cod, a dominant generalist predator, was overfished and collapsed in a northwest Atlantic food web, spatial β diversity increased; that is, the spatial structure of the fish assemblage became increasingly heterogeneous. If cod were to recover, would this situation be reversible, given the inherent complexity and non‐linear dynamics that typify such systems? A dramatic increase of cod in an ecologically similar large marine ecosystem may provide an answer. Here we show that spatial β diversity of fish assemblages in the Barents Sea decreased with increasing cod abundance, while decadal scale changes in temperature did not play a significant role. These findings indicate a reversibility of the fish assemblage structure in response to changing levels of an apex predator and highlight the frequently overlooked importance of trophic interactions in determining large‐scale biodiversity patterns. As increased cod abundance was largely driven by changes in fisheries management, our study also shows that management policies and practices, particularly those involving apex predators, can have a strong effect in shaping spatial diversity patterns, and one should not restrict the focus to effects of climate change alone.


| INTRODUC TI ON
Ecologists strive to understand the processes responsible for generating and modifying diversity in natural ecosystems Ricklefs, 1987). In particular, the study of the spatial variation in species composition, or beta (β) diversity, provides a link between processes operating at a local scale and affecting alpha (α) diversity and processes affecting changes in regional or gamma (γ) diversity.
The variability of local species richness trends (Blowes et al., 2019;Cardinale, Gonzalez, Allington, & Loreau, 2018;Vellend et al., 2013) has led to a call for analyses beyond species richness to achieve a more integrated understanding of processes acting at multiple scales Hillebrand et al., 2018;Socolar, Gilroy, Kunin, & Edwards, 2016). Teasing apart the relative importance of underlying processes, however, has been difficult, and our understanding of the relative importance of multiple fundamental mechanisms remains limited (Ohlmann et al., 2018;Xing & He, 2019).
Two types of processes have been emphasized: stochastic processes (associated with dispersal-colonization-extinction, ecological drift and priority effects), and deterministic processes (coupled to environmental factors, niche differences and interspecific interactions; Chase, Biro, Ryberg, & Smith, 2009). Stochastic processes will tend to increase β diversity whereas deterministic processes (if they vary little over space and time) will lead to lower β diversity.
Some processes, such as those linked to invasion or disturbance, are associated with lower β diversity, even if α diversity increases (Finderup Nielsen, Sand-Jensen, Dornelas, & Bruun, 2019;Iacarella et al., 2018;Stotz, Gianoli, & Cahill, 2019). Other processes, such as biotic interactions and climate change, may increase or decrease β diversity (Socolar et al., 2016). Predators could either strengthen stochastic processes and increase β diversity by reducing population size of all prey populations, or deterministically lead to extinction of rare species and therefore more homogenous communities (Chase et al., 2009). Chase et al. (2009), in one of the first empirical studies to assess the impact of predation on β diversity, found that predatory fish could lower prey β diversity in ponds. Using a model based on the equilibrium theory of island biogeography, Ryberg, Smith, and Chase (2012) suggested that predation effects may depend on feeding behaviour: generalist predators that feed indiscriminately may increase β diversity (and decrease α diversity), whereas specialist predators may decrease β and α diversity by selectively removing only particular prey species. In other empirical studies, the effect of predator-prey interactions on stochastic colonization and extinction processes were idiosyncratic and depended on prey and predator identities (Cirtwill & Stouffer, 2016), and on the relative importance of bottom-up versus top-down effects (Antiqueira, Petchey, dos Santos, de Oliveira, & Romero, 2018).
The effects of variation in climate on β diversity are even less understood (Socolar et al., 2016), and may be regionally specific, depending on species-specific dispersal abilities (Magurran, Dornelas, Moyes, Gotelli, & McGill, 2015;Qian & Ricklefs, 2007). If only a few species are able to alter and broaden their spatial distributions in response to climate change, a reduction in β diversity (biotic homogenization) might follow (Socolar et al., 2016). Because multiple drivers will impact diversity simultaneously, patterns of change in spatial β diversity may vary as much as those observed for species richness: a recent synthesis of diversity time series found large variation in temporal trends of both local diversity and turnover (Blowes et al., 2019).
Current research in marine ecosystems has focused on identifying linkages between large-scale diversity patterns and thermal tolerances associated with climate variability. For example, Magurran et al. (2015) argued that an increase in the turnover of groundfish assemblages west of Scotland over three decades was due to a spatially uneven pattern of rising ocean temperatures which, in turn, led to an increase in similarity of the fish assemblages over large geographic scales. At a larger spatial scale, Fossheim et al. (2015) suggested that recent warming in the Barents Sea has resulted in a rapid transformation of the fish assemblages, due primarily to the incursion of several large, migratory predator fish species which have displaced some of the resident Arctic fish species, potentially leading to more homogeneous fish communities at the scale of the Barents Sea as a whole.
While the focus on the effects of climate change is understandable, many large marine ecosystems have also been affected by changes in trophic interactions, particularly declines in many top predators (Baum & Worm, 2009). The idea that predation might be important in structuring large marine ecosystems has been slow to develop (Frank, Fisher, & Leggett, 2015;Verity & Smetacek, 1996), in contrast to the prominence of predation as a structuring force in benthic and intertidal marine communities (Estes & Palmisano, 1974;Paine, 1980). In freshwater ecosystems, comparative analyses of systems with and without predators, and whole-lake manipulations, have provided compelling evidence for the role of predation in structuring communities (Carpenter & Kitchell, 1993;Chase et al., 2009). In the north-western Atlantic, a natural experiment unfolded during the early 1990s as Atlantic cod (Gadus morhua), a large-bodied generalist predator, was overfished and driven to commercial extinction. The absence of cod revealed that it exerted a strong controlling influence on community structure, extending to the base of the food chain. The loss of cod resulted in the emergence of an ecosystem state dominated by a suite of its former prey and meso-predator species (Frank, Petrie, Choi, & Leggett, 2005). Moreover, the reduction in predation pressure associated with the depletion of northwest Atlantic cod also impacted the spatial pattern of co-occurring groundfish species (Shackell, Fisher, Frank, & Lawton, 2012). Furthermore, Ellingsen et al. (2015), using long-term, large-scale monitoring data from the Scotian Shelf, found that fish assemblages became more heterogeneous; that is, the assemblages became more dissimilar from place to place as measured by β diversity (Anderson, Ellingsen, & McArdle, 2006; when cod no longer dominated (see Figure 1). Notably, abiotic factors, including bottom water temperature, were unimportant in explaining the changing patterns of β diversity observed on the Scotian Shelf. However, changes in recorded bottom water temperature values over the time period examined for the Scotian Shelf were rather small, so it remains unclear how diversity in fish communities might be altered by concomitant effects of cod depletion alongside potential changes in climate.
Recently, the resident cod population in the Barents Sea, a large subarctic ecosystem, increased by more than sixfold since the historic low levels of the 1980s (ICES, 2018; see Figure 1). A major reduction in fishing mortality combined with the near abolishment of illegal, unreported and unregulated fishing practices from the late 1990s to the present were directly responsible for this increase (Kjesbu et al., 2014). The Barents Sea has many similarities to the Scotian Shelf in terms of species composition, thermal environment and availability of high-resolution monitoring data, and thus could provide a test of the generality of the structuring role of cod, in a system where changes in temperature were also substantial (Lind, Ingvaldsen, & Furevik, 2018). Many highly perturbed ecosystems shift from one state to another, but a reversal back to a pre-disturbance state rarely occurs following removal of the perturbation (Scheffer, Carpenter, & Young, 2005). If this process were reversible, we would expect that the recent increase in Barents Sea cod would result in a decline of β diversity, that is, a broad-scale spatial homogenization of resident fish assemblages.
Clearly, a robust test of this hypothesis requires consideration of concomitant changes through time in other potential drivers, such as temperature.
The situation in the Barents Sea permits an assessment of the relative importance of predation by cod versus temperature in structuring fish biodiversity, due to the fact that it is a large, heterogeneous system . Specifically, different areas within it have exhibited divergent temporal trends in temperature and abundance of large, mature cod.
We delineated large subregions (of approximately equal spatial extent) having contrasting patterns in cod abundance and temperature through time, and the patterns of β diversity within each of these subregions were assessed. This design-based approach (Bråthen et al., 2007;Butsic, Lewis, Radeloff, Baumann, & Kuemmerle, 2017;Shadish, Cook, & Campbell, 2002) to the analysis of observational data allowed us to disentangle the relative importance of an apex predator (cod) and a physical climate variable (bottom temperature) on the β diversity of the resident fish community. Ten years of monitoring data derived from a systematic, large-scale scientific sampling program conducted annually from 2004 to 2013 were used for the analysis.

F I G U R E 1
Is there a reversible dynamic relationship between apex predator abundance and fish community structure in large marine ecosystems? A recent study on the Scotian Shelf (SS), Canada demonstrated that fish assemblages became increasingly heterogeneous (β diversity increased) with decreasing cod populations . In contrast, the cod population in the Barents Sea (BS) has increased recently, and we predicted that the return of cod would cause increasing homogenization of fish assemblages (i.e. declining β diversity). Lower values of β diversity correspond to fish assemblages that are more similar from place to place. Inset figure: Y-axis denotes cod spawning stock biomass scaled to a maximum value of 1 for each time series: the SS (red) and the BS (blue). The dashed blue box identifies the time period investigated here with shallowest water in the south-east (mean depth in subregion 11:115 m) and deepest in the west (subregion 6:348 m). The data set from this survey is the most spatially extensive available from the Barents Sea (Johannesen, Høines, Dolgov, & Fossheim, 2012).
We have not included data from other surveys of the Barents Sea in our analyses due to data incompatibility, in particular with regard to species' identifications and the restricted spatial extent of other surveys.
The lengths of all cod from each trawl set were measured, and the age of one cod individual from each 5 cm length group was determined by sectioning otoliths. The number of cod per age group per trawl station was calculated using length-age keys. We did not take into account the uncertainty associated with the length-age keys. As the mean age at maturity is around 7 years (ICES, 2018), we counted abundance as the number of cod aged 7 years or older to a maximum of 17 years (here called mature) in our analysis, particularly because our focus was on the role of large, apex predators in structuring fish biodiversity. The mean length of 7 years old cod was 74 cm (Johannesen, Johnsen, Johansen, & Korsbrekke, 2019), and ontogenetic shifts in Barents Sea cod diet (Holt, Bogstad, Durant, Dolgov, & Ottersen, 2019) result in fish becoming a more important prey item for cod with increasing size of cod.
The data primarily consisted of benthic fish species, although other species, including pelagic fishes, were also routinely captured by the bottom trawl (and included in this study). All fish caught were identified to species level when possible, and for each species, the total number of individuals per trawl was recorded and then standardized per 15 min tow (as described above). Identification keys and taxonomic literature for fish species in the Barents Sea have been incomplete, although older Russian literature exists. There is ongoing work to improve and harmonize species identifications for data obtained from Russian and Norwegian vessels (Mecklenburg et al., 2018;Wienerroither et al., 2011). Difficulties remain regarding identification and unresolved taxonomy for several Arctic fishes, in particular for the families Zoarcidae and Liparidae (Christiansen & Reist, 2013). Juveniles (<10 cm) of the genus Sebastes are also difficult to identify. Thus, fish data were organized into two different data sets as follows: (Option 1) All taxa from the eelpout family were classed as 'unspecified' (Zoarcidae spp.), taxa from the snailfish family were classed as Liparidae spp. and juveniles of the genus Sebastes spp. were pooled with adults into the family Sebastidae spp., resulting in a regional count of 71 fish taxa (excluding cod) from 2004 to 2013; (Option 2) specimens from the Zoarcidae, Liparidae and Sebastidae families were kept at species level, and recordings at the family or genus level, including juvenile Sebastes, were removed, resulting in 95 fish taxa (excluding cod). Bottom temperature data were obtained from CTD (conductivity, temperature and depth) profiles taken adjacent to each trawl set during the ecosystem surveys.

| Design-based definition of subregions
We used a design-based approach (Bråthen et al., 2007;Shadish et al., 2002) to assess the relative roles of increases in mature cod abundance versus changes in bottom temperature on the β diversity of fishes. More specifically, we sought to partition the Barents Sea into spatially contiguous subregions having approximately similar areas (so as not to confound α diversity immediately with β diversity), of regular shapes (not snake-like or elongated, to avoid edge effects) and with consistent and contrasting patterns of change through time in (a) cod and (b) temperature. It was also important to identify regions of different types that were not geographically segregated (i.e. to avoid spatial pseudo-replication-it would not do to have all subregions of a similar type in the north and all regions of a different type in the south).
We defined subregions of the Barents Sea with different annual trends (over 10 years) in average bottom temperature and average abundance of large, mature cod (7+ years, i.e. mean length ≥74 cm; Johannesen et al., 2019; log(x + 1) transformed). As the sampling grid and survey routes were not identical from year to F I G U R E 2 Subregional trends in ocean bottom temperature and mature Atlantic cod abundance in the Barents Sea. Eleven contrasting subregions were identified (grey). Upper half of circles shows temporal trends in bottom temperature from negative (blue) to positive (red). Lower half of circles shows temporal trends in cod abundance from white (zero trend) through to dark green (positive trend). See Section 2 for definition of subregions, Figure S7 for region-specific trends in bottom temperature and abundance of mature cod, and Table S1 for summary statistics year, we did this in two steps. First, we interpolated values for cod and temperature at each point on a joint grid of 50 × 50 km for each year using an additive model with a bivariate tensor product smooth of the spatial coordinates (given in the Lambert Equal Area Projection), and a smooth term for depth (Wood, 2006). In addition, for temperature, we included a smooth effect of sampling date as temperature was still changing in the later part of summer. Second, we regressed interpolated values of cod and temperature (standardized to 200 m depth and median survey date 3 September) versus year to describe the temporal changes in these two variables separately at each point in the spatial grid ( Figure S1a). We then defined subregions of approximately similar sizes guided by maps of the regression coefficients and quantiles of the distribution of these coefficients across areas and years.
For temperature, we identified areas with either stable/decreasing (<0.02°C/year) or increasing temperature (>0.02°C/year), that is, temperature changes either below or above the median, respectively. For cod, we identified areas with stable (<0.19/year), medium increase (0.19-0.38/year), or large increase (>0.38/year) in cod abundance, that is, changes in (log) cod abundance within the first, second or third tertile of the ordered distribution of coefficients respectively. This resulted in 11 subregions ( Figure   S1b). We investigated spatially constrained classification methods, but they resulted in subregions of widely different sizes that were therefore unsuitable for comparative analyses of β diversity. To assess if temperature recorded during the survey in August-September reflected temperatures year-round, we compared the monthly temperature values obtained from an independent data set, along a transect located in the SW part of the Barents Sea, between Fugløya and Bear Island, which cover inflow of Atlantic and coastal water masses to the Barents Sea (Trofimov & Ingvaldsen, 2018). Correlations were always high, varying between 0.67 (for January-August) and 0.89 (for April-August; see Figure S2). Thus, we considered our measures of temperature to be a reasonable proxy for the annual average temperature regimes occurring within each subregion, but we acknowledge that correlations between seasons may vary over time and space (Boitsov, Karsakov, & Trofimov, 2012).
We also considered a classification of the trawls based mainly on topography, a classification based on temperature (Stige et al., 2014), and the Atlantis classification based on hydrography, depth and biology (Certain & Planque, 2015; Figure S3; Tables S1 and S2). All classifications performed substantially better in predicting temperature and mature cod abundance than a model without any spatial classification (Table S1), highlighting how spatially variable temporal trends have been within the Barents Sea. Our design-based classification based on cod and temperature had lower widely applicable information criterion (WAIC; Vehtari, Gelman, & Gabry, 2017;Watanabe, 2013), and explained a larger amount of variation in cod and temperature than divisions based on either topography or temperature (Table S1). The Atlantis division performed better in explaining temporal trends in cod and temperature (Table S1) but had 30 subregions that varied greatly in size (Table S2; Figure S3). Also, several of the Atlantis subregions at the fringe of our study area had fewer than three observations in some years ( Figure S3). In contrast, our cod-and-temperature design-based division consisted of 11 similar-sized subregions (Table S2). We therefore did ensuing analyses of β diversity within each subregion over time using these 11 design-based subregions only.

| Biodiversity measures
All measures of biodiversity were calculated excluding cod. This ensured that our models did not confound predictor variables which includes the effect of unseen species (Chao, Chazdon, Colwell, & Shen, 2005). We measured temporal changes in spatial β diversity within each subregion. The terms 'spatial β diversity', characterizing dissimilarity among assemblages over space for one time period, and 'temporal changes in spatial β diversity', that is, here: temporal change in β diversity within each subregion, should not be confused with the term 'temporal β diversity' which is a measure of change in community composition through time, often quantified as the dissimilarity between each time step and the time-series baseline (see e.g. McGill et al., 2015). We report 'spatial β diversity' as 'β diversity' throughout the article.
All biodiversity measures were calculated on each fish data set (Options 1 and 2) after excluding cod. Jaccard and Chao-Jaccard values were highly correlated with one another (Figure S4), suggesting that detection did not affect estimates of β diversity.
Because α diversity, measured as the annual mean number of species in all trawl sets within each subregion, varied temporally (see Figure S5c), we used Raup-Crick rather than either of the two Jaccard measures. In general, however, all three measures gave similar results.

| Statistical analysis
Prior to analyses and to define the design, bottom temperature data were standardized to 200 m depth and a common median to the species-area relationship (Crist & Veech, 2006), spatial extent was quantified as the area of the spatial convex hull of the group of trawl sets within each subregion for each year. Mean spatial extent within each subregion varied from 54,689 km 2 (subregion 4) to 115,975 km 2 (subregion 3) with the largest range between years in subregion 1 (58,244 km 2 ) and the smallest range between years in subregion 6 (10,584 km 2 ). Pairwise plots of predictor variables ( Figure S6) revealed no systematic collinearity between predictor variables among subregions, suggesting that our design was robust and valid for the purpose of disentangling potential effects of cod and temperature, as well as relevant covariates.
We used a hierarchical Bayesian framework (Gelman & Hill, 2007) to model how β diversity changed over time, in response to changes in mature cod abundance and bottom temperature within each subregion, using the 'rstanarm' library with default priors in R (Goodrich, Gabry, Ali, & Brilleman, 2019;R Core Team, 2018). Note that, although we naturally expect geospatial autocorrelation among subregions and potentially among trawls within subregions, quantitative inferences regarding estimated parameters (and associated errors) are indeed accurate and intact, provided they are restricted to the sampling extent of the Barents Sea only, and also to the specific 10 year time period of this investigation, which indeed is our objective here (de Gruijter & ter Braak, 1990;Pawley & McArdle, 2018). We explored models with a random intercept for each subregion, and models where both the intercept and the slope were allowed to vary for each subregion. We explored the performance of competing models using both WAIC and Leave-One-Out cross-validation Information Criterion (LOOIC) as the latter can be more robust in the case of weak priors or influential observations (Vehtari et al., 2017). We calculated the marginal (fixed effects only) and conditional (fixed plus random effects) r 2 for all models (Nakagawa & Schielzeth, 2013). We standardized the observations by subtracting the mean and dividing by the standard deviation. We present results using (a) unstandardized predictors; (b) standardized predictors; and (c) standardized response and predictor variables. This facilitates comparison of the magnitude of the effects in terms of the regression coefficients for predictor variables on their original measurement scales (unstandardized) and also the relative magnitude of effects of different predictor variables (standardized) on the different diversity measures (Gelman & Hill, 2007). We included both the spatial extent and the bottom depth as additional predictor variables (covariates) in our analyses. As expected, there was a negative effect of water depth on β diversity (more homogeneous assemblages occurred with increasing depth), and a positive effect of spatial extent (β diversity increased with increasing area; Table S3). We did not include measurement errors for cod and temperature values, as it would have required a complex statistical model of their spatial structure and local measurement variability, which is not known (Reeves, Cox, Darby, & Whitley, 1998). Hence, our models and inferences were conditional on the observed values for cod and temperature. When the 95% credible intervals of the posterior distribution of a regression coefficient broadly overlapped with 0, we considered there was no evidence for an effect of that regressor (Mysterud, Yoccoz, Langvatn, Pettorelli, & Stenseth, 2008). We did not use any hard threshold for differences in WAIC or LOOIC values, but combined these values with consideration of the parameter estimates from each model in an overall heuristic assessment (McShane, Gal, Gelman, Robert, & Tackett, 2019).
In addition to β diversity, we also examined effects of cod and temperature on changes in α diversity through time. However, there was no evidence for an overall temporal change in α diversity at the scale of the subregions (Tables S3 and S4; Figure S5c).
Neither mature cod abundance nor bottom temperature could explain variation in α diversity (Table S5).
Analyses were carried out using the two alternative fish data sets with different levels of taxonomic resolution (Options 1 and 2) to elucidate how uncertainty in species' identifications might impact our results and conclusions. Results based on Options 1 and 2 are presented in Tables S3-S5. As expected, the reduced data set (Option 1, where some taxa were lumped into groups) reduced variation among subregions, although the overall patterns of change in β diversity over time were similar to the results obtained by using the data set with more taxa (Option 2). In what follows, we present results based on Option 2.

| Subregional changes in ocean temperature and Atlantic cod
Bottom temperature increased most dramatically in the east (0.8-1.0°C per decade), particularly in subregions 5, 8 and 11 (Figure 2), whereas declining bottom temperatures typified the north (subregion 1) and south (subregions 9 and 10) at a rate of −0.6°C over 10 years (Figure 2; Figure S7). Overall, roughly half of the subregions either warmed or cooled. Due to the large geographic extent of the Barents Sea, average bottom temperatures differed by about 5°C between the south-west and the north-east ( Figure S7).
Temporal patterns of mature cod abundance ranged from almost no change in the south-west to large increases in the north, northeast (subregions 1, 2 and 5) and south-east (subregions 10 and 11; Figure 2; Figure S7). In subregion 1, where one of the largest decreases in bottom temperature occurred, cod abundance also exhibited a large increase. Conversely, some subregions showed increases in both cod and temperature (5, 8 and 11). Clearly, a mixture of opposing and matching temporal trends between cod abundance and temperature was apparent.

| Testing for homogenization of fish assemblages within subregions
There was an overall decrease through time in β diversity (i.e. increasing spatial homogenization) across 11 subregions (Figure 3; Table S3; β diversity was based on the modified Gower dissimilarity metric, see Section 2). The best overall model (according to WAIC; see Section 2) allowed both the slope and intercept to vary among subregions; thus, the rate of change in β diversity through time varied across the subregions (Figure 3; Table S4). The largest decrease occurred in the central part of the Barents Sea (subregions 4 and 7; Figure 3; Figure S5a), and there was no latitudinal gradient in β diversity across the subregions of the Barents Sea ( Figure S5). This model explained 36% of the marginal variation (fixed effects only) and 64% of the conditional variation (fixed and random effects; see Section 2; Table S4) in β diversity. In general, analyses of β diversity based on another dissimilarity metric (Raup-Crick; see Section 2) exhibited similar patterns (Tables S3 and S4; Figure S5b).

| Drivers of changes in β diversity
When comparing the effects of the two main drivers, bottom temperature and abundance of large, mature cod, the declining trend in spatial β diversity was best explained by cod abundance ( Figure 4; Table S5), supporting our original prediction (Figure 1).
There was variation in the magnitude of this effect at the subregional scale, with the largest impact of cod abundance on β diversity occurring in the central part of the Barents Sea (Figure 4).
The model allowing for variation in the effect of cod abundance among subregions performed slightly better than a model without such variation (based on WAIC: −315.6 vs. −315.1; Table S4). In contrast, no effect of bottom temperature on changes in β diversity was evident (Figure 4; Table S5).

| D ISCUSS I ON
Here, we document a decadal-scale pattern of spatial homogenization of fish assemblages in the Barents Sea. This decline in β diversity can be explained by an increase in abundance of cod, a generalist apex predator. Our study supported the prediction generated from a smaller but ecologically similar geographic area (i.e. the Scotian Shelf, one-tenth the size of the Barents Sea; Ellingsen et al., 2015), underscoring the significance of these earlier findings. The design-based approach (see Section 2) provided an effective means of disentangling the relative importance of potential drivers of biodiversity in this large and complex marine ecosystem (Bråthen et al., 2007;Butsic et al., 2017;Shadish et al., 2002). We infer from this collective work across the North Atlantic that a marine apex predator, whether increasing or declining in population size, performs a fundamental role in shaping the β diversity of fish assemblages in large marine systems.
Taken together, the results obtained from these two systems, which show that divergent trends in cod abundance yield predictable and contrasting biodiversity effects, demonstrate the potential for reversible community-level responses to large-scale fisheries management.
The absence of any evidence that bottom temperature affects β diversity may seem surprising, given recent studies highlighting the impacts of climate change on diversity (e.g. Blowes et al., 2019;Magurran et al., 2015;Thuiller, Lavorel, Araújo, Sykes, & Prentice, 2005). There are several possible reasons for the discrepancies between our results and those from previous studies. First, most studies have addressed changes in α diversity (and in particular species richness), and very few have addressed changes in β diversity, particularly in large marine ecosystems (Socolar et al., 2016). Second, our study observed relatively small changes in temperature (at most +1°C) over just 10 years. Temperature ranges for fish species in the Barents Sea are often wide (A.V. Dolgov, unpublished data;see e.g. Wienerroither et al., 2011), so we might not expect impacts on fish communities based on changes in temperatures of ca. 1°C. Such changes are reasonably large, however, compared to overall changes observed since 1950 for the Barents Sea (ca. +2°C; Trofimov & Ingvaldsen, 2018) and also compared to the significant impacts on fish communities reported from the North Sea in response to only a 0.5°C change in summer near-bottom temperatures (e.g. Magurran et al., 2015;Rutterford et al., 2015). Third, if temperatures become more homogeneous across a large area over time, then large-scale β diversity may well be reduced (as in Magurran et al., 2015). The heterogeneous changes in subregional average temperature values observed in our study through time may not correspond to an increased homogeneity of temperature at the broader scale of the Barents Sea as a whole (see also Årthun, Eldevik, & Smedsrud, 2019, for projected changes).
Interestingly, there was no latitudinal gradient in β diversity across the subregions of the Barents Sea, so replacing northern, arctic communities by southern, boreal communities may not affect β diversity.
Further studies should focus on species-specific dispersal rates and how they are related to local climate velocities.
How can an increasing apex predator population have a homogenizing effect on communities? Ryberg et al. (2012) suggested that generalist predators will increase β diversity because they prey randomly irrespective of the prey species, whereas specialist predators may decrease β diversity by removing species underlying species turnover. Our results seem to contradict their findings, and one reason may be that what matters is the spatial heterogeneity of the cod predation pattern compared to the spatial heterogeneity of the fish communities. This has been emphasized when studying the impacts of herbivory on plant communities (Adler, Raff, & Lauenroth, 2001), another form of predation. If predation by cod is spatially very homogeneous, as would be expected from a wide-ranging, opportunistic predator (Holt et al., 2019), then we would expect spatial heterogeneity to decrease (Adler et al., 2001). Both experimental and observational studies have also shown that when predation intensity is relatively high the interaction strengths among competing species are lessened, allowing for increased coexistence and therefore declining spatial segregation (Terborgh, 2015). In the absence of predation, however, competitive relationships are stronger, leading to reduced levels of coexistence and heightened levels of spatial segregation. The degree of spatial segregation is a measure of spatial species turnover and therefore of β diversity (Ulrich, Zalewski, & Uvarov, 2012). Furthermore, cod predation on small arctic fishes (liparids, cottids, etc.; Holt et al., 2019) may lead to their local extinctions and therefore be an additional factor leading to community homogenization (Chase et al., 2009).
Whether the observed changes in β diversity was the result of direct predation from cod, or indirectly through relaxed competition or predation from species heavily preyed upon by cod remain unresolved. Identifying the myriad of possible links among species in the cod's food web was beyond the scope of the current study, and we acknowledge that available time series may not hold sufficient information to disentangle these effects as interactions at the species-to-species level, although our β diversity approach, which examines changes in the spatial structure of the entire fish community concomitantly did so. What we have demonstrated, is that increased cod abundance in the Barents Sea resulted in changes in the entire fish community as shown by the decrease in β diversity in response to increased cod abundance.
The measurement of changes in spatial β diversity through time has added an important new dimension to our knowledge of marine biodiversity. Declines in marine β diversity (i.e. increasing spatial homogenization) found previously have been attributed solely to climate change, and have been considered to be a more important problem than the loss of local species richness . Our results highlight a contrasting view. Fisheries management policies and practices, particularly those involving apex predators, can have a strong effect in shaping overall biodiversity patterns, possibly more than climate change alone. The potential interactive effects of simultaneous changes in human interventions, biotic and abiotic drivers all need further investigation to achieve greater clarity for future decision-making that can enhance and sustain biodiversity in marine ecosystems. Such interactions among multiple drivers are, however, poorly understood (Socolar et al., 2016;Turner et al., 2020). Estimating interactions also requires considerably larger sample size than for estimating main effects (Gelman, 2018), particularly when drivers are correlated (Baum & Worm, 2009). Even if our design-based approach allowed us to disentangle the main effects of cod abundance and temperature changes on spatial β diversity, reliable estimates of interactions should be based on larger sample sizes.
The recent management measures directed towards cod in the Barents Sea, in conjunction with favourable environmental conditions, proved to be highly successful and resulted in a remarkable increase in their biomass. Equally dramatic was the associated change in the spatial structure of co-occurring fish species which may ultimately prove to be beneficial to the persistence of the current ecosystem state. Our study highlights the large impacts on marine ecosystems resulting from fisheries management practices, especially when involving large generalist predators such as cod. Clearly, future management decisions involving fisheries in the Barents Sea will require the inclusion of a dynamical trophic perspective, including climate variability, to better understand how fisheries affect trophic interactions. The utility of such an approach may be judged, in part, by the spatial configuration of its component parts which could potentially become part of a 'tool box' of indicators associated with an ecosystem-based approach to fisheries management.

ACK N OWLED G EM ENTS
We are grateful to many technicians, scientists and research vessel crew of IMR and PINRO who participated in surveys and otherwise contributed to the study. We acknowledge support from the Research Council of Norway (project no. 234359), the Norwegian Institute for Nature Research and the Fram Centre.

CO N FLI C T O F I NTE R E S T
The authors declare that they have no conflict of interest.

AUTH O R CO NTR I B UTI O N S
K.E.E., N.G.Y., T.T. and M.J.A. designed the study. A.V.D., E.J. and K.T.F. contributed data. E.J. and T.T. processed data. N.G.Y. and T.T.

DATA AVA I L A B I L I T Y S TAT E M E N T
Data at the subregional scale that support the findings of this study and R code used to produce all findings will be made available in the Dryad Digital Repository. Raw data from the Norwegian part of the Barents Sea, collected by the IMR, Norway, are available upon request. Raw data from the Russian part, sampled by the Polar Branch of Russian Federal Research Institute of Fisheries and Oceanography (PINRO), are not publicly available due to Russian Law, preventing publication of raw data collected on Russian territory. Access to Russian raw data is granted through contracted collaboration in joint projects with PINRO.