stresses on marine species, with potentially detrimental effects on ecosystems. Here, we evaluate the prospect of using climate velocities (CV) of the metabolic Frontiers in Marine Science OPEN ACCESS EDITED BY Stelios Katsanevakis, University of the Aegean, Greece REVIEWED BY Murray Duncan, Stanford University, United States Emily Slesinger, Alaska Fisheries Science Center, United States *CORRESPONDENCE Alexandra Parouffe alexandra.parouffe@gmail.com SPECIALTY SECTION This article was submitted to Marine Ecosystem Ecology, a section of the journal Frontiers in Marine Science RECEIVED 28 September 2022 ACCEPTED 12 December 2022 PUBLISHED 05 January 2023 CITATION Parouffe A, Garc¸on V, Dewitte B, Paulmier A, Montes I, Parada C, Mecho A and Veliz D (2023) Evaluating future climate change exposure of marine habitat in the South East Pacific based on metabolic constraints. Front. Mar. Sci. 9:1055875. doi: 10.3389/fmars.2022.1055875 COPYRIGHT © 2023 Parouffe, Garc¸on, Dewitte, Paulmier, Montes, Parada, Mecho and Veliz. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author (s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms. TYPE Original Research PUBLISHED 05 January 2023 DOI 10.3389/fmars.2022.1055875index (F) for assessing changes in habitat in the South East Pacific. Methods: Our approach is based on a species with mean ecophysiotype (i.e. model species) and the use of a global Earth System Model simulation (CESM- LE) under RCP 8.5 scenario. The SEP is chosen as a case study as it hosts an Oxygen Minimum Zone and seamounts systems sustaining local communities through artisanal fisheries. Results and Discussion: Our results indicate that CVF pattern is mainly constrained by the oxygen distribution and that its sign is affected by contrasting oxygen trends (including a re-oxygenation in the upper OMZ) and warming. We further show that CVF is weakly dependent on physiological traits composing F, which conveys to this metrics some value for inferring the projected mean displacement and potential changes in viability of metabolic habitat in a region where physiological data are scarce. Based on sensitivity experiments to physiological traits and natural variability, we propose a general method for inferring broad areas of climate change exposure regardless of species-specific F. We show in particular that for the model used here, the upper OMZ region can be considered a “safe” area for the species with ecophysiotype close to that of 71 species used to derive the model species. Limitations of the approach and perspectives of this work are also discussed.Evaluating future climate change exposure of marine habitat in the South East Pacific based on metabolic constraints Alexandra Parouffe1*, Ve´ronique Garc¸on1, Boris Dewitte2,3,4,5, Aure´lien Paulmier1, Ivonne Montes6, Carolina Parada2,7, Ariadna Mecho8 and David Veliz2,9 1Laboratoire d’Etudes en Géophysique et Océanographie Spatiales, Université de Toulouse, LEGOS (CNES/CNRS/IRD/UPS), Toulouse, France, 2Nu´cleo Milenio de Ecologı´a y Manejo Sustentable (ESMOI), Facultad de Ciencias del Mar, Departamento de Biologı´a Marina, Universidad Cato´lica del Norte, Coquimbo, Chile, 3Departamento de Biologı´a, Facultad de Ciencias del Mar, Universidad Cato´lica del Norte, Coquimbo, Chile, 4Centro de Estudios Avanzados en Zonas Aridas (CEAZA), Coquimbo, Chile, 5CECI, Universite´ de Toulouse, CERFACS/CNRS, Toulouse, France, 6Instituto Geofísico del Perú (IGP), Lima, Peru, 7Departamento de Geofı´sica, Universidad de Concepcio´n, Concepcio´n, Chile, 8Laboratoire des Sciences du Climat et de l’Environnement (LSCE/IPSL), Commissariat à l'Energie Atomique (CEA)/Saclay, Gif-sur-Yvette, France, 9Departamento de Ciencias Ecolo´gicas, Facultad de Ciencias, Universidad de Chile, Santiago, Chile Introduction:On-going climate change is now recognized to yield physiologicalKEYWORDS metabolic index, climate velocities, South East Pacific, oxygen minimum zone, deoxygenation frontiersin.org01 1 Introduction as the ratio of temperature dependent oxygen supply to animal oxygen demand rates. This ratio termed the metabolic index F Parouffe et al. 10.3389/fmars.2022.1055875With climate change, marine species face the combined effect of multiple threats: warming, acidification and oxygen decline (Levin, 2018; Sampaio et al. , 2021). Marine deoxygenation may be the most threatening of all stressors to marine life and has important detrimental consequences on marine ecosystems (Sampaio et al., 2021), particularly in the highly productive Eastern Boundary Upwelling Systems that host extended oxygen minimum zones (OMZ) (Garçon et al., 2019). The South East Pacific (SEP) hosts one of the largest OMZ resulting from the weak ventilation of thermocline waters combined with the relatively high rate of microbial decomposition of organic matter (Karstensen et al., 2008; Paulmier and Ruiz-Pino, 2009; Oschlies, 2018; Pitcher et al., 2021). Despite the current deoxygenation trend at global scale (e.g., Keeling et al., 2010; Bopp et al., 2013; Schmidtko et al., 2017; Oschlies, 2018), whether or not the SEP OMZ will expand or shrink remains uncertain. On the one hand, there is still currently a low consensus amongst the latest generation of the Global Earth System Models on the future trajectory of deoxygenation rate in the SEP (Bopp et al., 2013; Cabré et al., 2015) while available observations indicate regional discrepancies in the sign of the trend (Stramma et al., 2008; Graco et al., 2017; Ito et al., 2017). On the other hand, climate projections agree on significant warming along the western coast of South America with a pattern resembling the culminating phase of a strong El Niño event (Cai et al., 2021; Dewitte et al., 2021). Increased surface temperature yields an increased vertical stratification, a phenomenon already observed (Liao et al., 2021). Increased stratification not only alters ocean dynamics in a number of ways, and thus indirectly marine ecosystems, but it has also the potential to modify directly marine habitat through modulating thermal niches (Santana-Falcón and Séférian, 2022). As a matter of fact, shifts in species distribution have already been observed (Perry et al., 2005; Pinsky et al., 2013; Poloczanska et al., 2016) or projected (Cheung et al., 2009; Clarke et al., 2021a) as a result of habitat compression due to increasing temperature. Such changes can be captured by so-called climate velocities, a metrics that encapsulate both changes in space and time of a particular variable. Climate velocities are defined as the ratio of the long-term trend over the mean spatial gradient of a variable (Loarie et al., 2009). Previous oceanic studies have focused on temperature and showed a close relationship between species migration and speed of isotherm migrations (Burrows et al., 2014; Brito-Morales et al., 2020; Jorda et al., 2020; Thorne and Nye, 2021). However, for marine organisms, empirical data indicate that tolerance to hypoxia decreases with increasing temperature (Pörtner and Knust, 2007). This relationship can be modeledFrontiers in Marine Science 02(Deutsch et al., 2015) is a measure of aerobic scope in an oxygen- limited environment. It describes the capability of an environment to sustain aerobic metabolic activity relative to that at rest, and thus accounts for the non-linear interaction between thermal and oxygen stresses on a particular organism. Studies using the metabolic index highlighted a close match between the metabolic index and marine species distributions (Deutsch et al., 2015; Penn et al., 2018; Deutsch et al., 2020; Howard et al., 2020). In fact, the ratio necessary to sustain active aerobic metabolism defines an ecological trait Fcrit which has been found to correspond to the geographical limit of observed viable habitat (Deutsch et al., 2015). Below the threshold Fcrit, species are unable to perform active aerobic metabolism. This is of particular relevance in the SEP OMZ where marine species can undergo both types of stressors (warming and hypoxia) when the OMZ limits experience changes under environmental forcing, whether they originate from natural variability (e.g. El Niño/La Niña) or from anthropogenic activities. The SEP and its OMZ harbor a great number of species, many of them endemic (Clark et al., 2010; Friedlander et al., 2016), located in the vicinity of the OMZ upper limit along the continental slope or in the off-shore ocean where a chain of seamounts system stands along the Nazca, Salas y Gómez, and Juan Fernández ridges (Mecho et al., 2021; Tapia-Guerra et al., 2021). Although the OMZ hosts species that tolerate low oxygen conditions, it is also visited by species performing Diel Vertical Migrations (DVM), such as zooplankton (Wishner et al., 2008; Wishner et al., 2020), seeking punctual refugia from predators, or accumulating at the lower OMZ borders in micro-aerobic conditions attracted by the preserved organic matter sinking out of the OMZ (Paulmier et al., 2021). A great number of seamounts also shelters target species for artisanal fisheries (e.g., the golden crab, the Juan Fernandez lobster and the Juan Fernandez morwong (Arana and Ziller, 1985; Ernst et al., 2013; Friedlander and Gaymer, 2021) being the economic support of local communities. Artisanal fisheries supporting local Chilean and Peruvian economies (Gutiérrez et al., 2016) are also highly vulnerable to climate variability at interannual (e.g. El Niño, Carstensen et al., 2010) and decadal timescales (Bertrand et al., 2011). It is therefore useful to design resource management and adaptation strategies that account for organism responses to multiple stressors present in the SEP. In this context, we propose to investigate future marine habitat change in the SEP (horizon 2100) from a modelling perspective applying climate velocity to the metabolic index (CVF). Considering the scarcity of data to estimate physiological and ecological traits of species inhabiting the SEP (see Deutsch et al., 2015; Deutsch et al., 2020), our approach is to first explore the changes in metabolic index of a model species given climatefrontiersin.org change patterns for temperature and oxygen derived from an metabolic index to temperature and oxygen changes (e.g. differential climate sensitivity), considering current particular a warm bias in the eastern tropical Pacific, CESM match to observations in the tropical Pacific across models (Cabré et al., 2015). Although this model simulates a negative Parouffe et al. 10.3389/fmars.2022.1055875uncertainties in the trajectory of the SEP OMZ in the warmer climate (Oschlies, 2018), in order to explore where changes in habitat are more likely to be temperature or oxygen-driven. While most results are based on a model marine species with a mean ecophysiotype, we contrast them with those of a real species for which sufficient information is available to derive the metabolic index, the shrimp Oplophorus spinosus inhabiting oxygenated waters. Along with the results and sensitivity tests to values of the physiological traits, these are used as material for discussing implications of our work. The paper is structured as follow: Section 2 provides an overview of the dataset and methodology. Section 3 is devoted to the analysis and interpretation of climate velocities of the metabolic index and analysis on the representativeness of the results considering sensitivity analyses of climate velocities to physiological traits and uncertainties due to natural variability. Section 4 is the discussion section where we provide an assessment of exposure to climate change of broad areas based on climate velocities for species inhabiting the SEP in the horizon 2100 and discuss the limits of our approach. This section finally includes concluding remarks and perspectives. 2 Dataset and methods 2.1 Model data In order to diagnose the climate change patterns in the SEP (Figure S1), we use the long-term simulations of the NCAR Community Earth System Model Large Ensemble Project (CESM-LE) (Kay et al., 2015; Lovenduski et al., 2016), which is derived from the Community Earth System Model, version1, with the Community Atmosphere Model, version5 (CESM1 (CAM5); (Hurrell et al., 2013). Although presenting classical biases of current generation global climate models, inEarth System Model (Community Earth System Model Large Ensemble-Project). Secondly, we apply climate velocity (Loarie et al., 2009) to the metabolic index to infer future marine habitat change in the SEP, which, from a methodological point of view, expands former studies (Thorne and Nye, 2021; English et al., 2022) in a multi-factorial environmental stress (i.e. warming and oxygen changes) and for a variable (F) that considers species’ physiology. We show in particular that CVF conveys useful pieces of information to define potential marine metabolic habitat changes, providing in particular an overview on how broad metabolic habitats would evolve in space and time under climate change for a large range of species. Our goal is also to assess the relative sensitivities of climate velocities of theFrontiers in Marine Science 03trend in the global inventory of oxygen over the period 2006- 2100 consistently with other models participating to the Climate Model Intercomparison Project (Phase 6) (Kwiatkowski et al., 2020), it predicts a positive trend in O2 in the OMZ (Figures 1A, B, see also Figure S2 for an estimate of errors due to natural variability), which could be due to enhanced downward mixing of upper oxygenated waters in a warmer climate associated to intensified tropical variability (Carréric et al., 2020) in spite of the increased thermal stratification. Note that a few models participating to the Climate Model Intercomparison Project (Phase 6) also simulate such a positive tendency although disagreeing on its magnitude. Regarding temperature trends (Figures 1C, D), CESM is rather consistent with the ensemble mean of models participating to CMIP6 (Karamperidou et al., 2017; Cai et al., 2021). Besides its overall realism, this model resource is also perfectly suited for assessing uncertainties associated to natural variability (Figure S2), given that it provides simulations that solely differ in terms of initial conditions, thus accounting for fluctuations that originate from the internal dynamics of the earth system (Kay et al., 2015). We used temperature, salinity and oxygen outputs of the 34 members (i.e. simulations with unique initial conditions) covering the period 2006-2100 corresponding to the RCP 8.5 scenario. The oceanic model is the Parallel Ocean Program (POP, Smith et al., 2010) with a grid resolution of 1°x1° and 60 vertical layers. Hereafter, we systematically calculate the ensemble mean of specific quanti t ies (F and CVF , see below) using oxygen, temperature and salinity data for each of the 34 runs for two different mean climate conditions taken as “present” (2006- 2036) and “future” (2070-2100) conditions. 2.2 Metabolic index Ф The metabolic index is a measure of aerobic scope in an oxygen-limited environment. Following (Deutsch et al., 2015; Deutsch et al., 2020), the metabolic index F (non-dimensional) is calculated as the ratio of O2 (pO2 in atm, Figure S3) supplied by the environment to the oxygen demand of an organism atreasonably simulates the tropical variability, including ENSO (Deser et al., 2012; Karamperidou et al., 2017; Dewitte and Takahashi, 2019) and some key aspects of the sensitivity of the carbon cycle to climate variability (Long et al., 2013). While CESM tends to overestimate the volume and spatial extent of the OMZs in the tropical Pacific (Moore et al., 2013), it properly captures the large-scale spatial pattern distribution of dissolved oxygen (Long et al., 2016) and exhibits the closestfrontiersin.org BD 2100 ) for t ines r Parouffe et al. 10.3389/fmars.2022.1055875A C FIGURE 1 Ensemble mean of oxygen and temperature trends between 2006 and lines represent the projected oxygen isocontours (5 and 45 mmol.m-3 2070-2100, solid). For temperature trends (bottom panels), the black l (dashed) and future (plain) climates.Frontiers in Marine Science 04. At 200m (A, C) and 26°S (B, D). For oxygen (top panels), the blue he “present” (mean 2006-2036, dashed) and “future” climates (mean epresent the projected temperature isocontours (15°C) for presentrest, thus representing the coupled effect of oxygen and temperature on metabolic supply and demand (eq. 1 and 2; Deutsch et al., 2020): F = O2supplyO2demand = aS aD BϵpO2e E0 kB 1 T− 1 Tref   (1) Eq. 1 is simplified as: F = A0pO2e E0 kB 1 T− 1 Tref   ≥ Fcrit (2) where kB is the Boltzman constant and A0 and E0 the physiological traits of a given species (see Table 1). aD and aS are the resting metabolic rate and the efficacy of O2 supply per unit of body mass (B) respectively. A0 = aS aD , represents the hypoxic threshold so that 1/A0 is the minimum oxygen that can sustain metabolic activity at rest at a reference temperature Tref of 15°C (Penn et al., 2018; Penn and Deutsch, 2022). E0=ED- ES and represents the difference between the sensitivity to temperature of oxygen demand (ED) and oxygen supply (ES). ϵ is the allometric scaling of the supply to demand ratio. Bϵ is assumed equal to 1 (Penn et al., 2018; Deutsch et al., 2020; Penn and Deutsch, 2022). Below F = 1, oxygen demand exceeds oxygen supply capacity and organisms must suppress aerobic metabolic activities. Species must live in an environment capable of sustaining ecological activities, which is when F is superior to a threshold Fcrit. Fcrit is an ecological threshold as it represents the minimum energy requirements to support ecological activities such as swimming and feeding. In other terms, it is the threshold to maintain a population and as such, delimits geographical viable habitat. A0, E0 and Fcrit define an ecophysiotype unique to each species. From the most updated database (Deutsch et al., 2020), out of 71 species documented, we found estimates of A0, E0 andFcrit for only 6 species present in the SEP: 2 cephalopods (Dosidicus gigas and Octopus vulgaris) and 4 crustacea (Gaussia princeps, Notostomus elegans, Oplophorus gracilirostris and Oplophorus spinosus). As this represents a limited sample, following Penn and Deutsch (2022), we use a mean ecophysiotype: A0 = 20 atm - 1, E0 = 0.4 eV and Fcrit = 3 that is representative of a model species with mean traits of real species populating regions not necessarily in the SEP. A0 and E0 are the mean values of physiological traits estimated for 71 species (62 for Fcrit), benthic and pelagic, of 5 phyla (Annelida, Arthropoda, Chordata, Mollusca, and Cnidaria) in diverse ocean basins and biomes (Deutsch et al., 2020). In addition, we computed F for the shrimp, Oplophorus spinosus (A0 = 23.8 atm -1, E0 = 0.19 eV, Fcrit = 2.5). Metabolic traits for the shrimp were measured from animals off Hawaii (see database in Deutsch et al. (2020)) residing along a depth range of 140-750m with temperatures in the range 5-24°C (Cowles et al., 1991). We assume that thesefrontiersin.org difference in local conditions between adjacent cells. It is always positive and the direction of the gradient can be inferred from migrates or not. In addition, sinceF is linearly dependent on A0, (dT) and oxygen (dpO2) to CVF (Eq.6 and Supplementary Material eq. 4): TABLE 1 Physiological traits (A0 and E0) and Фcrit of the model species and O. spinosus. V) Parouffe et al. 10.3389/fmars.2022.1055875the angle between the latitudinal and longitudinal gradients. Climate velocities therefore represent the rate of change modulated by the spatial distribution of the climate variable. Hence, climate velocities of the Metabolic Index (CVF) indicate where (direction) and how fast (speed) the local and present metabolic index (F in 2006 in any cell (i, j)) is projected to displace by 2100. For example, if the metabolic index decreases (negative velocities) at some location, the species would need to move in the direction of increasing Ф, that is in the direction oftraits values are valid for the shrimps of the SEP as both environments are comparable over such a depth range. Lastly, this shrimp is found in oxygenated waters and near oxygen minimum layers and therefore may be used to represent how normoxic species may be impacted by warming and deoxygenation. Thus, we argue that the case of the shrimp can be used as a benchmark for comparing with the results obtained with the model species. 2.3 Climate velocities 2.3.1 Definition Climate velocities describe the speed and direction at which local climate conditions are changing (Loarie et al., 2009). Climate velocities (CV, km.yr-1) were calculated using the gradient based approach (eq. 3 and Supplementary Material eq. 3): CVF km : yr −1   = temporal trends  F : yr−1ð Þ spatial gradient  F : km−1ð Þ (3) The temporal gradient (F.yr-1) represents the local rate of change and corresponds here to the 2006-2100 local linear trend for F in every cell (i,j). The mean spatial gradient of F is calculated as the square root of the sum of the squared latitudinal and longitudinal gradients of F. The spatial gradients (F.km-1) were calculated using the same numerical scheme as that of the model code (Smith et al., 2010) so as to avoid artifacts associated to interpolation procedures. The spatial gradient indicates the Species A0 (atm -1) E0 (e Model 20 0.4 O. spinosus (shrimp) 23.8 0.19Frontiers in Marine Science 05CVFapprox = dFapprox ∇F = A0e E0 kB 1 T− 1 Tref     pO2 : −E0kBT2   : dT+dpO2 h i ∇F   (6) 2.4 Sensitivity analysis We perform several sensitivity analyses of F and CVF to physiological traits and natural variability.its gradient as well, so that CVF is independent of A0 (see Supplementary Text) and therefore is weakly sensitive to the species-specific F. This confers to CVF an added-value as a metrics for inferring change in metabolic habitat. 2.3.2 Drivers of climate velocities of the metabolic index We evaluate the contribution of oxygen and temperature to CVF using the differential of Ф as follows (eq. 4 and 5): dFapprox = ∂F∂T : dT + ∂F ∂ pO2 : dpO2 (4) dFapprox = A0e E0 kB 1 T− 1 Tref       pO2 : −E0kBT2   : dT + dpO2 h i (5) where dT and dpO2 are the long-term trends in temperature and oxygen and where T and pO2 are the mean (2006-2100) temperature and oxygen. dFapprox thus represents the linear approximation of the changes in F between the two climates. If we further divide it by the mean gradient of Ф ( ∇F ), we can access a linear estimate of CVF (CVF,approx) that informs on the relative contribution of the long-term trends in temperaturenot be sensitive to CVF. However CVF still informs on the distance between F and Fcrit in the future climate, which conveys information on the extent of the changes in metabolic habitat for that particular species independently of whether it Фcrit Reference 3 Deutsch et al., 2020; Penn and Deutsch, 2022 (Fcrit) 2.5 Deutsch et al., 2020frontiersin.orgthe gradient. Therefore, CVF indicates where species would migrate (Brito-Morales et al., 2018) in 2100 if they were tracking their present metabolic index (F in 2006). In practice, there are many instances when the value of F of some species is sufficiently larger than the value ofFcrit so that these species may 2.4.1 Sensitivity analysis on F and CVF to A0 and E0 Sensitivity of F to the traits is assessed using values of A0 in the range 5 to 60 atm-1 (A0 range, Supplementary S4A) and values of E0 in the range -0.2 to 1.3 eV (E0 range, Supplementary sensitivity to E0 is assessed since CVF is independent of A0. Climate velocities are calculated with F(E0) (hereafter CV(E0)) and vertical dimensions allows identifying changes in key areas documented range of traits (Figure S9). Therefore, given the Parouffe et al. 10.3389/fmars.2022.1055875and the associated standard deviation is then derived (Error_CV (E0)). The mean velocity of CV(E0) (hereafter CV(E0)mean) is compared to CVF of the model species (hereafter CVF,model) to emphasize the statistical significance of the 2D distribution of CVF,model in relation to CV(E0)mean. 2.4.2 Sensitivity of CVF to natural variability The CESM-LE offers the opportunity to estimate uncertainties associated to natural variability that accounts for variations in the long-term mean associated with natural modes of variability in the SEP such as ENSO or the Interdecadal Pacific Oscillation (Henley et al., 2015). This type of uncertainty can be derived from the estimate of the standard deviation of CVF amongst the 34 realizations (i.e. members) by the same model. We estimated the uncertainty associated to natural variability, oxygen trends and temperature in order to evaluate where they would most influence CVF. We also provide the ratio of these uncertainties relative to absolute value of CVF, which highlights themagnitude of the effect of natural variability on CVFmodel and thus on the confidence level of CVF. Uncertainty was calculated as the standard deviation of CVF amongst the 34 members, that is:ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi oNi=1 ½CVi∅−CV∅2 N q   where N=34 (hereafter referred to as Error_CVF_all). We also estimate uncertainties in CVF associated to that of the oxygen trend, i.e. oxygen trend-induced uncertainty (hereafter Error_CVF_O2) and the uncertainties in CVF associated to that of the temperature trend, i.e. temperature trend-induced uncertainty (hereafter Error_CVF_T). They (Error_CVF_T and Error_CVF_O2) correspond to the standard deviation of CVF amongst the 34 members when either oxygen and temperature is replaced by its ensemble mean. The relative error of CVF, model is also provided (i.e. the ratio of the absolute error to CVF, model). These estimates allow us to gain confidence in CVF with regards to the effects of natural variability. 3 Results We will focus our analysis along a vertical section at 26°S approximately aligned with the Nazca and Salas y Gomez ridges and for a horizontal section over the SEP at 200m (average depth of the upper oxycline). The analysis focused on the horizontalS4B). Hence F is calculated through the A0 range with E0 kept constant at 0.4 eV (hereafter F(A0)) and through the E0 range with A0 kept constant at 20 atm -1 (hereafter F(E0)). The sensitivity is then evaluated as the standard deviation across F (A0) (herea f ter Error_F (A0)) and across F (E0) (hereafter Error_F(E0)). The same is done for climate velocities, except that only theFrontiers in Marine Science 06lack of observed values in (A0, E0), mean Fcrit (=3) cannot be used here with confidence to project mean changes in viable metabolic habitat in the SEP. We show hereafter that CVF hassuch as: The OMZ mainly characterized by a reoxygenation, deep oxygenated waters characterized by a deoxygenation and the upper waters characterized by warming (Figure 1). 3.1 Changes in metabolic index of the model species between the present and future climates We present as a first step the mean metabolic index of the model species (Фmodel) for the future climate (2070-2100) and the associated change in F with respect to the present climate (Figures 2A, B). At 200m, values range from 0 to 5 (Figure 2A) and from 0 to 3 for the vertical section (Figure 2B). The distribution of F clearly reflects the influence of the mean oxygen distribution. The OMZ domain (see white contour Figures 2A, C) is non-viable (Fmodel< 1). A strong gradient is found along this oxygen isocontour marking the transition from hypoxic (45 mmol.m-3 white contour in Figure 2) to oxygenated waters (see Figure S5). Fmodel increases towards the south and above and below the OMZ as oxygen concentrations increase. However, Fmodel can be lower in surface waters than in deeper waters as a result of higher temperatures (Figure S5). We represented Fcrit (= 3) on Figure 2A (see dashed contours), which delimits the geographical viable habitat of the model species, for the two climates. The model species occupies well oxygenated waters (>200 mmol.m-3 Figure S5). A maximal 10% decrease causes a shift of viable habitat of approximately 5° south (~500km, Figures 2A, B) on the western side of the basin and a horizontal compression below 200m (Figures 2C, D) between 2006 and 2100 (Figure 2B) caused by the projected deoxygenation (Figures 1A, B). No gain of mean potential habitat between the present and future climates is simulated over the OMZ despite a consistent simulated reoxygenation in the upper OMZ (Figure 1B). In the OMZ, the distance between F and Fcrit remains significant for the model species. These small decreases inF cause a projected 3D reduction of habitat of 9.5% with a stronger loss in the upper 200m than below 200m Figure S6A). Changes in viable habitat are also highly dependent on whether Fcrit is positioned in a high gradient area or not. Further sensitivity tests show a high dependency of the magnitude of F to the trait A0 due to the linear relationship betweenF and A0 (Figures S7A, C) but a weaker sensitivity to E0 (Figures S7B, D). This can be illustrated with the shrimp O. spinosus (Figure S8) which shows higher values of F with a ratio Fmodel/Fshrimp ≈ A0,model/A0,shrimp. This also demonstrates that the position of the mean Fcrit strongly varies within thefrontiersin.org Parouffe et al. 10.3389/fmars.2022.1055875Frontiers in Marine Science 07much less dependence to A0 and E0, which makes it more convenient for addressing this issue given the scarcity of data in the SEP. 3.2 Climate velocities of the metabolic index (CVФ) 3.2.1 Pattern of CVФ We estimate here the speed and direction of the future displacement of the metabolic habitat described by F (see section 2.3) for the model species (Figure 3). The pattern of the spatial gradient of Fmodel indicates that it is constrained by pO2 (see eq.4). More specifically, it has bimodal spatial distribution constrained by the OMZ (Figure 3B, 45 mmol.m-3 blue isoline) and by oxygenated waters (Figure 3B, 150 mml.m-3 blue isoline). On the one hand, the OMZ is characterized by a high horizontal gradient area where F changes are in the higher range of 0.5 F/100km (in the outer edges of the OMZ) to 1.5 F/ 100km (between 30°S and 40°S) for the model species. This area characterizes a metabolic barrier as F varies abruptly there, separating two very distinct environments, i.e. separating hypoxic species (Vaquer-Sunyer and Duarte, 2008) and normoxic species. Juan Fernandez Archipelago (JFA) and the Desventuradas Islands Archipelago (DA) are within the area of high gradient with respectively 0.5 and 1.5 F/100km. On the other hand, oxygenated waters are characterized by a very weak gradient in the range of 0 to 0.5 F/100km with a mean spatial gradient there of ~0.1 F/100km indicating that the metabolic environment is stable and similar over long distances. The trends (Figure 3C) are positive on the eastern side of the basin and negative in the open ocean, covering the majority of the SEP, consistent with that of the oxygen trends at this depth (Figure 1). Positive trends represent locally an increase of up to 0.8 units (Figure 3C) by 2100 where F is lowest (F<2, FIGURE 2 Metabolic index of the model species (Фmodel). Mean metabolic Index (A, C) and change (DF) between the present and future climate (B, D) at 200m (A, B) and 26°S (C, D).The dashed isocontours indicate the position of the mean Fcrit (=3) in the present (black) and the future (red) climates. The white line is the mean 45 mmol.m-3 oxygen isocontour of the future climate. The dark grey shading is where F< 1. The dark stars represent Juan Fernandez (JFA), the Desventuradas Archipelagos (DA) and Easter Island (EI). The black dots (B, D) indicate where DF is not significant across the 34 members (Wilcoxon rank sum test p < 0.05).frontiersin.org c Inde patial ed wa hipela Parouffe et al. 10.3389/fmars.2022.1055875A B FIGURE 3 Climate velocities, spatial gradient and temporal trend of the Metaboli velocities (A), the dashed line is where CVFmodel = 0 km.yr -1. For the s 1 km.yr-1. The blue isolines represent the 45 (OMZ) and 150 (oxygenat The dark stars represent Juan Fernandez (JFA), the Desventuradas ArcFigures 2A, C). West of 80°W negative trends represent a Frontiers in Marine Science 08C x F of the model species at 200m over 2006-2100. For climate gradient (B), the dotted line delimits the region for which |CV|< ters) oxygen isocontours in mmol.m-3. (C) are the temporal trend. gos (DA) and Easter Island (EI).3.2.2 Drivers of CVF decrease of up to 0.1 units and up to 0.5 F by 2100. Trends in the key ecoregions of JFA and DA are weak. CVF,model are positive east of the basin consistent with temporal trends. They are mostly weak (i.e. |CVF|<1 km.yr -1, Figure 3B dotted line envelope) within the OMZ. There, the displacement of habitat represented by CVF,model is controlled by the spatial distribution ofFmodel (spatial gradient). In areas of slow positive/negative velocities, analogue metabolic environments are found “close” (within 100km) and in the direction of an increasing/decreasing gradient (Figure S10). This translates into a short displacement of species to find metabolic conditions analogous to those of their occupied habitat. CVF,model is locally high (>10 km.yr -1) and positive in the coastal region of Chile between 30°S and 40°S and off Peru at 20°S as the result of the combination of high temporal trends and a weak spatial gradient. There, marine species will benefit from a fast temporal F increase (high trend) in a stable environment (weak spatial gradient). Oxygenated waters are characterized by relatively high negative velocities (5km.yr-1) due to temporal trends in the higher range and weak gradients. Faster negative velocities (< -10 km.yr-1) are observed in the southernmost part of the basin and adjacent to the OMZ where the most negative trends are found. We now evaluate the relative contributions of temperature and oxygen trends to CVF (see section 2.3.2) in order to characterize the drivers behind the projected habitat dynamic and given current uncertainties in oxygen trend in particular to gain confidence in the assessment of areas of climate exposure (section 4.1). Figure 4 represents the relative contribution of pO2 and temperature to CVF. A negative sign implies a CVF decrease and conversely (in %, see section 2.3.2 and Supplementary Text). There are very contrasting behaviors between the OMZ and oxygenated waters and between the subsurface and deeper waters. At 200m (Figures 4A, B), positive oxygen trends are the main driver of CVF in the OMZ with an incursion into equatorial waters (>75%). Negative oxygen trends contribute at least to 50% of the change in CVF in the northernmost part of the SEP (south of 20°S) and along the OMZ border. On the contrary, temperature trends are the main driver in the south-end of the SEP with a contribution superior to 75%. At 26°S (Figures 4C, D), in the OMZ core, the trends in oxygen is the main driver (> 75%) and contributes to an increase in CVF. In oxygenated waters however ([02] > 150 mmol.m -3), oxygen trends contribute more (>50%) than temperature trends below 400m whereasfrontiersin.org BD f the ative/ % = 5 Parouffe et al. 10.3389/fmars.2022.1055875FIGURE 4 Contributions (in %) of the trends in oxygen and temperature to CVФ o 2.3.2 is used for this calculation. A negative/positive sign means a neg driver when % > 75, and O2 and temperature contribute equally when transition area between positive and negative velocities.A Ctemperature trends are the main driver in the thermocline (> Frontiers in Marine Science 09model species. At 200m (A, B) and 26°S (C, D). The Eq. 6 of section positive contribution of the driver. A driver is considered a main 0. The magenta contour is where [CV]= 0. km.yr-1 signaling the3.3.1 Representativity of the model species 75% above ~400m). There is an approximately equal contribution of both drivers to CVF around the transition zone between positive and negative oxygen trends (magenta contour Figure 4). CVF undergo either the additive negative effects of warming and deoxygenation, or, the compensatory effect of positive trends in oxygen and warming. Note that the effect of warming here will be modulated by E0 greatly affecting the sign and magnitude of the velocities. 3.3 Sensitivity analysis In this section we analyze how CVF varies as a function of E0 and natural variability. This allows us to assess the confidence level one can have in CVF,model to infer future displacement of local metabolic environment. We show in particular that CVF,model is statistically comparable to CV(E0) mean and an ensemble of real species. We also identify regions where natural variability is likely to contribute to the uncerta int ies in the est imate in CVF and habitat changes. This will serve as material for guiding our method to quantify changes in habitat and for discussing our results (section 4.1). In this section, we explore the sensitivity of CVF to the physiological trait E0 in order to infer if our results can be generalized to the case of real species. We first analyze the uncertainty associated to CVF across the E0 range (Figure 5). Error_CV(E0) (Figures 5A, C, see method section 2.4.1) indicates the weakest sensitivity of CV(E0) (s.d. < 1 km.yr -1) in the OMZ domain (white contours, Figures 5A, C) and a relatively stronger sensitivity (s.d. > 1 km.yr-1) in the thermocline (Figure 5C) beyond the OMZ limits ([O2] > 150 mmol.m -3, Figures 5A, C). Furthermore, CVF,model, CV(E0)Eo>0 and CV(E0)mean show statistically comparable bimodal distributions (Figure 6). CVF, model distribution slightly shifts to more negative velocities for CVF,model compared to CV(E0)mean. This can be verified with the velocities of the shrimp O. spinosus (Figure S11) whose velocities mirror that of CVF,model despite a weaker E0 (E0,shrimp = 0.19 eV vs. E0model = 0.4 eV). CV(E0) with negative E0 (light grey lines, Figure 6) are however not comparable as negative E0 tend to increase F under warming trends (see also Figure S12 and section 3.2.2) compensating deoxygenation. Lastly, the ratio CV (E0)mean/CVF,model is weak (<10%) in all the SEP except in the transition area (Figures 5B, D). In fact, the highest uncertainty occurs where the contribution of oxygen and temperature to CVF are equal (i.e. OMZ border and transition between positive and negativefrontiersin.org velocities, see magenta contour Figure 4). There, E0 affects the uncertainty are in fact found where CVF is relatively weak (i.e. gradient and temporal trends to illustrate how CVF can help define broad areas of exposure to climate change. Our approach pressure species face at any location. The magnitude of the velocity Parouffe et al. 10.3389/fmars.2022.1055875near balance between the contributions of oxygen and temperature to CVF so that any changes caused by E0 affect significantly the sign and magnitude of CVF, especially in this low velocity area. For example, increasing E0 strengthens the effect of warming over that of positive oxygen trends so that CVF tend to become positive, and conversely. Varying values of E0 cause ultimately a shift in the median of CVF (Figure S12) and hence displace the frontiers between positive and negative velocities. This effect is also visible when comparing CV(E0) mean to CVF,model where the relative difference in magnitude is >25% (Figures 5B, D). Hence, estimating the speed of habitat displacement with CVF is more trait sensitive in the transition area where the change in sign of velocities occurs. Hence, we argue that CVF,model is quite representative of the mean displacement of an ensemble of real species, except in the transition area highly sensitive to E0 because of the way it affects the equal contribution of oxygen and temperature trends to changes in F. This is particularly valuable and useful as very few species ecophysiotypes are documented. 3.3.2 Uncertainty in CVF due to natural variability The CESM-LE offers the opportunity to estimate uncertainties associated to natural variability that accounts for variations in the long-term mean associated with natural modes of variability in the SEP such as ENSO or the Interdecadal Pacific Oscillation (Henley et al., 2015). This type of uncertainty can be derived from the estimate of the standard deviation of CVF amongst the 34 realizations (i.e. members) of the same model. We estimated the absolute (Figures 7A–C and Figures S13A–C) and relative (Figures 7D–F; Figures S13D–F) uncertainty associated to natural variability, and whether it originates from the oxygen or temperature fields (see method in section 2.4.2) in order to provide guidance for determining areas of climate change exposure (see section 4.1). The results indicate a vast share of the domain where natural variability can induce fluctuations in the long-term mean of temperature and oxygen and thus CVF (Figures 7A–C; Figures S13A–C). Error_CVF_all appears to be mostly related to Error_CVF_O2 within the OMZ and slightly beyond the OMZ limit. Over the rest of the domain, Error_CVF_T is comparatively strong in oxygenated waters above 200m (Figure 7C; Figure S13C). Where interannual variations in oxygen and temperature are large, such as along the coast of South America where there is an efficient ENSO oceanic teleconnection (Sprintall et al., 2020) and over the subtropical pressure system where tropical and extra-tropical atmospheric teleconnections modulate the oceanic circulation, large uncertainties are expected. However since CVF tend to have large values (in absolute value) there, the relative uncertainty remains weak (Figures 7D–F; Figures S13D–F). The regions of large relativeFrontiers in Marine Science 10is used based on an arbitrary threshold equal to 1 km.yr-1 (Figure 8). This threshold value is convenient because it allows visualizing where the spatial gradient in F is larger or smaller than its long- term trend, i.e. where the local rate of change is larger than the spatial distribution of F. Therefore, slow velocities indicate that local changes can be “compensated” by close analog environments, whereas fast velocities indicate that local analogues will be found >100 km by 2100. Hence, we define four areas of incremental risk (Figure 8). Area I is the safest with fast positive velocities and therefore a fast increase in metabolic index locally. Area II isis based on the choice of threshold values for each of these three quantities. This choice is here somewhat subjective since we work with a model species, but it could be adjusted in case studies. We define areas of exposure to climate change based on the sign and magnitude of the velocities (Figure 8). Positive trends indicate a favorable evolution of the local metabolic environment (maintained to increased aerobic scope). However, negative trends indicate a reduction of aerobic scope and therefore an increased probability for F to fall below Fcrit. Hence local trends indicate the level of|CV|< 1 km.yr-1, see dotted contours Figure 3B), which corresponds to a zone approximately bordering the OMZ limits in its interior to the South of 20°S and expanding to the Peru upwelling system and equatorial region to the North. We define this area as the “transition area” (magenta contour Figures 7D; S13D) where uncertainty due to natural variability limits our confidence in the estimate of CVF. That said, confidence in estimate of CVF remains strong in oxygenated waters and the OMZ and can be relied on to infer future habitat change. 4 Discussion While our results illustrate the potential in using CVF for inferring changes in habitat, the characteristics of the SEP region (i.e. low consensus in deoxygenating trend in the OMZ region amongst climate models, and scarcity of physiological data) limit a quantitative inference of future habitat. Still, we argue that CVF can be useful for inferring broad regions where metabolic conditions change and therefore of potential metabolic refugia. 4.1 Inferring change in habitat from CVF: Areas of climate change exposure Taking advantage of the model species representativity, we now combine the information given by the velocities, spatialfrontiersin.org Parouffe et al. 10.3389/fmars.2022.1055875relative error exceeds 25%. This area encompasses the changes in Frontiers in Marine Science 11moderately safe; the metabolic environment improves at a slow rate and local analogue environments are found “close”. Area III is moderately at risk with a decrease of F but with close analogue environments. Finally, Area IV is most at risk with the fastest negative velocities. Locally, themetabolic environment degrades at a fast rate suggesting a short emergence time for negative conditions (e.g. Fcrit). Lastly, we define an area of high uncertainty which includes the uncertainty associated to traits (Figures 5B, D) and associated to natural variability (Figures 7D, S13D) and where the the contribution of drivers to CVFmodel, and allows to capture the coupled effect of temperature and E0 (Figure 5) and the change in sign of oxygen trends (Figure 7). Based on these definitions, we illustrate the potential consequences for marine species and habitat in the OMZ and its border, the open ocean and seamounts (Figure 9A) and along a vertical section at 26°S (Figure 9B, see also Figure S14 for a vertical section of the horizontal velocities at 26°S) which allows to visualize potential risks for the 2D horizontal marine FIGURE 5 Sensitivity of CVF to the physiological trait E0. At 200m (top) and 26°S (bottom). (A, C) provide the absolute uncertainty and (B, D) the relative uncertainty (see method section 2.4.1). The range of E0 used is provided in Figure S4. (B, D) represent the relative difference between CV(E0)mean and CVF,model in %. The white contours indicate the 45 and 150 oxygen isocontours (mmol.m -3). The magenta lines indicate where CVФ,model = 0 km.yr-1 (i.e. the transition area).metabolic habitats at depth. frontiersin.org Areas I and II spread along coastal waters and over the OMZ spreads across the OMZ border and represents a horizontal Parouffe et al. 10.3389/fmars.2022.1055875domain at 200m and in the OMZ core (Figures 9A, B). Consequently, species living and accumulating at the lower OMZ borders (Paulmier et al., 2021) or demersal or benthic (Levin and Gallo, 2019), may find less stressful and a more stable environment which should preserve species abundance. Area IIFrontiers in Marine Science 12compression of the OMZ (defined by the 45 mmol.m-3 isocontour) between 20°S and 40°S at 200m. There, hypoxia tolerant species inhabiting the outer edges of the OMZ (Birk et al., 2019; Roman et al., 2019) may extend their habitat distribution. Area I and II are further characterized by veryFIGURE 7 Uncertainty of CVF associated to natural variability at 200m. Error_CVF_all, Error_CVF_O2 and Error_CVF_T correspond to the uncertainty associated to natural variability and the oxygen-trend and temperature-trend induced uncertainty, respectively. (A–C) provide the absolute error and (D–F) the relative error (%, see method section 2.4.2). The white contours are where [O2] = 45 and 150 mmol.m -3. The magenta dashed contours indicate the transition area (|CVFmodel|= 0 km.yr -1).FIGURE 6 Pdf of climate velocities of the metabolic index at 200m. For CVFmodel (black) and CV(E0)mean (red). For the 71 species in grey with positive E0 (dark grey) and negative E0 (light grey).frontiersin.org nega nes I ncer en (O Parouffe et al. 10.3389/fmars.2022.1055875low F (F<2, Figure 2B) and by species adapted to hypoxic environment. Therefore, it is unlikely that hypoxia intolerant species migrate to area II. Since changes are mainly driven by oxygen (Figures 4A, C), species sensitive to temperature are more likely to sustain the weak changes in F attributed to increasing temperature. Hence, species composition and richness should be preserved. Area III is characterized by slow negative velocity and is located along the OMZ border. In area III, the temporal trends are in the lower range, meaning a slower decrease of F and therefore a limited compression of their habitat. Besides, being in a high gradient area, species would find analogue conditions a few kilometers away in the direction of the open ocean with more oxygenated waters. However, this needs to be taken with caution, as what matters for a species is the distance between F and its respectiveFcrit. Indeed, slow negative velocities could still be fast enough to cause deleterious changes as many species live close to their physiological limits (Tremblay et al., 2020; Wishner FIGURE 8 Definition of areas of climate change exposure. Areas with positive or area of slow velocities driven by the spatial gradient (|CV|< 1km.yr-1, zo CV| > 1 km.yr-1, zones I and IV). The transition area defined by a high u as uncertain. To each of these areas correspond main drivers (i.e. oxygviable as long as F remains above Fcrit. Equatorial waters and Frontiers in Marine Science 13coastal Peruvian waters are located in areas III and IV which may drive fisheries species away to the south and from the coasts. According to our results, artisanal and industrial fisheries as well as aquacultures from coastal waters off Peru and Chile, mostly located in areas III and IV, are at risk. Lastly, the transition area (Figure 9, dashed envelope) encompasses seamounts areas such as DA and JFA and Exclusive Economic Zones where fishing exploitations occur. These seamounts (Yanez et al., 2009) are known to act as stepping stones for species dispersion (i.e., Salas y Gomez Ridge), by diverse and highly vulnerable ecosystems and by the highest number of endemic species (Yanez et al., 2009; Dyer and Westneat, 2010; Friedlander et al., 2016; Mecho et al., 2021). Islanders’ local economy is based on the artisanal fisheries (Friedlander and Gaymer, 2021), such as the Juan Fernandez lobster. Environmental changes could be harmful to these communities if marine species were to migrate or perish. Besides, fishing grounds could be displaced beyond reasonable distances out of reach for these communities which use tive velocities are classified as “safe or “at risk”. We also distinguish I and III) and areas of fast velocities driven by the temporal gradient (| tainty related to natural variability and physiological traits is classified 2) and temperature (T) as analyzed in section 3.2.2.et al., 2020) and would therefore be weakly tolerant to environmental changes. Species at different ontogenetic stages are not responding equally against environmental changes. Organisms at embryonic, larval and spawning stages are less tolerant to elevated temperatures (Dahlke et al., 2020) and low oxygen conditions (Breitburg, 2002) than adults in marine environments. This is especially important as changes in area III are driven by both temperature increase and oxygen loss. Area IV covers the west of the SEP (Figure 9A) and most of the vertical range at depth (Figure 9B), from shallow waters above the OMZ to all the water column to the west. In area IV, changes are also driven by both drivers. Area IV is all the more at risk as analogue environments will be found over long distances (>100km away by 2100), suggesting long southward migrations (Figure S10). However, area IV is also characterized by the highest F, so that habitats may remain small boats and are not capable to sail long distances. Additionally, these islands are characterized by a high level of connectivity at larval stage. Changes there could reduce the migratory fluxes at the larval stage (Porobic et al., 2012; Porobic et al., 2013) of several species such as the sea urchin Centrostephanus sylviae (Veliz et al., 2021) or the Juan Fernandez lobster Jasus frontalis (Ernst et al., 2013), a fishery target species. Making use of species traits in these areas would reduce the uncertainty in the sign of the changes in these areas. 4.2 Cautionary consideration for the use of the metabolic index While the metabolic index can be used as an indicator of exposure to climate change, it is important to discuss its relevance for the OMZ waters and for some species. When environmentalpO2 is limiting the maximum metabolic rate (MMR), F is frontiersin.org P at he co s III a radien ainty tour Parouffe et al. 10.3389/fmars.2022.1055875equivalent to the measure of the factorial aerobic scope (FAS; (Seibel and Deutsch, 2020; Seibel and Birk, 2022) which is the factorial difference between maximummetabolic rate (MMR) and standard metabolic rate (SMR). The majority of species that have A B FIGURE 9 Areas of climate change exposure. Areas of climate exposure in the SE horizontal velocities. We defined areas of climate exposure based on t Areas I and II designate “safe” areas based on positive trends in Ф. Area and III are areas of slow velocities (|CV|< 1 km.yr-1) where the spatial g transition area (dashed envelope) is the area most subjected to uncert the sign of oxygen trends and the effect of E0. The black and solid con contour also coincides with F = 1.higher or lower than ambient pO2 at depth, but it will always be Frontiers in Marine Science 14much less than air-saturation in oxygen minimum zone species. When Pcmax < pO2, the metabolic index overestimates aerobic scope, which is, instead, dependent only on temperature. Species can migrate to keep within the boundary conditions defining the 200m (A) and the vertical representation at 26°S (B) derived from the mponents of the climate velocities of Ф using the model species. nd IV designate “at risk” areas based on negative trends in Ф. Areas II t controls the magnitude of CVF over the temporal trends. The in the sign and magnitude of the velocity due to high uncertainty in represents the OMZ limit at 45 mmol.m-3. Note that the OMZbeen studied in the literature have evolved in shallow, coastal waters where pO2 is routinely near air-saturation, regardless of temperature. The oxygen supply capacity for each of those species has evolved to meet maximum demand at the prevailing pO2, which for most of them is near air-saturation. As a result, the critical oxygen pressure for maximum metabolic rate (Pcmax) is also near air-saturation and both, MMR and FAS, decline in proportion to pO2 below Pcmax, again regardless of temperature (Seibel and Deutsch, 2020). For such species, the metabolic index (simplified as pO2/Pcrit at any given temperature), is a valuable estimate of the aerobic scope for activity. For vertical migrators, and likely other mesopelagic species in the upper water column, pO2, temperature and light are strongly correlated along a depth gradient. Metabolic demand for activity decreases with light and depth (Childress, 1995), such that MMR and FAS are highest in warm water and decline strongly with depth and temperature. Pcmax is also highly sensitive to temperature and tends to decline with depth at a similar rate as pO2 below the thermocline (Seibel and Birk, 2022). At the lower oxycline, oxygen begins to increase while temperature continues to decline with depth. For species living there, E0 is strongly negative and hypoxia tolerance is improved in warm water (Wishner et al., 2018). Pcmax is not readily predictable for those species. Therefore Pcmax may be ecological niches or can persist locally in new or changing habitats through both phenotypic plasticity and adaptive evolution (Anderson et al., 2012; Munday, 2014) as a response to rapidly changing environmental conditions (Baltar et al., 2019). Metabolic traits of SEP mesopelagic species have changed as a response to a unique environment and thus may not be captured by a mean ecophysiotype model representative of a broad range of species. Thus, Pcmax must be known for each non-coastal or shallow species before the metabolic index can be applied with confidence. 4.3 Concluding remarks We showed here that CVF conveys useful information to define potential marine habitat changes and distribution shifts, providing an overview on how broad environmental niches would evolve in space and time under climate change for a broad range of species inhabiting the SEP (see Figure 9). We also manage to identify regions where the confidence level in the estimate of CVF is likely to be the largest considering current regional uncertainties in deoxygenation trends in the SEP as simulated by current-generation climate models andconsidering the sensitivity of CVF to the physiological traits frontiersin.org E . Using species specific ecophysiotypes can additionally velocities and trajectories of critical habitat which inform on the displacement of boundaries of the metabolic viable have a lower order of magnitude (m.dec-1, Jorda et al., 2020) and they acknowledge support from ANID (Concurso de FutureMares project (Theme LC-CLA-06-2019, Grant We would like to thank Dr. K. von Schuckmann for early The authors declare that the research was conducted in the Parouffe et al. 10.3389/fmars.2022.1055875and therefore, are not directly comparable to horizontal velocities, placing them beyond the scope of this paper. Lastly, recent findings (Seibel and Birk, 2022) suggest that additional physiological information (Pcmax) may be necessary to determine whether the metabolic index can accurately be applied to vertical migrators and oceanic species (see section 4.2). In addition to the metabolic index, other indicators such as the Aerobic Growth Index (AGI, Clarke et al., 2021b) or simply the FAS (Seibel and Birk, 2022) may be valuable alternatives to consider. These represent directions for future investigations. Data availability statement Publicly available datasets were analyzed in this study. Thisdata can be found here: https://github.com/aparouffe/SEP_ habitability. The datasets used in this study are available at: https://www.cesm.ucar.edu/projects/community-projects/ DPLE/datasets.html The scripts for the climate velocity computations areavailable in a GitHub repository: https:// github.com/aparouffe/SEP_habitability. Author contributions VG and BD initiated the study. APar performed the computations and wrote the manuscript. All authors discussed and interpreted the results and provided input to the manuscript. All authors contributed to the article and approved the submitted version.habitat which varies from one species to another. Using CVF with species specific ecophysiotypes would also help reduce the uncertainty of CVF in key areas of the SEP. Other addit ions to the tools for designing MPA include information coming from vertical climate velocities. Indeed, vertical velocities for species performing vertical migrations or species occupying a range of depths would provide useful information regarding the likelihood of vertical migrations especially when vertical habitat seems at risk as suggested by our preliminary results. It was shown that vertical velocities0 support CVF for guiding management approaches to conservation (Brito-Morales et al., 2018). Indeed, defining the distance between the local F (at one location) and the limit of viable habitat (Fcrit) would enable: 1) emergence time estimations of negative conditions (F