# Imprint of intrinsic ocean variability on decadal trends of regional sea level and ocean heat content using synthetic profiles

William Llovel, Nicolas Kolodziejczyk, Sally Close, Thierry Penduff, Jean-Marc Molines, Laurent Terray

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

## To cite this version:

William Llovel, Nicolas Kolodziejczyk, Sally Close, Thierry Penduff, Jean-Marc Molines, et al. Imprint of intrinsic ocean variability on decadal trends of regional sea level and ocean heat content using synthetic profiles. Environmental Research Letters, 2022, 17, 10.1088/1748-9326/ac5f93. insu-03665078

HAL Id: insu-03665078

https://insu.hal.science/insu-03665078

Submitted on 11 May 2022

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
LETTER • OPEN ACCESS

# Imprint of intrinsic ocean variability on decadal trends of regional sea level and ocean heat content using synthetic profiles

To cite this article: William Llovel et al 2022 *Environ. Res. Lett.* **17** 044063

View the [article online](#) for updates and enhancements.

## You may also like

- Surface magnetization thermal fluctuation driven anomalous behaviour of ordinary Hall effect in Pt/YIG
Lili Lang, Shiming Zhou and Xuepeng Qiu

- ²⁷Al MAS NMR, EXAFS and X-ray induced absorption studies in Al-doped and Yb-Al-doped silica glasses fabricated by MCVD
Tomoya Okazaki, Shuichi Arakawa, Edson Haruhico Sekiya et al.

- Observed southern upper-ocean warming over 2005–2014 and associated mechanisms
William Llovel and Laurent Terray
ENVIRONMENTAL RESEARCH
LETTERS

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

CrossMark

OPEN ACCESS

RECEIVED

8 November 2021

REVISED

11 March 2022

ACCEPTED FOR PUBLICATION

21 March 2022

PUBLISHED

5 April 2022

Original content from
this work may be used
under the terms of the
Creative Commons
Attribution 4.0 licence.

Any further distribution
of this work must
maintain attribution to
the author(s) and the title
of the work, journal
citation and DOI.

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

LETTER

# Imprint of intrinsic ocean variability on decadal trends of regional sea level and ocean heat content using synthetic profiles

William Llovel<fnref n="1" />,*, Nicolas Kolodziejczyk<fnref n="1" />, Sally Close<fnref n="1" />, Thierry Penduff<fnref n="2" />, Jean-Marc Molines<fnref n="2" /> and Laurent Terray<fnref n="3" />

E-mail: wllovel@ifremer.fr

**Keywords:** ocean heat content, thermosteric sea level, Argo float, large ensemble simulations, atmospherically-forced variability, chaotic intrinsic variability, decadal trends

Supplementary material for this article is available online

Abstract

The global ocean is warming and has absorbed 90% of the Earth Energy Imbalance over 2010–2018 leading to global mean sea level rise. Both ocean heat content (OHC) and sea level trends show large regional deviations from their global means. Both quantities have been estimated from *in-situ* observations for years. However, *in-situ* profile coverage is spatially uneven, leading to uncertainties when assessing both OHC and sea level trends, especially at regional scale. Recently, a new possible driver of regional sea level and OHC trends has been highlighted using eddy-permitting ensemble ocean simulations over multiple decades: non-linear ocean processes produce chaotic fluctuations, which yield random contributions to regional decadal OHC and sea level trends. *In-situ* measurements capture a combination of the atmospherically-forced response and this intrinsic ocean variability. It is therefore important to understand the imprint of the chaotic ocean variability recorded by the *in-situ* measurement sampling in order to assess its impact and associated uncertainty on regional budgets. A possible approach to investigate this problem is to use a set of synthetic *in-situ*-like profiles extracted from an ensemble of forced ocean simulations started from different states and integrated with the same atmospheric forcing. Comparisons between the original ensemble outputs and the remapped, subsampled, *in-situ*-like profiles elucidate the contribution of chaotic ocean variability to OHC and regional sea level trends. Our results show that intrinsic variability may be large in eddy-active regions in the gridded model outputs, and remains substantial when using the *in-situ* sampling-based estimates. Using the latter, the same result is also found on large scales, for which atmospheric forcing has been identified as the main driver. Our results suggest accounting for this intrinsic ocean variability when assessing regional OHC and sea level trend budgets on decadal time scales.

1. Introduction

Planet Earth is warming at a global scale as more
radiation is entering than leaving at the top of the
atmosphere (Von Schuckman et al 2020). Directly
quantifying the Earth Energy Imbalance (EEI) is chal-
lenging and is best estimated by assessing changes
of the different climate reservoirs. The global ocean
has stored 90% of the excess heat associated with the

EEI, with a global heating rate of 0.87 ± 0.12 W m⁻<fnref n="2" />
over 2010–2018 (Von Schuckman *et al* 2020). Estim-
ating changes in the ocean heat content (hereafter
OHC) is therefore of great interest in accurately
assessing the EEI (Meyssignac *et al* 2019). Global
mean sea level (GMSL) rise is a direct response of
EEI increases and is mainly due to global ocean
warming (known as thermosteric sea level -TSL-)
and continental ice melting (from ice sheets and

<footnote n="1">¹ University of Brest/IFREMER/IRD/CNRS, LOPS, Brest, France</footnote>

<footnote n="2">² Université Grenoble Alpes, CNRS, IRD, Grenoble-INP, IGE, Grenoble, France</footnote>

<footnote n="3">³ CNRS/CERFACS, CECI, Toulouse, France</footnote>

<footnote>* Author to whom any correspondence should be addressed.</footnote>
mountain glaciers melting, WCRP Global Sea Level Budget Group 2018). Global mean TSL rise explains one-third of the net GMSL rise recorded by satellite altimetry data over 2005–2015 (Llovel et al 2019, Barnoud et al 2021). At a regional scale, OHC and TSL both show large deviations from their global trends since 2005 (Llovel and Lee 2015, Llovel and Terray 2016). Therefore, quantifying both OHC and TSL at regional scales is fundamental to better describing the impacts of the EEI.

OHC and TSL changes have been estimated based on *in-situ* observations for a number of years by various research groups (Boyer *et al* 2016). These estimates are based on *in-situ* temperature profiles with uneven spatial coverage. All estimates show that the upper 700 m of the global ocean has significantly warmed over past decades, but there are some discrepancies concerning the amplitudes of the interannual variability and decadal rates of warming amongst these various products (Lyman *et al* 2010). Differences arise from various sources of bias and error when estimating both the OHC and TSL changes. These include the quality control of data, the considered mapping methods, the baseline mean climatology to use (Lyman and Johnson 2014), the treatment of unsampled or undersampled regions (Durack *et al* 2014, Allison *et al* 2019) and the method for instrumental bias corrections (Cheng and Zhu 2015). The two main sources of discrepancies can be attributed to changing instrument technology and the uneven distribution of *in situ* records in time, and horizontal and vertical space (Boyer *et al* 2016).

Argo floats have provided near global coverage since 2005 and record temperature and salinity from the surface to 2000 m depth with a nominal temporal resolution of 10 days and a spatial resolution of about 3° × 3° (Roemmich et al 2019). This unprecedented data coverage has significantly reduced the uncertainties in OHC and TSL estimates from *in-situ*-based gridded products. However, there is still no perfect agreement when assessing both regional TSL and OHC trends, even over the well-sampled *in situ* period (Wang et al 2018).

Recently, a new possible driver of regional OHC trends has been highlighted using eddy-permitting ensemble ocean model simulations over 1980–2010 (Sérazin et al 2017). This multi-decadal variability spontaneously emerges from the numerical turbulent ocean without any low-frequency variability or trend in the atmospheric forcing (Penduff et al 2011). These chaotic fluctuations, resulting from non-linear ocean processes, may leave substantial random imprints on decadal regional OHC trends, in particular in the highly energetic western boundary currents and Antarctic Circumpolar Current (ACC). This intrinsic variability also leaves a large imprint on regional sea level trends (Llovel et al 2018) and interannual variability (Carret et al 2021) over 1993–2015.

Hydrographic measurements capture a combination of the ocean's response to atmospheric forcing and intrinsic ocean variability. Although the regional low-frequency variability may be reasonably captured by the design of the global Argo network, intrinsic variability may complicate the detection of the atmospherically-forced response and the attribution of observed signals to atmospheric drivers. Although the atmospherically-forced responses have been investigated for regional SL changes (Llovel et al 2010, Merrifield 2011, Forget and Ponte 2015) and regional OHC change (Llovel and Terray 2016, Thompson et al 2016), less attention has been paid to the intrinsic ocean variability contribution.

The present study aims to answer the following questions:

* Is the *in-situ* data coverage good enough to capture the spatial patterns simulated by the numerical model when assessing decadal TSL and OHC trends at the regional scale?

* How much uncertainty does the intrinsic ocean variability introduce in TSL and OHC decadal trends at the regional scale?

* In which regions does the atmospherically-forced response dominate the intrinsic ocean variability, and can the same results be obtained using unevenly-sampled *in-situ* data as with a numerical model?

Answering these questions is necessary to better assess the ability of hydrographic networks to estimate regional changes in atmospherically-forced OHC and TSL.

We approach this problem by using synthetic *in-situ* profiles (Allison et al 2019) derived from a large ensemble of numerical forced ocean simulations (Penduff et al 2014, Bessières et al 2017). We choose to focus on the well-sampled period since 2005, where the Argo data coverage is almost global, in order to assess the ability of the current *in-situ* observing system to capture ocean intrinsic variability. The 50 sets of synthetic profiles were mapped using the In-Situ Analysis System (ISAS), using its standard configuration for mapping real *in-situ* profiles (Gaillard et al 2016). This approach allows us to investigate the hydrographic coverage (down to 2000 m depth) needed to estimate the TSL and OHC changes and the imprint of intrinsic ocean variability as simulated by the ensemble numerical model and as seen by the *in-situ* network.

The paper is organized as follows. Section 2 describes the data and methods used in this study. Section 3 presents the imprints of atmospherically-forced and intrinsic variability on regional OHC and TSL trends. Finally, the results are summarized and discussed before addressing the implications of our findings in section 4.
## 2. Methods

### 2.1. OCCIPUT ensemble simulations

We exploit the Oceanic Chaos—ImPacts, struc- ture, predictability (OCCIPUT) ensemble of 1/4° ocean/sea-ice simulations (Penduff et al 2014, Bessières et al 2017). This ensemble consists of a set of 50 global hindcasts at 1/4° horizontal resolution integ- rated over 1960–2015. The numerical core is based on the NEMO3.5 model. Each member was initialized on January 1st 1960 from the final state of a 21 year single member spin-up. Then, a slight stochastic per- turbation in the seawater equation (Brankart 2013) was applied within each ensemble member during one year (i.e. 1960). This perturbation was switched off at the end of 1960 producing 50 different oceanic states. Each ensemble member was integrated from January 1st 1961 until the end of 2015 with the same atmospheric forcing based on Drakkar Forcing Set 5.2 (DFS5.2; Dussin et al 2016). These 50 ensemble members share the same numerical code and atmo- spheric forcing, but have different initial conditions. We will refer to this set of simulations as the OCCI- PUT ensemble simulations.

We also use a complementary one-member sim- ulation sharing the same code and setup but forced with a climatological forcing. This climatological sim- ulation has been integrated over 327 years and has been forced each year by the same annual atmospheric cycle derived from DFS5.2 (Penduff et al 2011). We consider the years 46–56 of this run, corresponding to the 2005–2015 period in the OCCIPUT simulations, to estimate any spurious model drift in the OCCI- PUT ensemble. This simulation will be referred to as OCCICLIM.

### 2.2. ISAS ensemble analyses

First, synthetic temperature and salinity profiles, cor- responding to the observed profiles contained in the EN4 database (excluding eXpendable BathyThermo- graphs profiles as the latter do not reach 2000 m depth and do not record salinity; Good et al 2013) have been extracted online from each member dur- ing the OCCIPUT ensemble integration. It is worth noting that the synthetic profiles are representative of the ocean state averaged over 1/4° model grid boxes rather than real point-wise longitude/latitude profiles as in real-world observations. We obtain 50 sets of in- situ-like synthetic profiles of temperature and salin- ity from the surface to 2000 m depth over 2005–2015. A corresponding set of in-situ-like profiles has been extracted from the OCCICLIM run to correct for any possible model drift in the synthetic OCCIPUT pro- files.

Each set of synthetic profiles is remapped on a regular 0.5° × 0.5° grid over 152 standard depths using the ISAS tool (Gaillard et al 2016). ISAS is based on optimal interpolation (OI; Bretherton et al 1976), and is designed to map in-situ profiles onto a regular

grid on a monthly basis to fill the spatio-temporal gaps between missing data and provide an associ- ated error estimate. We will refer to this remapped set of 50-member simulations as the ISAS OCCIPUT ensemble. We also remapped the in-situ-like profiles from OCCICLIM, and will refer to this as the ISAS OCCICLIM solution.

The ISAS configuration was the same as used for ISAS15 configuration (Kolodziejczyk et al 2019, Kolodziejczyk et al 2021). For mapping each member of ISAS OCCIPUT and ISAS OCCICLIM, we follow four steps described in Gaillard et al (2016):

* The synthetic temperature and salinity profiles are standardized onto 152 levels.

* A monthly mean and temporal standard devi- ation climatology of temperature and salinity are generated over the considered time period (i.e. 2005–2015). This climatology consists of a 5° × 5° binning (mean and variance) of temperature and salinity standardized synthetic profiles for each cal- endar month. The monthly mean climatology over the period 2005–2015 is used as a first guess for interpolating the monthly ISAS OCCIPUT temper- ature and salinity gridded fields.

* The annual variance climatology is used as the a priori variance for the OI.

* The mapping procedure is applied for each month between January 2005 to December 2015 using dif- ferent Gaussian weights on an a priori variance cor- responding to large (300 km) and meso-covariance scales (proportional to the Rossby radius), and unresolved scales.

This procedure is repeated for the in-situ-like synthetic profiles extracted from each member of the OCCIPUT ensemble simulations, and from the OCCICLIM simulation. We therefore obtain an ensemble of 50 gridded temperature and salinity fields (ISAS OCCIPUT) and one member gridded temperature and salinity fields (ISAS OCCICLIM).

### 2.3. OHC and TSL computations

OHC is computed by integrating the temperat- ure fields from the surface to 2000 m depth over 2005–2015 following the equation:

$$ \text{OHC}(x, y, t) = \rho_0 C_p \int_{2000}^{0} T(x, y, z, t) \, dz $$

where $\rho_0$ is the seawater density, taken here to be 1026 kg m$^{-3}$, and $C_p$ is the heat capacity of sea water, taken here to be 4000 J kg$^{-1}$ K$^{-1}$. $T(x, y, z, t)$ cor- responds to the potential temperature anomaly with respect to 0 °C.

TSL is computed by estimating the density of seawater at each standard depth with the temperat- ure fields and the salinity climatology computed as
$$TSL(x, y, t) = -\frac{1}{\rho_0} \int_{2000}^{0} (\rho(T(x, y, z, t), \bar{S}(x, y, z), p) - \rho(0, 35, p)) dz$$

where p is the pressure at the standard depth. $\bar{S}$ is the monthly mean salinity climatology estimated over 2005–2015. We compute OHC and TSL fields for both OCCIPUT and ISAS OCCIPUT ensembles over 2005–2015.

We also computed the OHC and TSL fields from OCCICLIM and ISAS OCCICLIM simulations to assess any spurious model drift. At each grid point of these two data sets, we estimated the decadal OHC and TSL linear trends. We remove these decadal trends from each member of the OCCIPUT ensemble simulation and from the ISAS OCCIPUT ensemble gridded products for each field. We thus obtain two 50-member ensembles for OHC and TSL, both with gridded outputs, corrected for any spurious model drift.

## 2.4. Ensemble statistics

The processing steps presented above yield two ensembles of 50-member TSL and OHC simulations. For each grid point, we compute the decadal linear trends from each member of these ensembles. We then consider the decadal ensemble mean trend to be an estimate of the atmospherically-forced response, and the ensemble standard deviation $\sigma_{Vi}$ to be an estimate of the uncertainty associated with the intrinsic ocean variability for V. V will be either the regional trends of OHC(x,y,N) or TSL(x,y,N) from OCCIPUT or ISAS OCCIPUT, where (x,y) denotes the spatial coordinates and N represents the number of members (50 simulations).

The uncertainty is defined using the ensemble standard deviation as follows:

$$\sigma_{Vi} = \sqrt{\frac{1}{N-1} \sum_{N=1}^{50} (V_n - <V>)^2}$$

where $<V>$ represents the ensemble mean of V and $V_n$ the nth member of the considered ensemble (Llovel et al 2018).

We will also investigate the signal-to-noise ratio (hereafter SNR) defined as $|<V>|/\sigma_{Vi}$ to compare the relative magnitudes of the atmospherically-forced and intrinsic ocean variability, and to identify regions where the atmospherically-forced response dominates over the intrinsic variability and vice versa.

# 3. Results

## 3.1. The impact of spatial distribution of in-situ observations on regional trends over 2005–2015

We first investigate the ability of the hydrographic network sampling to accurately capture regional TSL trends over 2005–2015. We consider one member (member #1) from the OCCIPUT ensemble simulation and its corresponding member from ISAS OCCIPUT ensemble simulation. This member corresponds to one possible realization of the simulated ocean. Overall, we find similar large-scale spatial patterns for regional TSL trends of member #1 from OCCIPUT (figure 1(a)) and ISAS OCCIPUT (figure 1(b)). Positive trends are found in the Indian Ocean, the Central East Pacific and the mid-latitude North Atlantic Oceans. Negative values are also found in the Western Pacific and the North East Atlantic Oceans. These patterns are similar to the regional sea level trend maps simulated by OCCIPUT over the same period (Llovel et al 2018). This highlights that temperature change is the main driver of regional sea level trends over this decade. Overall, the amplitudes of regional TSL trends are smaller in member #1 ISAS OCCIPUT than in member #1 OCCIPUT. Furthermore, western boundary currents and the ACC present smoother patterns with less prominent mesoscale features.

Figure 1(c) shows the differences between figures 1(a) and (b). Large differences are found in the western boundary currents. This indicates that the ISAS OI method tends to slightly underestimate the trend amplitudes, and smooths out regional patterns associated with mesoscale variability and front displacements. This is due to a combination of two factors: first, trends may be underestimated by the OI method due to the relaxation toward the climatology. Second, as most of the mesoscale features cannot be resolved by the sampling of the in-situ network, they are smoothed out by the OI method. In the subtropical gyres, zonal bands of opposite trends are simulated by the OCCIPUT ensemble simulations. Similar features have been already identified in ocean simulations (Penduff et al 2011) and observed by both Argo floats and satellite altimetry and attributed to slow-varying zonal jet structures (van Sebille et al 2011). These patterns are not resolved by the OI analyses. We thus find that the in-situ gridded product has the ability to capture the large-scale patterns of the regional TSL trends that are simulated by the NEMO model, but not mesoscale features.
<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 1.** TSL trend maps computed over January 2005–December 2015 for member #1 of (a) OCCIPUT (b) ISAS OCCIPUT and (c) the difference between (a) and (b). Values are in mm yr⁻¹.</figcaption></figure>

Similar results are found for OHC trend maps as resolved by member #1 for both OCCIPUT and ISAS OCCIPUT (see figure S1 in supplementary material available online at stacks.iop.org/ERL/17/044063/mmedia). Interestingly, these regional OHC trend patterns are similar to the in-situ-based OHC trend map computed over 2006–2015 described in

Kolodziejczyk et al. (2019). Indeed, positive OHC trends are found in the subtropical North Atlantic, North Indian, eastern tropical and South Pacific Oceans. Negative trends are found in the subpolar North Atlantic and western tropical Pacific Oceans as simulated by member #1 for both OCCIPUT and ISAS OCCIPUT. This comparison confirms
c)

<figure><img src="image_7.png" /><figcaption>**Figure 1.** (Continued.)</figcaption></figure>

that the model has skills in reproducing decadal OHC trends. In addition, this suggests that the *in-situ* spatial coverage is dense enough for assessing regional decadal TSL and OHC trends over 2005–2015.

## 3.2. Atmospherically-forced regional trends over 2005–2015

The atmospherically-forced (i.e. ensemble mean) regional TSL trends from OCCIPUT simulations are shown in figure 2(a). Large scale regional patterns are broadly similar to those of member #1 from OCCIPUT ensemble (figure 1(a)), especially for the tropical Pacific Ocean where positive trends are found in the east and negative trends in the west. Positive trends are also identified for the North Indian Ocean and the subtropical North Atlantic Ocean. The latter basin exhibits negative trends in the subpolar gyre, especially in the eastern part.

However, the atmospherically-forced regional TSL trends are not identical to those obtained for member #1 OCCIPUT, as shown by the difference trend map (figure S3(a) in supplementary material). The differences are locally large, especially in the western boundary currents and the ACC, with absolute values exceeding 10 mm yr⁻¹. The subtropical gyres also exhibit large differences between the ensemble mean and member #1, with values exceeding 5 mm yr⁻¹ and zonal extensions corresponding to oceanic jets. We find similar results for regional OHC trend maps based on OCCIPUT (see figures

S2(a) for the atmospherically-forced trend map and S1(a) for member #1 in supplementary material).

Figure 2(b) shows the atmospherically-forced regional TSL trends estimated from ISAS OCCIPUT. As previously mentioned, this trend map looks like member #1 presented in figure 1(b). Again, large-scale decadal trends of regional TSL are captured by the *in-situ* data distribution, albeit with lower amplitude. Differences are nevertheless found in both the amplitude and spatial patterns of the highly energetic western boundary currents (Gulf Stream, Kuroshio, Zapiola and Agulhas currents; figure S3(b) in supplementary material). The ISAS OI tends to smooth large scale spatial patterns and trend amplitudes. We find similar results for regional OHC trend maps based on ISAS OCCIPUT (figures S2(b) for the atmospherically-forced trend map and S1(b) for member #1 in supplementary material).

## 3.3. Imprint of intrinsic variability on regional trends

Figure 2(c) depicts the imprint of ocean intrinsic variability on regional TSL trends (ensemble spread from OCCIPUT simulations). Spread values larger than 1 mm yr⁻¹ are found for the subtropical gyres (both south and north) of the Pacific and the Atlantic Oceans. The Indian Ocean also experiences trend values larger than 1 mm yr⁻¹ for the southern subtropical gyre and near the coasts of India and the Somalia-Yemen-Oman coasts. Spread values ranging from 5 mm yr⁻¹ to more than 12 mm yr⁻¹ are
a)

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

b)

<figure><img src="image_9.png" /><figcaption>**Figure 2.** Atmospherically forced TSL trend maps (ensemble means from the 50 members) over January 2005–December 2015 for (a) OCCIPUT and (b) ISAS OCCIPUT. Imprint of chaotic ocean variability on TSL trend maps (ensemble spreads from the 50 members) for (c) OCCIPUT and (d) ISAS OCCIPUT. Values are in mm yr⁻¹.</figcaption></figure>

found in the western boundary currents and their vicinity.

Figure 2(d) shows the imprint of ocean intrinsic variability on regional TSL trends (ensemble spread from the ISAS OCCIPUT ensemble runs) over 2005–2015. The ensemble spread of regional TSL

trends shows large values of intrinsic variability in highly energetic currents and the ACC. The ensemble spread spatial patterns are clearly reduced (compared to the OCCIPUT ensemble) almost everywhere in terms of both extent and amplitude. The ACC region shows spread values larger than 5 mm yr⁻¹ in the
<figure><img src="image_10.png" /><figcaption>d)</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Figure 2.** (Continued.)</figcaption></figure>

Agulhas region only. The large spread values identified in the OCCIPUT ensemble runs ranging from 5 mm yr⁻¹ to more than 12 mm yr⁻¹ south of Australia and between east Australia and New Zealand are now smaller than 4 mm yr⁻¹. The South Atlantic Ocean also exhibits weaker spread in the trend when compared to the OCCIPUT ensemble

runs (figure 2(c)). Large values of spread in the TSL trend are still present in the Kuroshio and Gulf stream regions. However, the spatial extents of these regions are significantly reduced, especially in the Gulf stream region. Both *in-situ* sampling and the ISAS OI method thus tend to filter out scales smaller than about 300 km. However, chaotic variability is
<figure><img src="image_12.png" /></figure>

<figure><img src="image_13.png" /><figcaption>**Figure 3.** SNR maps. The signal is the decadal atmospherically-forced TSL trends (in absolute values) and the noise is considered to be the decadal intrinsic ocean variability imprint on decadal TSL trends (ensemble trend spread) for (a) OCCIPUT ensemble simulations and for (b) ISAS OCCIPUT ensemble simulations. Values are unit less.</figcaption></figure>

still present, especially in the subtropical gyre and energetic western boundary currents, and in the ACC. Both spread maps show the limited influence of intrinsic variability in the tropics. This is consistent with the major contribution of wind forcing to inter-annual and decadal regional sea level variability in these regions (Merrifield 2011, England et al 2014, Llovel et al 2018).

Figures S2(c) and (d) show the ocean intrinsic variability imprint on regional OHC trends over

2005–2015 for the OCCIPUT and ISAS OCCIPUT ensemble simulations. Similar results are found compared to the TSL, and intrinsic ocean variability is noted to have an imprint on regional OHC trend over 2005–2015 that is substantial for highly energetic western boundary currents.

### 3.4. Signal-to-noise ratio maps

To evaluate the origin (intrinsic ocean variability versus atmospherically-forced variability) of regional
TSL trends, we consider the SNR maps for OCCIPUT (figure 3(a)) and ISAS OCCIPUT (figure 3(b)). We recall the reader that the SNR is defined as the ensemble mean trend (absolute value) divided by the ensemble spread trend at each grid point. SNR values greater than 2 (at 95% confidence interval) denote trends that can be attributed to the external forcing (natural variability and anthropogenic drivers). Figure 3(a) shows that SNR > 2 over a large part of the intertropical band (15°S–15°N) and of continental shelves (i.e. Argentina, European and Arctic regions). On the other hand, SNR < 2 means that TSL trends cannot be unambiguously attributed to the atmospheric forcing: this is the case over 44% of the global ocean area (i.e. SNR < 2).

Very interestingly, figure 3(b) shows similar spatial patterns, but with coarser spatial resolution. The degraded spatial resolution comes from the OI mapping. SNR > 2 values are found in the tropical oceans, the center of subtropical gyres of the South Pacific Ocean and the North Atlantic Ocean, and the rim of subtropical gyres of the North Atlantic Ocean. Therefore, the *in-situ* spatial distribution remapped using the ISAS tool preserves -with a coarser spatial resolution- the SNR spatial patterns found in the OCCIPUT ensemble simulation. SNR patterns are however not well reproduced in the Arctic and Southern Oceans, mainly due to the limited number of *in-situ* measurements in these regions. We thus find that, outside the high latitudes, the *in-situ* data coverage is dense enough to preserve the atmospherically-forced and intrinsic ocean variability ratios after binning the individual profiles with the ISAS tool. We also find that TSL trends cannot be unambiguously attributed to the atmospheric forcing over 40% of the global ocean area. This value remains large even after the remapping method.

Figures S4(a) and (b) show the SNR maps for the regional OHC trends as resolved by both OCCIPUT and ISAS OCCIPUT ensemble simulations, respectively. We find overall similar spatial patterns for OHC SNR maps compared to TSL SNR maps for both OCCIPUT and ISAS OCCIPUT ensemble simulations. We find values larger than 2 in the equatorial band, the south Pacific subtropical gyre, the north Atlantic subtropical gyre for both ensemble simulations. As for regional TSL trends, the ISAS tool preserves the SNR between the atmospherically-forced versus the intrinsic ocean variability for the regional OHC trends. In other words, OHC trends cannot be unambiguously attributed to the atmospheric forcing over 45% and 39% over the global ocean area for the OCCIPUT and ISAS OCCIPUT ensemble simulations, respectively.

## 4. Conclusions and discussion

*In-situ* data and, more recently, Argo floats have remodeled our knowledge of ocean circulation and

dynamics, and revealed large regional variability in decadal trends of TSL and OHC. Using synthetic profiles extracted from climate model output, Allison et al (2019) showed that the introduction of Argo profiles reduces the mapping method uncertainty in estimating the variability and spatial structures of OHC anomalies. We adopt the same strategy by assessing TSL and OHC trends using synthetic profiles extracted from an ensemble of forced ocean simulations sharing the same atmospheric forcing but with different initial conditions. By remapping an ensemble of *in-situ*-like synthetic profiles onto a regular grid using the ISAS tool, we show that the current *in-situ* coverage is dense enough to capture the regional decadal TSL and OHC trends simulated by the NEMO ocean model over 2005–2015. This is in line with recent studies showing that the international Argo program allows us to mitigate spatial and temporal sampling bias and errors (Durack et al 2014, Allison et al 2019). We note, however, coarser spatial resolution and smaller amplitudes, due to both mapping procedure (using Gaussian weights) and climatology relaxation in undersampled regions. Decadal trends in high latitude regions are not well captured as the *in-situ* coverage is still sparse. This is a challenge to be overcome in the near future as these regions are experiencing faster warming than the mid latitudes or the tropics (Meredith et al 2019). Comparing a single ocean simulation and its associated synthetic profiles is an effective approach to give confidence in the OI configuration.

Our approach allows us not only to assess the reliability of the mapping method, but also to quantify for the first time the intrinsic ocean variability recorded by the observing system. The atmospherically-forced response has been investigated using the OCCIPUT and ISAS OCCIPUT ensemble datasets. Atmospherically-forced trends are rather similar in the OCCIPUT and ISAS OCCIPUT ensembles; however, differences exist, especially in western boundary currents and basin interiors. This general resemblance supports the hypothesis that large scale TSL and OHC regional trends are a direct response to the atmospheric forcing.

We also find that the OCCIPUT ensemble indicates that intrinsic ocean variability has a significant influence on both TSL and OHC regional trends over 2005–2015, especially in the western boundary currents and the ACC. These patterns are consistent with the imprints of intrinsic ocean variability on sea level trends over the same period (Llovel et al 2018). Again, this confirms that temperature changes are the main driver of regional sea level decadal trends, as was previously shown with ocean models (Fukumori and Wang 2013, Griffies et al 2014) and with *in-situ*-based gridded products (Llovel and Lee 2015).

More interestingly and suggested here for the first time, substantial intrinsic ocean variability remains significant after remapping *in-situ*-like profiles with
the ISAS tool -with contributions to the TSL trend exceeding 10 mm yr⁻¹ locally. Although the mapping of sparse *in-situ* profiles filters out part of the intrinsic variability, hot spots still remain, especially in the Kuroshio current and its extension, the Gulf Stream and the Agulhas current. This finding highlights the need to consider this uncertainty when assessing regional budgets of quantities such as sea level or OHC when analyzing *in-situ* gridded products. However, this assessment of intrinsic variability after remapping onto a regular grid should be further analyzed using other remapping tools, as there is no consensus on the best method to interpolate unevenly-spaced observations and the mapping methods are a key source of uncertainty (Allison *et al* 2019).

We also show the presence of regions where the atmospherically-forced response dominates the intrinsic variability and vice versa for regional TSL and OHC trends over 2005–2015. Unsurprisingly, the tropics show a strong direct response to the atmosphere, along with the continental shelves. We find that TSL trends cannot be unambiguously attributed to the atmospheric forcing over 44% and 40% of the global ocean area in the OCCIPUT and ISAS OCCIPUT datasets, respectively. Concerning the OHC trends, these fractions reach 45% and 39% from OCCIPUT and ISAS OCCIPUT datasets. This means that the intrinsic variability introduces uncertainty in attributing these decadal trends to an atmospheric influence. These fractions are large and consistent with decadal sea level trends over 2005–2015 (see Llovel *et al* 2018).

We find that the mapped *in-situ* products show skill in reproducing the respective contributions of atmospherically-forced and intrinsic ocean variability as simulated by the OCCIPUT ensemble simulations when considering regional TSL and OHC trends. This suggests that no or little spurious variability has been introduced when regridding the *in-situ* profiles onto a regular grid for assessing the decadal trends of regional TSL and OHC. Indeed, mapping methods can introduce spurious variability. For instance, mesoscale features that cannot be resolved by the Argo network sampling may introduce distortion or spurious patterns in time and/or space in gridded products. Therefore, the covariance scales for the OI methods are critical to reproduce accurately large-scale patterns and variability in order to filter out unresolved oceanic processes (Allison *et al* 2019). Undeniably, the current configuration of ISAS mapping method does not introduce any such spurious variability.

Our results should be interpreted with care, as intrinsic variability is model dependent. To verify this, similar analyses should be conducted with other ensemble simulations, using different models, atmospheric forcings and spatial resolutions. However, we believe that it is unlikely that intrinsic variability has

been overestimated in this work, as it has been shown to increase with increasing model spatial resolution for interannual sea level variability (Sérazin *et al* 2015). Our results highlight the substantial imprint of intrinsic ocean variability on regional TSL and OHC trends, both in the NEMO numerical model and in *in-situ*-based remapped products over 2005–2015.

## Data availability statement

The data that support the findings of this study are available upon reasonable request at therry.penduff@cnrs.fr.

## Acknowledgments

This work is a contribution of the CRATERE project supported by the French national program LEFE/INSU. This work is a contribution to the OCCIPUT and PIRATE projects. PIRATE (https://sealevel.jpl.nasa.gov/science/ostscienceteam/scientistlinks/scientificinvestigations2017/penduff/) is funded by CNES through the Ocean Surface Topography Science Team (OSTST). OCCIPUT (http://meom-group.github.io/projects/occiput) has been funded by ANR through contract ANR-13-BS06-0007-01. We acknowledge that the results of this research have been achieved using the PRACE Research Infrastructure resource CURIE based in France at TGCC.

## Conflict of interest

The authors declare no conflict of interest.

## ORCID iDs

William Llovel https://orcid.org/0000-0002-0798-7595

Nicolas Kolodziejczyk https://orcid.org/0000-0002-0751-1351

Sally Close https://orcid.org/0000-0001-9665-5316

Thierry Penduff https://orcid.org/0000-0002-0407-8564

Jean-Marc Molines https://orcid.org/0000-0003-1665-6816

Laurent Terray https://orcid.org/0000-0001-5512-7074

## References

* Allison L C *et al* 2019 *Environ. Res. Lett.* 14 084037

* Barnoud A *et al* 2021 *Contributions of altimetry and Argo to non-closure of the global mean sea level budget since 2016* *Geophys. Res. Lett.* 48 e2021GL092824

* Bessières L, Leroux S, Brankart J-M, Molines J-M, Moine M-P, Bouttier P-A, Penduff T, Terray L, Barnier B and Sérazin G 2017 *Development of a probabilistic ocean modelling system based on NEMO 3.5: application at eddying resolution* *Geosci. Model Dev.* 10 1091–106
Boyer T *et al* 2016 Sensitivity of global upper ocean heat content estimates to mapping methods, XBT bias corrections, and baseline climatologies *J. Clim.* **29** 4817–42

Brankart J-M 2013 Impact of uncertainties in the horizontal density gradient upon low resolution global ocean modelling *Ocean Model.* **66** 64–76

Bretherton F, Davis R and Fandry C 1976 Technique for objective analysis and design of oceanographic experiment applied to MODE-73 *Deep Sea Res. Oceanogr. Abstr.* **23** 559–82

Carret A, Llovel W, Penduff T and Molines J-M 2021 Atmospherically forced and chaotic interannual variability of regional sea level and its components over 1993–2015 *J. Geophys. Res.* **126** e2020JC017123

Cheng L and Zhu J 2015 Influences of the choice of climatology on ocean heat content estimation *J. Atmos. Ocean Technol.* **32** 388–94

Durack P *et al* 2014 Quantifying underestimates of long-term upper-ocean warming *Nat. Clim. Change* **4** 999–1005

Dussin R, Barnier B, Brodeau L and Molines J-M 2016 The making of Drakkar forcing set DFS5 *DRAKKAR/MyOcean Report* 01 04 2016 (Grenoble: LGGE)

England M H *et al* 2014 Recent intensification of wind-driven circulation in the Pacific and the ongoing warming hiatus *Nat. Clim. Change* **4** 222–7

Forget G and Ponte R M 2015 The partition of regional sea level variability *Prog. Oceanogr.* **137** 173–95

Fukumori I and Wang O 2013 Origins of heat and freshwater anomalies underlying regional decadal sea level trends *Geophys. Res. Lett.* **40** 563–7

Gaillard F, Reynaud T, Thierry V, Kolodziejczyk N and von Schuckmann K 2016 *In situ*-based reanalysis of the global ocean temperature and salinity with ISAS: variability of the heat content and steric height *J. Clim.* **29** 1305–23

Good S A, Martin M J and Rayner N A 2013 EN4: quality controlled ocean temperature and salinity profiles and monthly objective analyses with uncertainty estimates *J. Geophys. Res. Oceans* **118** 6704–16

Griffies S M *et al* 2014 An assessment of global and regional sea level for years 1993–2007 in a suite of interannual CORE-II simulations *Ocean Modell.* **78** 35–89

Kolodziejczyk N, Llovel W and Portela E 2019 Interannual variability of upper ocean water masses as inferred from Argo array *J. Geophys. Res.* **124** 6067–85

Llovel W, Fukumori I and Meyssignac B 2013 Depth-dependent temperature change contributions to global mean thermosteric sea level rise from 1960 to 2010 *Glob. Planet. Change* **101** 113–8

Llovel W, Guinehut S and Cazenave A 2010 Regional and interannual variability in sea level over 2002–2009 based on satellite altimetry, Argo float data and GRACE ocean mass *Ocean Dyn.* **60** 1193–204

Llovel W and Lee T 2015 Importance and origin of halosteric contribution to sea level change in the southeast Indian Ocean during 2005–2013 *Geophys. Res. Lett.* **42** 1148–57

Llovel W, Penduff T, Meyssignac B, Molines J-M, Terray L, Bessières L and Barnier B 2018 Contributions of atmospheric forcing and chaotic ocean variability to regional sea level trends over 1993–2015 *Geophys. Res. Lett.* **45** 405–13

Llovel W, Purkey S, Meyssignac B, Blazquez A, Kolodziejczyk N and Bamber J 2019 Global ocean freshening, ocean mass increase and global mean sea level rise over 2005–2015 *Sci. Rep.* **9** 17717

Llovel W and Terray L 2016 Observed southern upper-ocean warming over 2005–2014 and associated mechanisms *Environ. Res. Lett.* **11** 1240

Lyman J M and Johnson G C 2014 Estimating global ocean heat content changes in the upper 1800 m since 1950 and the influence of climatology choice *J. Clim.* **27** 5

Lyman J *et al* 2010 Robust warming of the global upper ocean *Nature* **465** 334–7

Meredith M *et al* 2019 Polar regions *IPCC Special Report on the Ocean and Cryosphere in a Changing Climate* ed H-O Pörtner, D C Roberts, V Masson-Delmotte, P Zhai, M Tignor, E Poloczanska, K Mintenbeck, A Alegría, M Nicolai, A Okem, J Petzold, B Rama and N M Weyer (Cambridge: Cambridge University Press) pp 203–20

Merrifield M A 2011 A shift in western tropical Pacific sea level trends during the 1990s *J. Clim.* **24** 4126–38

Meyssignac B *et al* 2019 Measuring global ocean heat content to estimate the earth energy imbalance *Front. Mar. Sci.* **6** 432

Nicolas K, Annaig P-M and Fabienne G 2021 ISAS temperature and salinity gridded fields *SEANOE* (https://doi.org/10.17882/52367)

Penduff T, Barnier B, Terray L, Bessières L, Sérazin G, Gregorio S, Brankart J, Moine M, Molines J and Brasseur P 2014 Ensembles of eddying ocean simulations for climate, CLIVAR exchanges, special issue on high resolution ocean climate modelling p 19

Penduff T, Juza M, Barnier B, Zika J, Dewar W K, Treguier A-M, Molines J M and Audiffren N 2011 Sea-level expression of intrinsic and forced ocean variabilities at interannual time scales *J. Clim.* **24** 5652–70

Roemmich D *et al* 2019 On the future of Argo: a global, full-depth, multi-disciplinary array *Front. Mar. Sci.* **6** 439

Sérazin G, Jaymond A, Leroux S, Penduff T, Bessières L, Llovel W, Barnier B, Molines J-M and Terray L 2017 A global probabilistic study of the ocean heat content low-frequency variability: atmospheric forcing versus oceanic chaos *Geophys. Res. Lett.* **44** 5580–9

Sérazin G, Penduff T, Grégorio S, Barnier B, Molines J-M and Terray L 2015 Intrinsic variability of sea-level from global 1/12° ocean simulations: spatio-temporal scales *J. Clim.* **28** 4279–92

Thompson P R, Piecuch C G, Merrifield M A, McCreary J P and Firing E 2016 Forcing of recent decadal variability in the equatorial and North Indian Ocean *J. Geophys. Res.* **121** 6762–78

van Sebille E, Kamenkovich I and Willis J K 2011 Quasi-zonal jets in 3D Argo data of the northeast Atlantic *Geophys. Res. Lett.* **38** L02606

von Schuckman K *et al* 2020 Heat stored in the Earth system: where does the energy go? *Earth Syst. Sci. Data* **12** 2013–41

Wang G *et al* 2018 Consenses and discrepancies of basin-scale ocean heat content changes in different ocean analyses *Clim. Dyn.* **50** 2471–87

WCRP Global Sea Level Budget Group 2018 Global sea-level budget 1993–present *Earth Syst. Sci. Data* **10** 1551–90