Simulated Abrupt Shifts in Aerobic Habitats of Marine Species in the Past, Present, and Future

The physiological tolerances of marine species toward ambient temperature and oxygen can jointly be evaluated in a single metric: the metabolic index. Changes therein characterize a changing aerobic habitat tailored to species‐specific thermal and hypoxia sensitivity traits. If the geographical limits of marine species as indicated by critical thresholds of the metabolic index shift abruptly in response to ocean warming and deoxygenation, aerobic habitat could potentially be lost abruptly. Here, we assess the spatio‐temporal detectability of abrupt shifts in potential habitats for selected marine species within the Shared Socioeconomic Pathway 5–8.5 (SSP5‐8.5) scenario run with the fully coupled Norwegian Earth System Model version 2 (NorESM2‐LM). We use an environmental time series changepoint detection routine and analyze the number and timing of these abrupt changes over the past, present and future. We construct nine ecophysiotypes with low, medium, and high resting vulnerability to hypoxia and sensitivity of hypoxia vulnerability to temperature, respectively, with six different thresholds for minimal oxygen demand. For all ecophysiotypes with positive temperature sensitivity to hypoxia, the volume of non‐viable habitat in the upper ocean expands between 1850 and 2100. Changepoints in the metabolic index are detected in 49.0 ± 9.2% of the volume that eventually becomes non‐viable for all ecophysiotypes over the course of the 21st century. More than 75% of these abrupt shifts occur in response to warming close to the surface, while at depth, the abrupt shifts driven by changes in oxygen partial pressure become more important, with potentially severe consequences for marine species, populations, and ecosystems.


Introduction
Marine habitats can be modified by gradual or shock-induced changes that reduce the strength of stabilizing feedback mechanisms until a critical threshold is crossed, that is, a so called tipping point, which forces a stable regime to reorganize into a new, different stable regime (Biggs et al., 2018).Ocean warming, deoxygenation, or ocean acidification, but also other anthropogenic stressors such as overfishing, plastic contamination, pollution, or eutrophication, can eventually trigger ecosystem-wide regime shifts (Heinze et al., 2021).Such regime shifts within marine ecosystems are increasingly observed in all ocean basins and marginal seas (Levin & Möllmann, 2015;Möllmann & Diekmann, 2012).Different definitions of regime shifts are available (Lees et al., 2006), but in general they are associated with an abrupt, persistent, and substantial reorganisation of ecosystem structure and functioning (e.g., Wernberg et al., 2016), and they are expected to occur more frequently and may be irreversible, if the anthropogenic perturbation remains unabated (Rocha et al., 2015).
The detectability of regime shifts in marine habitats depends on the magnitude and the spatial extent of the shift itself, and is increasingly complicated due to nonlinear dynamics of ecological systems, complex interactions between the physical-chemical environments and biota, and species-dependent physiological tolerances to change (Andersen et al., 2009;de Young et al., 2008;Pelletier et al., 2020).Moreover, large internal variability can obscure gradual changes from abrupt shifts, or, patterns of internal variability may be misinterpreted as trends or shifts (Beaulieu & Killick, 2018).Analyzing changepoints in time series can be used as a first step toward objectively identifying the correct timing of abrupt shifts of statistical properties (Beaulieu et al., 2012), that is, a measure of the abruptness as one of the indicators for a regime shift.For example, the frequency of changepoints in phytoplankton biomass were used to identify abrupt shifts in the global plankton community composition in response to changing ambient conditions (Cael et al., 2021(Cael et al., , 2022)).Further, the presence of abrupt shifts in sea surface temperature and mixed layer depth, but also nutrient concentration, chlorophyll, or primary production, as identified using changepoint analysis, were indicative of regime shifts in the Gulf of Alaska (Beaulieu et al., 2016) or in the North-East Atlantic (Allen et al., 2020).
The habitat boundaries of marine species have been shown to closely align with the species' tolerance toward ambient temperature and oxygen (O 2 ) (Deutsch et al., 2015(Deutsch et al., , 2020)).These physiological tolerances can be combined in a single metric: the metabolic index (Φ).The metabolic index is a species-specific measure of the ratio of oxygen supply in relation to an organisms' resting oxygen demand (Deutsch et al., 2015), and hence can be used to describe the extent of aerobic habitats of marine species but also changes therein.For example, the rapid spreading of temperature-dependent deoxygenation has led to a widespread loss of aerobic habitat and thereby contributed to the marine mass extinction of marine species at the end of the Permian Period (Penn et al., 2018).Likewise, in response to climate change, ongoing warming and deoxygenation will not only affect the future distribution of marine species, but could lead to extirpation of species with small tolerance toward deoxygenation (Howard et al., 2020;Parouffe et al., 2023;Penn & Deutsch, 2022).Depending on the future emission pathway, the projected habitat loss may even be irreversible for centuries (Santana-Falcón et al., 2023).
If the rate of future ocean warming and deoxygenation is slow, marine species could adapt or acclimatize to changing aerobic habitat conditions (Fox et al., 2019;Seebacher et al., 2015).The physiological plasticity in the aerobic energy metabolism, that is, the ability of individuals to produce a different phenotype, allows species to adjust standard metabolic rates to different temperatures (e.g., Oellermann et al., 2022).Beyond that, because metabolic requirements change with body size; that is, the oxygen demand varies across the size spectrum (Deutsch et al., 2015;Penn et al., 2018), ongoing warming can lead to a reduction in the body size of marine species (Deutsch et al., 2022).This thermally driven decrease in oxygen demand would change local aerobic habitat requirements.However, especially for species that are only adapted to low fluctuations in available aerobic habitat today (Mora et al., 2013), the capacity to adapt to warming in such way may be limited and is likely rate-dependent.So, if future warming and deoxygenation occurs rapidly, this could lead to abrupt changes in the metabolic index, and the consequences for marine ecosystems would be even more severe (Turner et al., 2020).
In this study, the detectability of abrupt shifts in the metabolic index is assessed using the environmental time series changepoint detection routine developed by Beaulieu and Killick (2018) for a range of conceptualized ecophysiotypes that represent different aerobic habitats.The metabolic index describes ambient conditions of a potential habitat and not the actual distribution of marine species.Therefore, abrupt changes in the metabolic index can only be indicative of regime shifts in the aerobic habitat, and henceforth, only abrupt shifts will be discussed.Abruptness needs to be defined in relation to a reference period.Here, annual averages of the metabolic index are analyzed and abrupt shifts describe shifts between decadal to multi-decadal mean states in potential marine habitats.Using simulations over the historical period and in a high CO 2 future conceptualized with the Shared Socioeconomic Pathway 5-8.5 (SSP5-8.5),which are both run with the fully coupled Norwegian Earth System Model version 2 (NorESM2-LM), this study aims to answer the following question: Does the number and the timing of abrupt shifts in the metabolic index and their driving mechanisms change?

Model and Experiment Description
All experiments in this study were run using NorESM2-LM output (Seland et al., 2019(Seland et al., , 2020)).NorESM2 is based on the Community Earth System Model (CESM2.1,Danabasoglu et al., 2020), but uses the Bergen Layered Ocean Model (BLOM, Bentsen et al., 2013), which is a modified version of the Miami Isopycnic Coordinate Ocean Model (MICOM, Bleck & Smith, 1990); it further uses the isopycnic coordinate Hamburg Ocean Carbon Cycle model for ocean biogeochemistry (iHAMOCC, Tjiputra et al., 2020), and the atmospheric aerosol module OsloAero6 (Kirkevåg et al., 2018), as well as tuning of the atmospheric component following Toniazzo et al. (2020).NorESM2-LM has an atmosphere-land horizontal resolution of nominal 2°, while the ocean and ice components are run with 1°resolution.BLOM deploys 51 isopycnic model layers and two non-isopycnic surface layers representing the mixed layer.Here, annual output from a single member (r1i1p1f1) is used from two standard Coupled Model Intercomparison Project Phase 6 (CMIP6, Eyring et al., 2016) scenarios: historical andSSP5-8.5 (ScenarioMIP, O'Neill et al., 2016).All model outputs used in this study have been regridded to a 1°× 1°-resolution regular grid.

Metabolic Index
Historically, hypoxia is often defined by a single oxygen concentration threshold applied to a wide range of temperature and depth levels in the global ocean, while for the respiration of marine species the partial pressure of oxygen (pO 2 ) is the critical parameter (Hofmann et al., 2011).Therefore, pO 2 of the ambient water (see Supporting Information S1), and its temperature T are used to calculate the metabolic index following Deutsch et al. (2015Deutsch et al. ( , 2020)): Further, k B is the Boltzmann constant, α d is the resting metabolic rate per unit body mass (B) at a reference temperature (T ref ), and α s the efficacy of O 2 supply per unit body mass.The metabolic index is estimated assuming that the allometric scaling of the supply-to-demand ratio ϵ (Equation 1) approaches zero (Nilsson & Östlund-Nilsson, 2008), such that B ϵ approaches one.Although body size affects the level of metabolism, and hence the response of individual species to warming and deoxygenation (Duncan et al., 2023), here, a constant organism mass is assumed.
The first key physiological trait of an organism is its resting vulnerability to hypoxia at the reference temperature, V h = α d /α s .The second key physiological trait is the temperature sensitivity to that hypoxic vulnerability described by E 0 .The third trait describes the minimum pO 2 demand for essential ecological activities, that is, feeding, defense, growth, or reproduction, which is a factor of Φ above resting levels (Φ crit ).As such, Φ crit is a species-dependent measure for hypoxia, because below a metabolic index of 1 (Φ < 1) and below Φ crit , not even the resting metabolic demand of an organism is met, leading to physiologically unsustainable conditions where either aerobic activity is effectively inhibited or anaerobic metabolism is initiated (Deutsch et al., 2015).

Decomposition of Changes in the Metabolic Index
The mechanisms responsible for spatio-temporal changes in metabolic habitats are further investigated by considering the total derivative of the metabolic index.Here, the relative contributions from pO 2 and T are taken into account.The contribution to dpO 2 is further decomposed into a solubility and a combined disequlibrium, circulation, and remineralization component (Frölicher et al., 2009), that is, the equivalent to oxygen solubility (O * 2 ) and apparent oxygen utilization (AOU) expressed in units of pressure: As pO * 2 is essentially a function of temperature and only to a lesser extent salinity in seawater (Garcia & Gordon, 1992), a thermally driven regime change in the metabolic index is constructed by adding ∂Φ In order to decompose changes in the metabolic index into thermal and nonthermal contributions, time series are divided into segments using the changepoints (see Section 2.5) to avoid fitting a linear trend in presence of shifts (Beaulieu & Killick, 2018).The period prior to the changepoint is used to evaluate driving mechanisms by considering dpO 2 dt and dT dt (see Equation 2), further, the contributions are classified according to the direction of change in Φ, that is, into positive and negative changes in the metabolic index.Only negative changes in Φ are considered that indicate a loss in metabolic habitat, that is, changepoints after an increase in Φ are not accounted for.Deutsch et al. (2020) compiled the key metabolic traits for 61 species including their critical thresholds, that is, Φ crit (see Figure 1).Based on this data set, nine ecophysiotypes are constructed with low, medium and high resting vulnerability to hypoxia (V h ∼1.2 kPa, 4.5 and 10.4 kPa, respectively) and sensitivity of hypoxia vulnerability to temperature (E 0 ∼ 0.03, 0.34 and 0.86 eV, respectively), representing the 5th, 50th, and 95th percentile of the key metabolic traits.These different ecophysiotypes cover a range of potential metabolic habitats over the same temperature and pO 2 limits as illustrated in Figure 2. Here, exemplary critical thresholds for each ecophysiotype are shown: Φ crit = 1, 2, 3, or 5.

Construction of Ecophysiotypes
While type-specific hypoxia vulnerability V h and their temperature sensitivities E 0 were constructed with the help of percentiles, simply assuming a single Φ crit , that is, the third key physiological trait of a species, for each ecophysiotype would lead to large errors (see Figure 1), because individual species have specific Φ crit .Some of the species used for the construction of ecophysiotypes have values of Φ crit as large as seven or eight (Deutsch et al., 2020), however, most values of Φ crit are smaller than five.Therefore, henceforth six values of Φ crit for each ecophysiotype are depicted: Φ crit = 1, 1.5, 2, 2.5, 3, and 5.In other words, each value of Φ crit in each ecophysiotype represents a different, albeit constructed, species (a total of 54 species).To illustrate this, Table S1 and Figure S1 in Supporting Information S1 show an exemplary species for each ecophysiotype based on the data provided by Deutsch et al. (2020).

Changepoint Analysis
Changepoint detection has been widely used to identify abrupt shifts in various statistical properties of time series data (Beaulieu et al., 2012;Beaulieu & Killick, 2018).Here, the EnvCpt changepoint package in R by Killick et al. (2020) has been employed, which fits 12 different models to individual time series: six models without changepoints and six models with (multiple) changepoints in all model parameters.These models are (1) a constant mean, (2-3) a constant mean with first-or second-order autocorrelation, (4) a linear trend, (5-6) a linear trend with first-and second-order autocorrelation, which characterizes the memory inherited in the system.The best fit to the individual time-series is identified with the Akaike information criterion (AIC; Akaike, 1974), which considers the number of fitted parameters to penalize the model likelihood as a guard against potential overfit (Beaulieu & Killick, 2018).If a model with changepoints is considered to be the best fit, a pruned exact linear time (PELT) algorithm in combination with a modified Bayesian information criterion (MBIC) as the penalty function is used to estimate the optimal number of changepoints (Killick et al., 2012).Note that the choice of the information criterion has been shown to influence the timing and structure of the detected changepoints.For example, the Bayesian Information Criterion, BIC, has been shown to select a changepoint model less often than the AIC and to find different changepoints in the time series (Kaderli et al., 2020).
Changepoints in the metabolic index for all ecophysiotypes are calculated for timeseries from 1850 to 2100 for all grid points and depth-levels between 0 and 700 m depth.Choosing the length of the regarded time series, that is, the number of data points, affects the number and timing of detected changepoints as a different internal variability structure is recorded in the time series.Changepoint detection has been shown to be more robust within longer timeseries (Killick et al., 2012), justifying the analysis between 1850 and 2099.The lower depth boundary was chosen to include most of the observed occurrences of the species  1 from Deutsch et al. (2020).The constructed ecophysiotypes LH, MH, HH, LM, Median, HM, LL, ML, and HL (the first letter describes vulnerability to hypoxia V h , the second letter describes their temperature sensitivity E 0 ) are based on the 5th (Low, "L"), 50th (Median, "M"), and 90th (High, "H") percentile of all species and are plotted on top (white diamonds).
Earth's Future 10.1029/2023EF004141 FRÖB ET AL. used in Deutsch et al. (2020).Although the PELT algorithm already ensures that temporal segments between changepoints are not too small (Killick et al., 2012), the minimum segment length is set to 10 years, that is, the period between two changepoints in a single time series with multiple changepoints is forced to exceed 10 years.Thus, changepoints illustrate shifts in decadal to multi-decadal mean states in potential marine habitats.The results show whether a changepoint was recognized, but not which of the selected changepoint models, that is, either trend-based or mean-based with or without autocorrelation, identified this changepoint.Note that changepoints found in the first or last 3 years of each time series are discarded in order to minimize boundary effects.

Spatial Distribution of Metabolic Habitats and Their Change in Time
For the Median ecophysiotype, the metabolic index Φ increases toward the cooler, oxygenated waters of the high latitudes (Figure 3, upper left panel; see also Figures S2 and S3 in Supporting Information S1 for temperature and oxygen distribution).With increasing depth and lower temperature and lower pO 2 , the potential habitat of the Median ecophysiotype decreases, in particular in the low latitudes around the oxygen minimum zones (Figure 3, center and lower left panels).However, the pattern of Φ depends strongly on the metabolic traits (see Figure 2).The more tolerant species are to hypoxia (i.e., the lower V h ), the higher Φ (ecophysiotypes LL, LM, and LH).Therefore, the meridional gradient between the high and low latitudes is larger for ecophysiotypes LL, LM, and LH when compared to the species with median metabolic traits.Vice versa, large areas in the surface ocean in the mid and low latitudes show values below Φ = 1 for species that are more vulnerable to hypoxia (types HM and HH).For ecophysiotypes with negative temperature sensitivities to hypoxia vulnerability E 0 (existent in types LL, ML, and HL, shown for ML in Figure S4 in Supporting Information S1), the meridional gradient in Φ is very small, but reversed compared to the pattern in Φ for positive values in E 0 .For ecophysiotypes with high E 0 (types LH, MH, and HH, shown for MH in Figure S5 in Supporting Information S1), Φ in the high latitudes can be larger than 25, but is low in the low latitudes.Here, in the tropics, the colder temperatures below surface and the decline in pO 2 with depth lead to a strong gradient in Φ for species with low E 0 (types LM, Median, and HM), which is less pronounced if E 0 is high (because Φ is already low in the entire water column).With ongoing ocean warming and deoxygenation until the end of the 21st century, the extent of the potential habitat for the Median ecophysiotype decreases in general more rapidly in the surface ocean, except in the North Atlantic subpolar gyre where the metabolic conditions improve in response to surface cooling (right panels of Figure 3, see also Figures S2 and S3 in Supporting Information S1).The expansion of non-viable habitat (Φ < Φ crit ) in, for example, the warm and low-oxygenated tropical surface waters illustrates a particularly vulnerable area to such warming and deoxygenation (right panels of Figure 3).At the same time, in response to the warming, the viable habitat of tropical species shifts polewards, enforcing there tropicalization of marine species (not shown).At depth, changes in the non-viable habitat are minor, and the nonlinear response to changes in both temperature and oxygen (potentially caused by changes in organic particle fluces and different remineralization rates) leads to both positive and negative changes in Φ.
The change of metabolic habitats over the course of the 21st century depends strongly on the metabolic traits (see Figure 2).The volume of the global potential habitat may decline or increase depending on the rate of change and individual temperature-dependent hypoxia traits.Most species have positive values for E 0 , the temperature sensitivity to hypoxia vulnerability, meaning that the O 2 supply over demand ratio decreases with temperature, that is, ocean warming leads to a loss of aerobic capacity: Deutsch et al., 2020).For ecophysiotypes with negative temperature sensitivities to hypoxia vulnerability E 0 (existent in types LL, ML and HL, shown for ML in Figure S4 in Supporting Information S1), the anomalies in Φ over the course of the 21st century are minimal, but positive, because the warming actually leads to an increase in Φ (see Figure 2).While the overall pattern of Φ anomalies between the end of the 21st century and a preindustrial climate is preserved for ecophysiotypes with higher E 0 (types LH, MH, HH) compared to those with median E 0 (types LM, Median, HM), the overall magnitude of change is larger; that is, losses and gains in potential habitats are intensified.

Barrier for Potential Habitats: Φ crit
Not all ecophysiotypes are projected to experience the same loss in their potential habitat by the end of the 21st century, but certain combinations of metabolic traits show a very different response to ongoing warming and deoxygenation between 1850 and 2100 (see Figure 4).The boundary between viable and non-viable habitat of a species is indicated by Φ = 1, while the threshold between resting and active metabolic rates is described by the critical value of Φ, that is, where the minimum O 2 demand for essential ecological activities is met.Assuming that the boundaries of the geographical ranges of a species' habitat strongly align with Φ crit (Deutsch et al., 2020), the volume of the potential habitat is calculated by adding the volume of all grid cells where Φ > Φ crit (see Supporting Information S1 for further details).For the nine ecophysiotypes and the different considered critical thresholds, the lost viable volume generally follows the geographical limits of the metabolic index that allows an active metabolism (Figure 4 and Figure S6 in Supporting Information S1).These limits are located further north with higher Φ crit , but they also shift northwards over time if the temperature sensitivity to hypoxia is positive (types LM, Median, HM, LH, MH and HH).If all ecophysiotypes across different levels for minimum O 2 demand were put on top of each other, the resulting map would reveal that all ocean basins are affected by loss of potential metabolic habitat by the year 2100.
The relative change in potential habitat volume driven by changes in oxygen and temperature shows that species with negative temperature sensitivities to hypoxia (types LL, ML, HL) or species with weak hypoxia vulnerability experience only small changes in habitat volume in the surface ocean and temporally even a small expansion of their potential habitat at depth, independent of Φ crit (see Figure 5).For the other ecophysiotypes, around 10% of the potential habitat is lost in the upper 200 m with higher losses for higher values of Φ crit .At deeper depths, ecophysiotypes with the highest vulnerability to hypoxia (types HL, HM, HH) show the largest loss in potential habitat, while ecophysiotypes with weak hypoxia vulnerability (types LL, LM and LH) can potentially temporarily experience a small gain in potential habitat.Note that the anomaly in potential habitat does not consistently increase with higher values of Φ crit for all species with positive temperature sensitivity to hypoxia or for species with higher hypoxia vulnerability at the surface or below (see Figure 5). .Spatial pattern of lost habitat volume (Φ falls below Φ crit in an area that was viable in the temporal mean over 1850-1899) integrated over 0-200 m and between 1850 and 2100 using six different values of potential Φ crit for nine ecophysiotypes (L is low, M is median, and H is high; the first letter describes vulnerability to hypoxia V h , the second letter describes its temperature sensitivity E 0 ).By removing the temporal and the depth dimensions (by integrating over the preindustrial viable habitat of the first 200 m that eventually drops below Φ crit during the course of the 250 years time span), particularly vulnerable areas become visible.Note that white areas indicate either that over the water column no potential habitat is lost between 1850 and 2100 or that there never was potential habitat in the first place.

Global Pattern of Detected Changepoints
For the Median ecophysiotype and the 6 values for Φ crit , the absolute volume that experiences abrupt shifts in Φ within the boundaries of the lost habitat are shown in Figure 6.The panels show this volume over time (top panel), per latitude (left panel), depth (right panel), and longitude (bottom panel), while the center panel maps out the regional pattern of the lost habitat including abrupt shifts.The total cumulative volume that experiences shits in Φ for the Median ecophysiotype ranges between 0.24 and 5.8 × 10 6 km 3 (for comparison, the volume of the Mediterranean Sea is 3.7 × 10 6 km 3 ).Over time, abrupt changes occur more frequently after 1950, with peaks in the early 2000s, and after 2025, especially visible for Φ crit = 5.The higher Φ crit is, the more changepoints are detected in higher latitudes, which simply reflects the predefined boundaries of lost habitat relative to the viable area in 1850-1899.Over depth, there is a tendency that more abrupt changes are detected, except for Φ crit = 2.5, where most changepoints are detected in the upper 200 m of the water column.Also for the critical threshold Φ crit = 2.5, most changepoints per longitude are detected in the Somali Basin north of Madagascar, whereas for all other thresholds, the most changepoints per longitude are found in the Pacific Ocean.All of these signals across different axes reveal that within the habitat that is eventually lost, abrupt changes in Φ potentially occur in all regions of the global ocean and over all depth levels.Over all the considered nine ecophysiotypes and six values of Φ crit , 49.0 ± 9.2% of the viable habitat that is eventually lost between 1850 and 2100 between 0 and 200 m experiences abrupt shifts in Φ (also see Section 4.2).Prior to 1950, almost no changepoints are detected, that is, abrupt shifts in Φ occur, except for the ecophysiotypes with high vulnerability to hypoxia and a negative temperature sensitivity to hypoxia (see Figure 7).The volume where abrupt shifts are detected rapidly increases after 1950 to approximately 0.5% per year of the total lost habitat.Most of the abrupt shifts are detected between 2030 and 2040, that is, approximately 1% per year of the total lost habitat.This pattern is generally reflected over the entire water column as well (see Figure S7 in Supporting Information S1).Some ecophysiotypes, especially those with negative temperature sensitivities to hypoxia, may experience a gain of potential habitat (see Figure 2), but changepoints therein are not considered here.The peaks in the lost volume that experiences abrupt shifts in Φ in the 1960s, early 2000s, and a maximum around 2030 are independent of the ecophysiotypes.This relates well to a large-scale reorganization of sea surface temperature patters centered around the late 1960s (Baines & Folland, 2007).Di Cecco and Gouhier (2018) also find increases in temporal and spatial autocorrelation patterns in temperature in the 1950s and around 2030 (in the Figure 7. Fraction of lost volume that experiences changepoints within the lost aerobic habitat (Φ falls below Φ crit in an area that was viable in the temporal mean over 1850-1899) between 1850 and 2100 integrated over latitude and longitude and over 0-200 m depth for nine ecophysiotypes (L is low, M is median, and H is high; the first letter describes vulnerability to hypoxia V h , the second letter describes its temperature sensitivity E 0 ) for the SSP5-8.5 scenario.The absolute volume that experiences abrupt shifts in Φ differs for all ecophysiotypes and for different values of experiences abrupt shifts in Φ differs for all ecophysiotypes and for different values of Φ crit , which is why the fraction to the absolute lost volume is depicted.Note that for a better visibility a 10-year running mean is shown.representative concentration pathway RCP8.5), indicating shifts in statistical properties within global and regional temperature time series.As the application of the changepoint routine consistently detects most changes at similar time periods, this points to one potentially important driver of abrupt metabolic index changes, that is, temperature.Note that changepoints are detected relative to the dynamics of the underlying system, that is, a high variance in Φ in the surface ocean due to high variance in temperature does not immediately equate to a larger number in changepoints, instead, the parameters of the underlying distribution need to change abruptly.Finally, the decline in the fraction of lost habitat volume including changepoints after the 2030s is simply because Φ falls below Φ crit more frequently, that is, the potential habitat is indeed lost and changepoints hereafter are no longer counted.
Changepoints in Φ that occur well within the viable habitat, that is, Φ ≫ Φ crit may be less severe than changepoints that occur closer to the boundaries between viable and non-viable habitats, that is, when Φ Φ crit tends to zero.The distance between Φ and Φ crit at the moment when the changepoint occurs (see Figure 8 and Figure S8 in Supporting Information S1), can be interpreted as a measure of the severity of the abrupt shift in Φ.This distance exceeds 0.75 for most changepoints that occur for ecophysiotypes with a low vulnerability to hypoxia V h (types LL, LM, and LH), but is smaller than 0.25 for most changepoints detected for ecophysiotypes with a high vulnerability to hypoxia V h (types HH, HM, and HL).This is to some extent related to the curvature of the lines of the same Φ in T-pO 2 space (see Figure 2), but the fact that for all ecophysiotypes the fraction of changepoints increases with decreasing distances between Φ and Φ crit until 2100, indicates that abrupt shifts in Φ become more severe.In other words, the probability that a habitat will shift abruptly from viable to non-viable aerobic conditions, increases over the course of the 21st century.This is especially relevant as there is uncertainty associated with the estimated values of Φ crit (Deutsch et al., 2020), some species not having a clear Φ crit in the first place, and within this study, only constructed values of Φ crit were chosen.

Drivers of Abrupt Metabolic Habitat Changes
For all ecophysiotypes, there is a tendency that the more habitat is lost, the higher Φ crit is (see Figure 9).Approximately 50% of that habitat between surface and 200 m that is eventually lost over the course of the 21st century experiences abrupt shifts in Φ (see Figure 9 and Figure S8 in Supporting Information S1 for 0-700 m).This fraction as integrated over depth and time is independent of Φ crit and the ecophysiotype.The ecophysiotypes with positive temperature sensitivity to hypoxia vulnerability (E 0 ) show that more than 75% of the habitat that experiences abrupt shifts in Φ is driven by changes in temperature, that is, the thermal regime, except for type LM.Fraction of number of changepoints relative to the total number of changepoints binned into 4 classes according to the distance of Φ to Φ crit during a changepoint (ΔΦ ≤ 0.25, 0.25 < ΔΦ ≤ 0.5, 0.5 < ΔΦ ≤ 0.75, ΔΦ > 0.75) in each time step over latitude and longitude for nine ecophysiotypes (L is low, M is median, and H is high; the first letter describes vulnerability to hypoxia V h , the second letter describes its temperature sensitivity E 0 ).All six Φ crit are jointly evaluated between 0 and 200 m depth.Note that for a better visibility a 10-year running mean is shown.
Earth's Future This fraction is slightly lower if the upper 0-700 m are considered, that is, the non-thermal contribution to abrupt changes is larger below 200 m depth (see Figure S9 in Supporting Information S1).For types with negative E 0 , temperature has negligible impact on abrupt shifts, because for these types, warming generally leads to an increase in Φ, not a decrease.The more vulnerable the species are to hypoxia (V h ), the larger is the fraction of lost habitat, and the larger is the fraction of lost habitat where changepoints are detected.

Discussion
Changes in marine habitat of diverse taxa can be described using the metabolic index which takes both environmental and physiological parameters into account.Despite the obvious advantages, it can be challenging to get more generalized results with this approach due to the variety of taxa living in the ocean.Previous studies using the metabolic index therefore have assessed habitability in the ocean by calculating the metabolic index for all species provided by Deutsch et al. (2020) and then constructing a mean habitat (Santana-Falcón et al., 2023) or calculate mean metabolic traits first to calculate the metabolic index for a mean ecophysiotype (e.g., Parouffe et al., 2023;Penn & Deutsch, 2022), or describe only individual species (e.g.,Howard et al., 2020;Parouffe et al., 2023).Compared to these very valid approaches, here, the construction of the nine ecophysiotypes allows to generalize the patterns of habitability and changes therein, while still keeping a species-specific complexity with regards to hypoxia vulnerability V h , its temperature sensitivity E 0 , and critical oxygen demand.Most species reported by Deutsch et al. (2020) are accounted for within this framework without excluding species with unreported metabolic traits, making this approach highly valuable for future assessments.
A recent assessment of the vulnerability of the marine ecosystem to climate change finds that almost 90% of the analyzed species are at critical risk under the SSP5-8.5 scenario (Boyce et al., 2022).In their multi-index analysis, the authors use thermal tolerance margins of marine organisms (e.g., Sunday et al., 2012) to quantify exposure to hazardous climate conditions.When comparing the habitat loss between species with the same temperature sensitivity but different hypoxia vulnerabilities, the habitat loss at 0-200 m depth is significantly higher for species with a higher hypoxia sensitivity (see Figure 4).This applies also to the vertical structure of the habitats.While vertically migrating species are less vulnerable to exposure, that is, migration in the vertical may be an adaptive strategy to overcome potential loss in habitat (Boyce et al., 2022), the thermal habitat in the vertical is projected to contract for most marine species in response to future warming (Santana-Falcón & Séférian, 2022).And as shown here, high hypoxia vulnerability even increases the habitat loss that occurs at depth (see Figure 5).The inclusion of hypoxia vulnerability in habitability assessments beyond the use of thermal tolerance traits is likely to exacerbate habitat loss at depth and hence habitat fragmentation.Therefore, the risk for marine life in a scenario with very high end-of-the-century CO 2 emissions (Riahi et al., 2017) could even be higher than quantified by Boyce et al. (2022).
While it has been shown before, that the loss in habitability is evident for most species in all ocean basins on decadal to multi-decadal timescales (Boyce et al., 2022), here, the potential abruptness of this loss is highlighted.The frequency of changepoints that are detected every year increases sharply after 1950, that is, between 0.5 and 1% of the eventually lost habitat experiences abrupt changes every year (see Figure 7), which may be indicative of regime shifts (Scheffer et al., 2001).Due to the ongoing warming and deoxygenation over the course of the 21st century the metabolic index decreases for most species such that these abrupt changes become more severe, because the probability that the habitat shifts abruptly from viable to non-viable aerobic conditions increases (see Figure 8).Species that are currently exposed to large fluctuations in temperature, for example, on seasonal scales often have a broader thermal tolerance (Nadeau et al., 2017), which probably makes them less sensitive to climate change (e.g., Carilli et al., 2012).However, how abrupt changes erode resilience of such species toward temporal changes in the ambient conditions remains an open question (Turner et al., 2020).The same applies to species that inhabit a larger, less fragmented habitat and therefore have access to potentially favorable habitats (Boyce et al., 2022).As spatio-temporal autocorrelation has been shown to increase for temperature over the course of the 21st century (Di Cecco & Gouhier, 2018), the erosion of spatial variation in the composition of ecological communities would likely be reinforced (Cheung et al., 2013;Vergés et al., 2014), and probably even more so if environmental changes were abrupt instead of gradual.
How vulnerable species are to being exposed to abrupt changes in their aerobic habitats depends also on their ability to migrate or their adaptivity, for example, acclimatization of physiological processes to variations in the environmental conditions of habitats.Potential shifts in the distribution of marine species (Poloczanska et al., 2016), that is, migration toward the poles, in the east-west, or toward the equator in response to the complex patterns of shifting isotherms, could overcome local, abrupt, changes in the aerobic habitat conditions.However, not all species are physically able to migrate, and migration may be inhibited by bathymetry and biogeographic barriers (Pinsky et al., 2020).At the same time, abrupt changes may reduce resilience of local native species and favor thereby invasive species (Diez et al., 2012).Modeled abrupt shifts in the aerobic habitat are detected mostly before 2030, therefore the adaptive potential for marine species is likely limited, especially for those species with generation times measured in decades (Peck, 2011).Tropical species have narrower thermal ranges compared to temperate species and have even less capacity for acclimation through phenotypic plasticity (Tewksbury et al., 2008).Therefore, they may be particularly vulnerable to abrupt habitat change, specifically as biological richness is projected to decline most in the tropics (Penn & Deutsch, 2022).Overall, there is a need for more comprehensive studies on the potential effects of gradual and abrupt shifts to hypoxic conditions on marine ecosystems (Borges et al., 2022).It remains hence unclear what the here identified abrupt changes in the potential habitats would mean for the species inhabiting them, but it is clear that the abruptness in the loss of habitat is another stressor that makes adaptation more difficult for the concerned species.
The results of this study rely on simulations with a single ESM: NorESM2-LM.There is considerable model uncertainty in projecting the metabolic index due to different projections of ocean warming and deoxygenation (Bopp et al., 2013;Kwiatkowski et al., 2020), and a different ESM will likely yield a different pattern of the number and timing of detected abrupt shifts therein.Inter-model uncertainty of the projected warming between 1850 and 2100 has been linked to a different climate sensitivity within CMIP6 models (Kwiatkowski et al., 2020); so that NorESM2-LM potentially shows later loss of thermal habitat in the future scenario due to a lower climate sensitivity.Related to this, Santana-Falcón et al. (2023) show that the mean metabolic indices are slightly higher in NorESM2-LM (lower climate sensitivity) than in CNRM-ESM2-1 (higher climate sensitivity).Inter-model uncertainty in oxygen projections is even larger than that of projections of ocean warming (Kwiatkowski et al., 2020), but due to its relatively slow ocean warming, it can be assumed that the deoxygenation in NorESM2-LM is also relatively slow.Overall, it seems likely NorESM2-LM will show a later loss of habitat than other CMIP6 models.Additionally, the overall aerobic habitat loss by 2100 will be smaller in a lower emission scenario.However, most changepoints are detected prior to 2030 and hence in a period of time where differences between the different SSP scenarios are still considerably small (Kwiatkowski et al., 2020).Analyzing nine ecophysiotypes with 6 different critical thresholds for different ESMs would go beyond the scope of this study, the same applies to the analysis of more than one scenario.

Summary and Conclusion
The impact of ocean warming and deoxygenation on marine ecosystems can be jointly evaluated using a single metric: the metabolic index.The metabolic index describes the ratio of oxygen supply over the organisms resting oxygen demand.Changes in the metabolic index describe specific changes in the aerobic habitat based on key thermal and hypoxia sensitivity traits.Here, a framework of nine ecophysiotypes is constructed that represent various levels of resting vulnerability to hypoxia, sensitivity of hypoxia vulnerability to temperature, and critical thresholds of the metabolic index in order to quantify the loss in potential aerobic habitat.The long-term trend in the metabolic index, driven by warming and deoxygenation, is superimposed by rapid and abrupt shifts, the timing of which can be described by changepoints.Using a changepoint detection routine, the spatio-temporal detectability of abrupt shifts in potential habitats is assessed in the NorESM2 SSP5-8.5 scenario run.
For all ecophysiotypes with positive temperature sensitivity to hypoxia, the volume of non-viable habitat expands from 1850 to 2100.The overall potential for loss is evident in all major ocean basins across all depth levels.Ecophysiotypes with a negative temperature sensitivity to hypoxia experience only small changes in habitat volume in the surface ocean and even a small expansion of their potential habitat at depth.Changepoints in the metabolic index are detected in 49.0 ± 9.2% of the volume that eventually becomes nonviable for all ecophysiotypes with different critical thresholds over the course of the 21st century, that is, approximately half of the eventually lost potential habitat experiences abrupt shifts.For most ecophysiotypes, more than 75% of these abrupt changes occur in response to warming, in particular close to the surface, while at depth, the abrupt shifts driven by changes in oxygen partial pressure may not become dominant, but at least more important.
In this study, only potential aerobic habitats were discussed in a framework of different ecophysiotypes.The actual distribution of individual species is likely more localized and more structured.To what extent individual species would experience abrupt changes in their in situ metabolic niche in response to continued anthropogenic forcing, how species respond to environmental change, and how communities are impacted in detail, how resilient marine species truly are to inter-annual to decadal variability in the metabolic index and abrupt shifts therein, remains to be analyzed.Nevertheless, incorporating ecosystem metrics such as the metabolic index into risk assessment studies and analyzing gradual and abrupt changes to assess exposure to hazardous conditions will advance identification of species most vulnerable to climate change.Only then, the loss in habitat can be quantified, potential shifts in species distribution be predicted, and potential protective measures be adapted to conserve and manage marine ecosystems in a changing climate.

Figure 1 .
Figure1.Species-specific metabolic traits based on their hypoxia vulnerability V h and their temperature sensitivity E 0 , as indicated by filled circles.Colors indicate the critical threshold, Φ crit .Data from the Supplementary Table1fromDeutsch et al. (2020).The constructed ecophysiotypes LH, MH, HH, LM, Median, HM, LL, ML, and HL (the first letter describes vulnerability to hypoxia V h , the second letter describes their temperature sensitivity E 0 ) are based on the 5th (Low, "L"), 50th (Median, "M"), and 90th (High, "H") percentile of all species and are plotted on top (white diamonds).

Figure 2 .
Figure2.Range of metabolic indices for nine ecophysiotypes (L is low, M is median, and H is high; the first letter describes vulnerability to hypoxia V h , the second letter describes their temperature sensitivity E 0 ).Gray contour lines show Φ = 1, 2, 3, and 5 to illustrate various species-specific Φ crit .

Figure 3 .
Figure3.Distribution of the metabolic index for the Median ecophysiotype for the surface ocean (top panels), at 200 m (center panels), and 700m depth (lower panels) for the historical and SSP5-8.5 scenario.Left panels show the mean between the years 2070-2099 (left panels), the contour lines show Φ = 2.5 (light gray) and Φ = 1 (dark gray).Right panels show the anomaly with respect to 1850-1899, the contours show Φ = 2.5 between 2070 and 2099 (dark gray) and between 1850 and 1899 (light pink).

Figure 4
Figure 4. Spatial pattern of lost habitat volume (Φ falls below Φ crit in an area that was viable in the temporal mean over 1850-1899) integrated over 0-200 m and between 1850 and 2100 using six different values of potential Φ crit for nine ecophysiotypes (L is low, M is median, and H is high; the first letter describes vulnerability to hypoxia V h , the second letter describes its temperature sensitivity E 0 ).By removing the temporal and the depth dimensions (by integrating over the preindustrial viable habitat of the first 200 m that eventually drops below Φ crit during the course of the 250 years time span), particularly vulnerable areas become visible.Note that white areas indicate either that over the water column no potential habitat is lost between 1850 and 2100 or that there never was potential habitat in the first place.

Figure 5 .
Figure5.Global mean lost habitat volume (Φ > Φ crit ) relative to the total volume for the first 700 m for nine ecophysiotypes (L is low, M is median, and H is high; the first letter describes vulnerability to hypoxia V h , the second letter describes its temperature sensitivity E 0 ) and six different values of Φ crit .The solid lines illustrate the mean habitat volume temporally averaged over 2070-2099 relative to the extent of the habitat volume temporally averaged over 1850 to 1899.

Figure 6 .
Figure 6.Habitat volume that experiences abrupt shifts in 10 13 m 3 for the Median ecophysiotype over time (top panel), per latitude (left panel), per depth (right panel), and per longitude (bottom panel).Map shows the spatial pattern of abruptly lost habitat volume integrated between 0 and 700 m depth and between 1850 and 2100 for six different values of Φ crit .

Figure 8 .
Figure8.Fraction of number of changepoints relative to the total number of changepoints binned into 4 classes according to the distance of Φ to Φ crit during a changepoint (ΔΦ ≤ 0.25, 0.25 < ΔΦ ≤ 0.5, 0.5 < ΔΦ ≤ 0.75, ΔΦ > 0.75) in each time step over latitude and longitude for nine ecophysiotypes (L is low, M is median, and H is high; the first letter describes vulnerability to hypoxia V h , the second letter describes its temperature sensitivity E 0 ).All six Φ crit are jointly evaluated between 0 and 200 m depth.Note that for a better visibility a 10-year running mean is shown.

Figure 9 .
Figure9.Fraction of lost habitat with no changepoints (dark blue bars) as well as thermal and non-thermal contribution of the abruptly lost habitat (purple and light blue bars).Results are shown for 0-200 m according to six different Φ crit and nine ecophysiotypes (L is low, M is median, and H is high; the first letter describes vulnerability to hypoxia V h , the second letter describes its temperature sensitivity E 0 ).Note that the y-axis have different increments.