# Chaotic oceanic excitation of low-frequency polar motion variability

Lara Börger, Michael Schindelegger, Mengnan Zhao, Rui Ponte, Anno Löcher, Bernd Uebbing, Jean-Marc Molines, Thierry Penduff

► To cite this version:

Lara Börger, Michael Schindelegger, Mengnan Zhao, Rui Ponte, Anno Löcher, et al. Chaotic oceanic excitation of low-frequency polar motion variability. Earth System Dynamics, 2024, 16 (1), pp.75-90. 10.5194/esd-2024-21. hal-04765246v2

HAL Id: hal-04765246

https://hal.science/hal-04765246v2

Submitted on 27 Nov 2025

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_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

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

# Chaotic oceanic excitation of low-frequency polar motion variability

Lara Börger<fnref n="1" />, Michael Schindelegger<fnref n="1" />, Mengnan Zhao<fnref n="2" />, Rui M. Ponte<fnref n="2" />, Anno Löcher<fnref n="1" />, Bernd Uebbing<fnref n="1" />, Jean-Marc Molines<fnref n="3" />, and Thierry Penduff<fnref n="3" />

**Correspondence:** Lara Börger (lboerger@igg.uni-bonn.de)

Received: 8 July 2024 – Discussion started: 29 July 2024

Revised: 24 October 2024 – Accepted: 6 November 2024 – Published: 9 January 2025

**Abstract.** Studies of Earth rotation variations generally assume that changes in non-tidal oceanic angular momentum (OAM) manifest the ocean's direct response to atmospheric forces. However, fluctuations in OAM may also arise from chaotic intrinsic ocean processes that originate in local nonlinear (e.g., mesoscale) dynamics and can map into motions and mass variations at basin scales. To examine whether such random mass redistributions effectively excite polar motion, we compute monthly OAM anomalies from a 50-member ensemble of eddy-permitting global ocean/sea ice simulations that sample intrinsic variability through a perturbation approach on model initial conditions. The resulting OAM (i.e., excitation) functions, $\hat{\chi}^O$, are examined for their spread, spectral content, and role in the polar motion excitation budget from 1995 to 2015. We find that intrinsic $\hat{\chi}^O$ signals are comparable in magnitude to the forced component at all resolved periods except the seasonal band, amounting to ~46% of the total oceanic excitation (in terms of standard deviation) on interannual timescales. More than half of the variance in the intrinsic mass term contribution to $\hat{\chi}^O$ is associated with a single global mode of random bottom pressure variability, likely generated by nonlinear dynamics in the Drake Passage. Comparisons of observed interannual polar motion excitation against the sum of known surficial mass redistribution effects are sensitive to the representation of intrinsic $\hat{\chi}^O$ signals: reductions in the observed excitation variance can be as high as 68% or as low as 50% depending on the choice of the ensemble member. Chaotic oceanic excitation thus emerges as a new factor to consider when interpreting low-frequency polar motion changes in terms of core–mantle interactions or employing forward-modeled OAM estimates for Earth rotation predictions.

# 1 Introduction

Changes in ocean currents and the global distribution of ocean mass imply variability in oceanic angular momentum (OAM). Transfer of portions of this momentum to the solid Earth causes measurable fluctuations in the latter's rotation, comprising both changes in angular velocity and the movement of a conventional pole relative to Earth's crust – commonly referred to as wobble or polar motion (Gross, 1992). Past studies (e.g., Ponte et al., 1998; Gross et al., 2003; Bizouard and Seoane, 2010; Harker et al., 2021; Afroosa et al., 2021; Börger et al., 2023) have shown that variability in non-

tidal OAM is indeed a prominent source for wobble excitation on timescales from days to several years, modulating or even exceeding the rotational contributions of other geophysical fluids. An assumption in all of these works is that the path to accurate OAM estimates is to encode known physics in a numerical ocean model and to integrate that very model in time under the necessary external (mostly atmospheric) forcings, either with or without data constraints. However, such an approach overlooks the possibility of a significant random component in the evolving ocean state, emerging internally within the fluid irrespective of atmospheric fluctua-

<footnote n="1">¹Institute of Geodesy and Geoinformation, University of Bonn, Bonn, Germany</footnote>

<footnote n="2">²Atmospheric and Environmental Research (AER), Lexington, MA, USA</footnote>

<footnote n="3">³Université Grenoble Alpes, CNRS, INRAE, IRD, Grenoble INP, Institut des Géosciences de l'Environnement (IGE), Grenoble, France</footnote>
tions. Evidently, any such intrinsic chaotic variability would be inaccessible – or at least challenging to capture – with the usual means of forward or inverse modeling.

A well-known example of the ocean's intrinsic variability is the nonlinear, turbulent flow as manifested in mesoscale eddies. While these eddies operate at typical scales of O(10–100) km and O(10–100) d (Rhines, 2019), their kinetic energy is known to cascade toward much larger spatial and longer temporal scales (e.g., Venaille et al., 2011; Arbic et al., 2014; Sérazin et al., 2018). In detail, ocean simulations performed in the eddying regime (at grid spacings ≤ 1/4°) have demonstrated that interactions between relatively short-lived mesoscale eddies and large-scale baroclinic instability can generate chaotic interannual to secular variability in regional sea level (e.g., Llovel et al., 2018; Carret et al., 2021), large-scale volume transports (Leroux et al., 2018; Cravatte et al., 2021), or the eddy field itself (Hogg et al., 2022). Moreover, the same simulations suggest that nonlinear dynamics feed random fluctuations in ocean mass between entire basins (Zhao et al., 2023). Given this large-scale chaotic imprint on dynamical variables, intriguing questions arise as to whether oceanic intrinsic variability also affects OAM and what such random variability potentially means for the interpretation of polar motion observations.

We address these open questions using output from a 1/4° large-ensemble ocean hindcast experiment, as analyzed in numerous previous studies of oceanic chaos and its inverse cascades (e.g., Leroux et al., 2018; Sérazin et al., 2018; Cravatte et al., 2021; Hogg et al., 2022). Similarly to ensemble methods in climate science (e.g., Maher et al., 2021), the member simulations are generated from a single model under identical external forcing but perturbed initial conditions, allowing one to disentangle intrinsic from forced (i.e., atmosphere-driven) variability. Emphasis is lain on the period from 1995 to 2015 and particularly the interannual band, where oceanic effects are a leading cause for wobble excitation, along with continental hydrology (Gross et al., 2003; Wińska, 2022; Börger et al., 2023). We separate the ensemble OAM estimates into forced and intrinsic components and assess their contributions to the observed wobble excitation in the presence of mass change and motion elsewhere in the Earth system. This in turn calls for a consideration of all known climatic excitations of polar motion, as described below.

## 2 Methods

### 2.1 Basic framework and data

Polar motion is a two-dimensional quantity, conventionally expressed as complex-valued angular distance $\hat{p} = \hat{p}(t) = p_1 + ip_2$, where t indicates time and subscripts 1 and 2 refer to longitudes of 0 and 90° E in the terrestrial reference frame. Under the conservation of angular momentum (Munk and MacDonald, 1960; Moritz and Mueller, 1987), the equa-

tion governing the excitation of observed variations in $\hat{p}$ is (e.g., Brzeziński, 1992)

$$ \hat{p} + i\hat{\sigma}_c^{-1}\dot{\hat{p}} = \hat{\chi}, \quad (1) $$

where $\hat{\chi} = \chi_1 + i\chi_2$ is the polar motion excitation function (Barnes et al., 1983; Gross, 2007), the dot denotes the time derivative, and $\hat{\sigma}_c = 2\pi(1 + i/2Q_c)/T_c$ is the complex-valued Chandler frequency with standard values for period $T_c = 433.0$ d and quality factor $Q_c = 179$ (Gross et al., 2003). Given proper observations of $\hat{p}$, the deconvolution on the left side of Eq. (1) yields the “geodetic” (or observed) excitation of polar motion, labeled as $\hat{\chi}^*$ throughout this work. We specifically use rotation data from the daily sampled SPACE2018 series (Ratcliff and Gross, 2019) that employs a Kalman filter to adjust measurements from the primary space geodetic techniques. We remove long-period tidal effects in polar motion using the conventional model of Peitit and Luzum (2010), deconvolve $\hat{p}$, and finally average the daily $\hat{\chi}^*$ values to monthly estimates.

The practical problem defined by Eq. (1) lies in identifying the causes for the observed excitation from angular momentum changes in geophysical fluids. As is customary (e.g., Wahr, 1982), we consider excitations both due to the redistribution of mass and particle movement; that is,

$$ \hat{\chi} = \hat{\chi}^m + \hat{\chi}^v, \quad (2) $$

where the mass or pressure term (superscript m) accounts for perturbations in inertia products and the motion term (superscript v for velocity) represents the fluid's relative angular momentum. Formulae for computing $\hat{\chi}^{m,v}$ as area integrals over globally distributed mass changes and horizontal velocities are given in Appendix A. Note that these geophysical excitation functions are also called “effective angular momentum functions” (Barnes et al., 1983; Dickman, 2003), where the term “effective” pertains to the non-rigid component of the Earth's response to the excitation process. Below we occasionally refer to the oceanic $\hat{\chi}$ estimate as the OAM function, thus omitting the prefix “effective” for convenience.

In comparisons with $\hat{\chi}^*$, we consider excitations from four main climate system components, comprising the atmosphere (denoted to as $\hat{\chi}^A$), the ocean ($\hat{\chi}^O$), the continental hydrosphere ($\hat{\chi}^H$), and the cryosphere as represented by the Greenland and Antarctica ice sheets ($\hat{\chi}^C$). The satellite data underlying $\hat{\chi}^H$ and $\hat{\chi}^C$, along with the ensemble approach to oceanic excitation (Sect. 2.2), constrain us to have monthly sampling. For consistency with the ocean model output analyzed here, we use atmospheric angular momentum time series from the ERA-Interim reanalysis of the European Centre for Medium-Range Weather Forecasts (Dee et al., 2011). The $\hat{\chi}^A$ series is an extension of a previous calculation by us (Schindelegger et al., 2013) and accounts for excitation signals associated with the ocean's inverted barometer (IB) response to atmospheric pressure loading. Figure 1 shows that,
on interannual timescales (here periods ≥ 14 months<fnref n="1" />), the atmosphere is somewhat less effective than the ocean in exciting Earth’s wobbles, particularly through $\chi_1$ (see Gross et al., 2003). Peak-to-peak amplitudes in $\chi_2^A$ occasionally attain 15 mas (milliarcseconds) and are largely governed by mass redistribution over Northern Hemisphere landmasses (Nastula and Salstein, 1999).

Contributions to $\hat{\chi}^m$ from continental hydrology – or, equivalently, terrestrial water storage (TWS) – and the cryosphere are deduced from time-variable gravity fields, as practiced in previous works (e.g., Chen et al., 2013; Adhikari and Ivins, 2016; Nastula et al., 2019; Göttl et al., 2021). To cover the same 21-year period (1995–2015) as for our oceanic excitation series, we use a long-term reconstruction of monthly gravity fields based on satellite tracking data (Löcher and Kusche, 2021, updated version). The approach pursued for this reconstruction is to expand surface mass changes from the Gravity Recovery and Climate Experiment/Follow On (GRACE/-FO; Tapley et al., 2019) mission into empirical orthogonal functions (EOFs) and use the leading EOFs as basis functions for gravity field modeling in the processing of the tracking observations. The solution has no temporal gaps and preserves the spatial resolution of GRACE/-FO gravity fields (i.e., wavelengths of ~300 km) also for times prior to GRACE. We synthesize the monthly varying Stokes coefficients as equivalent water heights across Earth’s surface, mask out the ocean, and partition the mass anomalies into changes over the Greenland Ice Sheet, the Antarctic Ice Sheet, and land. Spatial leakage effects near continental boundaries and ice sheet margins are mitigated using an iterative forward modeling algorithm (Chen et al., 2015); see Appendix B for details. Figure 1 illustrates that hydrological mass changes drive large (~10–15 mas) interannual to intradecadal variations in polar motion excitation (Adhikari and Ivins, 2016; Meyrath and Van Dam, 2016). Detrended time series for cryospheric sources ($\hat{\chi}^C$) reveal residual long-term drifts that bear on the accelerated melting of polar ice sheets (Chen et al., 2013; Deng et al., 2021).

In addition to $\hat{\chi}^{A,O,H,C}$, we account for non-uniform mass changes in the ocean due to gravitational attraction and loading (GAL), defining the excitation function $\hat{\chi}^L$. Processes described by GAL are the continent-ocean net mass transfer, the gravitational attraction of water toward mass anomalies anywhere in the Earth system, the crustal deformation under these loads, and the ensuing changes in Earth’s gravitational field (e.g., Clarke et al., 2005). The surface mass variations arising from these effects are long-wavelength in nature, leading to small but non-negligible angular momen-

tum changes (Quinn et al., 2015; Adhikari and Ivins, 2016, Fig. 1). We calculate the GAL effects on ocean mass and distribution by solving the sea level equation (Farrell and Clark, 1976), assuming elastic (i.e., instantaneous) deformation of a rotating Earth and a purely static adjustment in the ocean. The considered loads are gridded monthly anomalies of TWS, cryospheric mass, and IB-corrected atmospheric pressure, taken from the same sources as those adopted for our $\hat{\chi}^{A,H,C}$ estimates; see Appendix C.

## 2.2 Oceanic excitation

We study the externally forced and chaotic intrinsic components of oceanic excitation based on a 50-member ensemble of eddy-permitting ocean/sea ice simulations from the OceaniC Chaos – ImPacts, strUcture, predicTability (OC-CIPUT) project (Penduff et al., 2014). The OCCIPUT large ensemble was generated using a global 1/4° tri-polar grid setup of the Nucleus for European Modelling of the Ocean engine for Boussinesq ocean dynamics and thermodynamics. The full ensemble hindcast spans 56 years (1960–2015). Following a common 21-year spin-up simulation, a small stochastic perturbation was applied in the density equation of each member throughout the year 1960 (Brankart et al., 2015), triggering the dispersion of the individual simulations and subsequent saturation of the ensemble spread (which occurs more quickly in turbulent, eddy-rich areas; Penduff et al., 2014; Bessières et al., 2017). These 50 realizations of the time-evolving ocean state were then integrated forward from their respective, slightly perturbed initial conditions and driven by the same 6-hourly atmospheric data from the DRAKKAR forcing set DFS 5.2 (Dussin et al., 2016), which itself is based on ERA-Interim. As the ensemble members are identical in terms of the numerical model, boundary conditions, and atmospheric forcing, their spread can be solely attributed to nonlinear ocean dynamics. Note also that the existence, origin, and spatiotemporal scales of this intrinsic variability is consistent across eddy-permitting and eddy-resolving ocean models (e.g., Sérazin et al., 2015).

To keep the involved data volumes manageable, we analyze the OCCIPUT ensemble output from 1995 to 2015 (cf. Zhao et al., 2021) and focus on monthly mean values. Specifically, we evaluate the mass and motion terms of $\hat{\chi}^O$ (see Appendix A) from monthly averages of ocean bottom pressure, $p_b = p_b(\mathbf{x}, t)$, and vertically averaged northward and eastward currents, $\bar{u}(\mathbf{x}, t)$ and $\bar{v}(\mathbf{x}, t)$, where $\mathbf{x}$ represents horizontal space. Each of the three variables $\xi(\mathbf{x}, t) \in (p_b, \bar{u}, \bar{v})$ consists of externally forced and intrinsic signals; i.e.,

$$ \xi(\mathbf{x}, t, \mu) = \xi^f(\mathbf{x}, t) + \xi^i(\mathbf{x}, t, \mu), \qquad (3) $$

where $\mu$ is a running index over the 50 ensemble members, superscript f indicates the forced variability common to all members, and superscript i identifies the intrinsic signal. As in previous studies (see, in particular, Hogg et al., 2022), we

<footnote n="1">¹This cutoff period was chosen to effectively suppress the seasonal components in the excitation time series. For studying the excitation of the Chandler wobble at $T_c = 433.0$ d, which is not attempted here, a different approach to analysis would be required (see Gross, 2000).</footnote>
<figure><img src="image_4.png" /></figure>

**Figure 1.** Modeled excitation $\hat{\chi} = \chi_1 + i\chi_2$ of interannual polar motion (in milliarcseconds) from 1995 through 2015, separated into contributions from atmosphere ($\hat{\chi}^A$, light blue), ocean ($\hat{\chi}^O$, blue), terrestrial hydrology ($\hat{\chi}^H$, dark blue), the Antarctic Ice Sheet ($\hat{\chi}^{C,Ant}$, red), the Greenland Ice Sheet ($\hat{\chi}^{C,Gre}$, orange), and gravitational attraction and loading ($\hat{\chi}^L$, cream). Interannual signals were obtained by removing each time series' trend and seasonality and low-pass filtering the residual series with a cutoff period of 14 months. We clipped 4 months both at the beginning and at the end of the time series to avoid filter artifacts. The $\hat{\chi}^O$ estimate is from an arbitrary OCCIPUT ensemble member. For traceability, we compare our estimates of $\hat{\chi}^H$, both unfiltered and low-pass filtered, with GRACE-only $\hat{\chi}^H$ time series in the Supplement (Figs. S1 and S2).

estimate the forced component at a given month and location by the ensemble mean $\xi^f(\mathbf{x}, t) = \langle \xi(\mathbf{x}, t, \mu) \rangle$, adopting $\langle \cdot \rangle$ as the ensemble-mean operator. Upon separating signals as per Eq. (3), we compute angular momentum time series (Appendix A) for forced and intrinsic contributions to ($p_b, \bar{u}, \bar{v}$) and for their sum. Hence, intrinsic variability in oceanic excitation (abbreviated as $\hat{\chi}^{O,i}$) is distinguished from forced variability ($\hat{\chi}^{O,f}$) at the level of spatially distributed dynamical variables and is not based on excitation time series.

Throughout this work, $p_b$ denotes ocean bottom pressure anomalies, obtained by removing for each ensemble member a common time-mean reference pressure field and a spatial mean $p_b$ value per month and ensemble member. The latter correction prevents small biases or drifts in $\hat{\chi}^O$ caused by spurious global ocean mass changes in OCCIPUT under the Boussinesq approximation. Legitimate contributions to the globally averaged $p_b$, which would be removed by the correction, pertain to the freshwater flux into the ocean from atmospheric and continental sources. However, this term is accounted for through the framework of GAL and thus in $\hat{\chi}^L$. Note that the OCCIPUT simulations do not consider static effects, such as long-period tides or the pole tide; thus these phenomena are absent from $\hat{\chi}^{O,f}$.

# 3 Results and discussion

## 3.1 Overview in spectral domain

A natural first step is to assess the magnitude and importance of the intrinsic oceanic excitation in relation to the externally forced component across a range of frequencies. In Figs. 2 and 3, we show amplitude spectra (symbol A) for $\hat{\chi}^{O,i}$ (50 realizations) and $\hat{\chi}^{O,f}$ (1 realization), separated into mass and motion terms and the two coordinate directions. For comparison purposes, we form chaotic amplitude fractions, $R = A_i/A_f$, representing the ratio of intrinsic over forced amplitude per frequency (not to be mixed up with time-invariant fractions of variance considered in previous OCCIPUT analyses, e.g., Carret et al., 2021; Hogg et al., 2022). Values of R are computed separately for (a) annual, semiannual, and terannual harmonics extracted from $\hat{\chi}^{O,i}$ and $\hat{\chi}^{O,f}$ and (b) the deseasonalized excitation functions. We smooth the estimated fractions under item (b) using a moving window, with window sizes decreasing (essentially linearly) from 21 frequency bins at the high end of the spectrum to 1 bin at the longest period.

As apparent in Figs. 2 and 3 (panels a and c), the intrinsic oceanic excitation is embedded in a relatively flat spectrum, bounded by amplitudes of ~1.5–2.0 mas in the mass term and ~0.7–1.0 mas in the motion term. This spectral characteristic renders the chaotic excitation comparable in magni-
<figure><img src="image_5.png" /></figure>

**Figure 2.** Spectral representation of the OCCIPUT-based estimates of intrinsic and forced oceanic excitation of polar motion, 1995–2015. Shown are amplitudes $A_i$ and $A_f$ (in milliarcseconds) of the intrinsic and forced oceanic mass term, $\hat{\chi}^m$, in (a) $\chi_1$ and (c) $\chi_2$, deduced from a 512-point fast Fourier transform using Welch's method (Welch, 1967). Black vertical bars indicate 68% confidence limits for selected frequencies. Panels (b) and (d) illustrate spectrally smoothed versions of the chaotic amplitude fraction $R = A_i/A_f$ (dimensionless), estimated from the excitation time series with seasonal terms removed. Solid black lines mark the respective medians across the 50 realizations of $R$. Statistics for the annual, semiannual, and terannual terms – computed separately based on seasonal cycle fits to $\hat{\chi}^{O,i}$ and $\hat{\chi}^{O,f}$ – are summarized by the purple whiskers; box edges indicate the 25th and 75th percentiles, and vertical bars extend to the minimum/maximum values of seasonal $R$.

tude to the forced variability at most but not all frequencies. Notable exceptions occur in $\chi_2^v$ and wherever $\hat{\chi}^{O,f}$ features peak in seasonality. Values of $R$ are close to 0.15 for the dominant ($\sim 3.5$ mas) annual harmonic in both $\chi_1^v$ and $\chi_2^v$, along with similarly small fractions at periods of 1 and $1/3$ year in $\chi_1^m$. We can tentatively interpret the results for $\chi_1^m$ in terms of intrinsic vs. forced ratios of $p_b$ signals from previous OCCIPUT analyses (Zhao et al., 2021). In particular, Fig. 3 of Zhao et al. (2021) places large intrinsic contributions to the seasonal $p_b$ cycle over western boundary currents and the eddying Southern Ocean, with little sign of long-wavelength effects that would drive strong variability in $\hat{\chi}^{O,i}$.

Switching from seasonal harmonics to the broadband spectral behavior in Figs. 2 and 3, we observe a gradual increase in the chaotic amplitude fraction from $\sim 0.3$ to $\sim 0.6$ as periods lengthen from 60 to 120 d. This is followed at intra-annual frequencies by somewhat of a plateau of the median $R$ near 0.9 in the mass term, along with a broadening of the ensemble spread to $R \approx 0.3-2.0$ in all OAM components. The prominence of $\hat{\chi}^{O,i}$ at these periods (roughly 120 to 300 d) coincides with the dominant timescale for the emergence of mesoscale eddies (Chelton et al., 2011) and the ensuing nonlinear inverse cascades (e.g., Arbic et al., 2014; O'Rourke et al., 2018; Sérazin et al., 2018). Concurrently, the forced oceanic excitation is relatively weak at intra-annual

periods (especially in $\chi_1^O$; Figs. 2a and 3a), making this specific band another attractive (future) target for exploring the polar motion imprint of oceanic chaos. We note, though, that the intra-annual wobble excitation is largely dominated by atmospheric mass redistribution effects (Gross et al., 2003).

For the interannual signals central to this work, the spectral representations in Figs. 2 and 3 indicate appreciable power in $\hat{\chi}^{O,i}$ relative to $\hat{\chi}^{O,f}$, even at sub-decadal periods and particularly in $\chi_1^O$. Some low-frequency, intrinsic oceanic process must thus be at work at scales large enough to perturb the globally integrated OAM quantities (see Sect. 3.2). More specific to Figs. 2 and 3, we observe a number of troughs in the forced amplitude $A_f$ that in turn afford peaks in the chaotic amplitude fractions. Such maxima in $R$ are seen, for example, for periods $T \approx 1.2-1.3$ years in all four $\hat{\chi}^O$ components or near 3 years in the two mass terms (spectral uncertainty is generally large, though). By contrast, the forced component exceeds the intrinsic amplitude by a factor of about 2 at $T \approx 1.6$ years in both mass and motion terms of $\chi_2^O$. One can therefore expect any time series of interannual $\chi_2^O$ signals to feature stronger deterministic variability than $\chi_1^O$.

For a condensed view of Figs. 2 and 3, we have additionally formed ratios between intrinsic and forced standard deviations, i.e., $\sigma_i/\sigma_f$, calculated from $\hat{\chi}^{O,i}$ and $\hat{\chi}^{O,f}$ time series.
<figure><img src="image_6.png" /><figcaption>**Figure 3.** As in Fig. 2 but for the two components of the motion term $\hat{\chi}^v$.</figcaption></figure>

Summing up mass and motion terms, and filtering the partic-
ipating time series to periods ≥ 14 months, we obtain σᵢ/σf
values of 0.50–0.83 (median 0.69) for χ₁º and 0.28–0.50
(median 0.36) for χ₂º. Similar calculations for the complex-
valued quantity χº yield σᵢ/σf ratios in the range of 0.38–
0.60 and a median at 0.47. When taking σᵢ relative to the vari-
ability of the total oceanic excitation, the median is 0.46 in
the case of χº (range 0.43–0.50). These numbers emphasize
the importance of chaotic contributions to equatorial OAM
changes on timescales much longer than typically associated
with the life cycles of mesoscale eddies.

3.2 Modal decomposition

Analysis of the space–time variability in modeled intrinsic
$p_b$ fields, denoted as $p_b^i$, can reveal more about the prove-
nance of the chaotic OAM signals apparent in Figs. 2 and 3.
Following Zhao et al. (2023), we have decomposed each
of the 50 realizations of $p_b^i$ into real-valued EOF patterns
and the associated principal components. We find a common
leading EOF (mode 1 hereafter), accounting for 16 % to 22 %
of the global $p_b^i$ variance. This particular mode, displayed
exemplarily for one ensemble member in Fig. 4a, exhibits a
distinct global structure, characterized by covariability across
the Atlantic, the Arctic, and partly the Indian Ocean, while
much of the Southern Ocean oscillates with opposite sign.
The suggested origin for this mode lies in the Drake Passage
(Zhao et al., 2023). Here, strong eddies and random current
fluctuations, possibly involving interaction with topography
(Provost et al., 2011), are thought to induce sharp bottom
pressure gradients, evident in Fig. 4a between Cape Horn

and the South Shetland Islands. The basin-wide $p_b^i$ anomalies seen elsewhere in the ocean would then be the result of fast barotropic adjustment to these local gradients along $H/f$ contours (where $H$ is the local water depth and $f$ is the Coriolis parameter; Zhao et al., 2023).

The described processes, while different in nature to the
nonlinear inverse energy cascades noted in Sect. 1, are
unique in that they produce a global pattern of mass redis-
tribution with opposite polarity between ocean basins. Such
large-scale see-saws of mass have previously been shown to
be an effective source of polar motion excitation (cf. Afroosa
et al., 2021). We have thus synthesized monthly $p_b^i$ grids
from mode 1 alone and deduced the corresponding oceanic
excitation (that is, the oceanic mass term; Appendix A). The
derived estimates of $\hat{\chi}^{O,i}$ are compared in Table 1 against the
full $\hat{\chi}^{O,i}$ functions from OCCIPUT in terms of the percentage
of variance explained (PVE). This simple evaluation demon-
strates that mode 1 is indeed the leading driver for the mass
component of intrinsic oceanic excitation across timescales
and for both $\chi_1$ and $\chi_2$. Specifically, at the interannual pe-
riods of interest in this work, the OAM changes associated
with $p_b^i$ mode 1 explain 39 % to 97 % of the variance in the
$\hat{\chi}^{O,i}$ mass term, with a median PVE of 63 % across the 50
ensemble members. These numbers encapsulate the idea that
localized mesoscale processes somewhere in the world ocean
govern variability in a truly global and geophysically relevant
quantity.
<table><caption>**Table 1.** Contribution of mode 1 $p_b^i$ variability to the intrinsic component of the oceanic mass term, 1995–2015.ª</caption><thead><tr><th></th><th></th><th>χ<sub>1</sub></th><th>χ<sub>2</sub></th><th>χ&#770; = χ<sub>1</sub> + iχ<sub>2</sub></th></tr></thead><tbody><tr><td rowspan="3">PVE in χ<sub>j</sub><sup>O,i</sup> by EOF mode 1</td><td>Min</td><td>38.4</td><td>35.0</td><td>40.4</td></tr><tr><td>Median</td><td>56.6</td><td>49.6</td><td>52.3</td></tr><tr><td>Max</td><td>97.4</td><td>97.9</td><td>97.7</td></tr><tr><td rowspan="3">PVE in low-frequency χ<sub>j</sub><sup>O,i</sup> by EOF mode 1<sup>b</sup></td><td>Min</td><td>35.6</td><td>28.3</td><td>39.0</td></tr><tr><td>Median</td><td>67.8</td><td>63.9</td><td>63.3</td></tr><tr><td>Max</td><td>96.7</td><td>97.1</td><td>96.8</td></tr></tbody></table>

<figure><img src="image_7.png" /><figcaption>**Figure 4.** Attribution of the mass term of $\hat{\chi}^{O,i}$ to global-scale intrinsic bottom pressure variability, $p_b^i$. Shown are (a) the leading EOF mode of $p_b^i$ variability (in centimeters) for ensemble member 47 and (b, c) the corresponding equatorial oceanic excitation signal (dashed black curve), obtained from synthesizing mode 1 in time and space and computing the associated OAM time series. The full intrinsic oceanic excitation, derived from $p_b^i$ without modal partitioning, is plotted in gray. We selected member 47 for this comparison, as its statistics in Table 1 are close to the median percentage of variance explained by mode 1 in $\hat{\chi}^{O,i}$. Note that other ensemble members suggest a mode 1 structure very similar to the one shown in panel (a); see also Zhao et al. (2023).</figcaption></figure>

## 3.3 Interannual variability

We now examine the role of oceanic chaos in the low-frequency equatorial excitation budget over 1995–2015. To that end, Fig. 5 shows the residual geodetic excitation, obtained by subtracting from $\hat{\chi}^*$ modeled contributions from the atmosphere, terrestrial hydrology, the two ice sheets, and GAL (cf. the respective time series in Fig. 1). Together, these non-oceanic sources account for 29.7% of the variance in the full equatorial geodetic excitation (Table 2), governed in large parts by the PVE in $\chi_1$ (51.6%). The residual is to be compared with 51 variants of oceanic excitation from OCCIPUT, comprising both the forced component ($\hat{\chi}^{O,f}$) and 50 realizations of $\hat{\chi}^O$ (i.e., $\hat{\chi}^{O,f} + \hat{\chi}^{O,i}$). Intrinsic processes clearly add extra variability to estimates of $\hat{\chi}^O$ relative to $\hat{\chi}^{O,f}$, as conveyed by an ensemble spread of nearly ±5 mas throughout the analysis period. This is an important result that implies some leeway in how to evaluate and interpret the observed excitation of Earth’s wobbles on interannual timescales.

Adopting the forced component from OCCIPUT as a middle ground for an estimate of $\hat{\chi}^O$ reduces the root mean square (RMS) of the residual equatorial geodetic excitation from 9.6 to 6.7 mas (PVE = 50.8%). Similarly, we can select from any of the 50 ensemble members as valid representations of the time-evolving ocean state and thus $\hat{\chi}^O$. This approach to evaluation leads to a considerable spread in the statistics in Table 2, involving a few ensemble members that exhibit slightly better skill than the forced component (minimum $\hat{\chi}$ RMS of 6.5 mas, PVE = 53.7%) but also a number of cases that reduce little variance (only 1.5 mas in RMS) in the residual geodetic excitation. A total of 15 members even yield negative PVE values in $\chi_1$. However, these statistics are not to be taken to define absolute bounds on the efficacy of OAM in exciting low-frequency polar motion, as the intrinsic variability in any of the ensemble members is merely a random realization of the potential intrinsic variability in the real world. In fact, oceanic chaos may act to produce yet more time series of $\hat{\chi}^O$ that are also consistently close to the

<footnote>ª Values are percentage of variance explained (PVE) by a monthly varying $\hat{\chi}^{O,i}$ mass term computed from mode 1 vs. the full intrinsic mass term. Subscript *j* is used to discriminate between $\chi_1$, $\chi_2$, and $\hat{\chi}$.</footnote>

<footnote>ᵇ Participating time series filtered to periods ≥ 14 months.</footnote>
<figure><img src="image_8.png" /><figcaption>**Figure 5.** The low-frequency polar motion excitation budget from 1995 to 2015 in the presence of oceanic chaos. Shown are (as dashed black lines) the (a) $\chi_1$ and (b) $\chi_2$ components of geodetic excitation ($\hat{\chi}^*$) – corrected for contributions from the atmosphere ($\hat{\chi}^A$), terrestrial hydrology ($\hat{\chi}^H$), ice sheets ($\hat{\chi}^C$), and GAL ($\hat{\chi}^L$) – compared with 50 estimates of oceanic excitation from OCCIPUT ($\hat{\chi}^O$, gray curves), the OCCIPUT ensemble mean ($\hat{\chi}^{O,f}$, blue curve), and a hypothetical ensemble member (purple curve) producing the best match with the corrected geodetic excitation. All time series have been detrended, filtered to periods $T \ge 14$ months, and cut back by 4 months at the two end points.</figcaption></figure>

<table><caption>**Table 2.** Observed vs. modeled excitation of interannual polar motion, 1995–2015. a, b</caption><thead><tr><th></th><th></th><th>χ<sub>1</sub></th><th>χ<sub>2</sub></th><th>χ̂ = χ<sub>1</sub> + iχ<sub>2</sub></th></tr></thead><tbody><tr><td>RMS of χ<sub>j</sub><sup>*</sup></td><td></td><td>7.5</td><td>8.6</td><td>11.4</td></tr><tr><td>PVE in χ<sub>j</sub><sup>*</sup> by non-oceanic sources,<br>χ<sub>j</sub><sup>A</sup> + χ<sub>j</sub><sup>H</sup> + χ<sub>j</sub><sup>C</sup> + χ<sub>j</sub><sup>L</sup></td><td></td><td>51.6 (5.2)</td><td>13.3 (8.0)</td><td>29.7 (9.6)</td></tr><tr><td rowspan="2">PVE in χ<sub>j</sub><sup>*</sup> by OCCIPUT, χ<sub>j</sub><sup>O</sup></td><td>Forced</td><td>18.5 (6.8)</td><td>35.4 (6.9)</td><td>28.2 (9.7)</td></tr><tr><td>Ensemble members</td><td>[-17.5–27.4] ([6.4–8.1])</td><td>[10.2–44.7] ([6.4–8.2])</td><td>[5.3–34.2] ([9.3–11.1])</td></tr><tr><td rowspan="3">PVE by OCCIPUT χ<sub>j</sub><sup>O</sup> in residual series,<br>χ<sub>j</sub><sup>*</sup> – χ<sub>j</sub><sup>A</sup> – χ<sub>j</sub><sup>H</sup> – χ<sub>j</sub><sup>C</sup> – χ<sub>j</sub><sup>L</sup></td><td>Forced</td><td>23.9 (4.5)</td><td>62.1 (5.0)</td><td>50.8 (6.7)</td></tr><tr><td>Ensemble members</td><td>[-24.4–35.1] ([4.2–5.8])</td><td>[44.8–64.6] ([4.8–6.0])</td><td>[29.1–53.7] ([6.5–8.1])</td></tr><tr><td>Hypothetical member<sup>c</sup></td><td>77.4 (2.5)</td><td>87.4 (2.9)</td><td>84.0 (3.8)</td></tr><tr><td rowspan="3">PVE in χ<sub>j</sub><sup>*</sup> by sum of excitation processes,<br>χ<sub>j</sub><sup>A</sup> + χ<sub>j</sub><sup>O</sup> + χ<sub>j</sub><sup>H</sup> + χ<sub>j</sub><sup>C</sup> + χ<sub>j</sub><sup>L</sup></td><td>Forced</td><td>63.1</td><td>67.1</td><td>65.4</td></tr><tr><td>Ensemble members</td><td>[39.8–68.6]</td><td>[52.1–69.3]</td><td>[50.2–67.5]</td></tr><tr><td>Hypothetical member<sup>c</sup></td><td>89.0</td><td>89.1</td><td>89.0</td></tr></tbody></table>

residual geodetic excitation while still lying within the ensemble spread in Fig. 5.

We have constructed such a hypothetical – but still legitimate – equatorial OAM function by identifying, month per month and in the complex plane, the particular OCCIPUT

$\hat{\chi}^O$ sample closest to the reduced observation (i.e., $\hat{\chi}^*$ minus non-oceanic sources). As one switches from one low-pass filtered member to the other, small jumps are incurred, requiring an additional round of smoothing to periods $\ge 14$ months. As per design, the resulting time series follow the residual

<footnote>a Values are percentage of variance explained (PVE) and the corresponding RMS of residuals (milliarcseconds) in parentheses, except for the first line. Entries for ensemble members apply to the sum of the forced and intrinsic signals. All time series have been detrended and filtered to periods $T \ge 14$ months. Subscript *j* in the intermediate headers is used to discriminate between $\chi_1$, $\chi_2$, and $\dot{\chi}$.</footnote>

<footnote>b Replacing the SPACE2018 series in the analysis by other, more recent Earth rotation determinations changes the tabulated PVE values by 0.2% at most.</footnote>

<footnote>c Hypothetical member drawn from the $\hat{\chi}^O$ ensemble based on a constraint of minimum distance to the residual observed excitation; see the main text.</footnote>
geodetic excitation more tightly than $\hat{\chi}^{O,f}$ or any of the en-
semble members (Fig. 5). The improved correspondence is
also reflected in the PVE by the hypothetical $\hat{\chi}^O$ (84.0% of
the reduced geodetic excitation), leaving only 3.8 mas of in-
terannual RMS variability unexplained. Likewise, using the
hypothetical member instead of the ensemble mean in evalu-
ating the sum of all geophysical fluid excitations against $\hat{\chi}^*$
increases the PVE from 65.4% to 89.0%. For comparison,
selecting from the original OCCIPUT members to form the
total surface layer excitation yields PVE values from 50.2% to
67.5%. Overall, these statistics illustrate the sensitivity of
such budget analysis to the representation of chaotic OAM
signals.

As evident in Fig. 5, the 3.8 mas RMS signal unaccounted for in the observations is mostly made up by discrepancies between geodetic and geophysical excitations before early 2002 (i.e., the launch of GRACE), large differences throughout the year 2005 in $\chi_1$ (to a lesser extent also 2009/2010), and a ~7 mas cusp around the turn of 2011/2012 in $\chi_2$. It is possible that part of the discrepancies arise from errors in OCCIPUT, related to, e.g., imperfectly modeled energetic circulations in the Southern Ocean (Ponte and Piecuch, 2014; Harker et al., 2021). A more likely source, though, is residual errors in the adopted gravity fields that project onto the inter-annual excitation budget, particularly through TWS changes and the associated fluctuations in $\hat{\chi}^H$. Comparing our hydro-logical time series in Fig. 1 with GRACE-only estimates of $\hat{\chi}^H$ in Fig. S2 does indeed point to differences over some of the aforementioned periods (2009/2010 in $\chi_1$; 2011/2012 in $\chi_2$). It therefore appears that the specifics of the gravity field reconstruction – including the limited spectral sensi-tivity of the involved satellites and the temporal parameteri-zation in terms of GRACE-/FO EOFs (Löcher and Kusche, 2021) – occasionally compromise the accuracy of the derived hydrological angular momentum changes. More detailed ex-ploration of these factors is left for future studies.

Although the geophysical models and data used in this
study are necessarily imperfect, it is worth inspecting the
excitation budget from yet another angle. In Fig. 6, we
show amplitude spectra for the two components of $\hat{\chi}^*$ mi-
nus all climate-system-related excitations, also including,
for present purposes, the forced oceanic component ($\hat{\chi}^{O,f}$).
Spectral power in these residuals – while not well resolved
within our 21-year analysis window – may be due to the
said errors in the adopted angular momentum series. Alter-
natively, residual power may point to secondary excitation
processes left out of the considerations and coming from,
for example, the Earth’s core (Chen et al., 2019; Ding et al.,
2019; Kuang et al., 2019). Inferences along the latter lines
are, however, complicated by the difficulty in accounting for
the intrinsic component of oceanic excitation, $\hat{\chi}^{O,i}$, as ex-
pressed by the gray background spectra in Fig. 6. In detail,
the OCCIPUT ensemble suggests a likelihood for $\hat{\chi}^{O,i}$ to
attain amplitudes between 0 and 1.5 mas across all interan-
nual and sub-decadal frequencies (see Figs. 2 and 3). Always

selecting the most energetic ensemble member at each fre-
quency (blue lines in Fig. 6) raises the upper bound for in-
trinsic effects even further, to about 2 mas.

A 1.0–1.5 mas uncertainty in modeled wobble excitation due to oceanic chaos is generally comparable in size to the residual geodetic excitation across most frequencies (dashed black curves in Fig. 6). This prompts the question as to what extent it is possible to conclusively link such residuals to the dynamics of Earth’s deep interior. For instance, Chen et al. (2019) suggested that an ~2 mas oscillation with a period of T = 5.9 years (left after correcting the observed excitation over 1962–2018 for modeled atmospheric, oceanic, and hydrological effects) may be caused by the exchange of equatorial angular momentum between the outer core and the mantle via the electromagnetic torque (Kuang et al., 2019). However, analyses of satellite-derived degree-2 Stokes coefficients (Rosat et al., 2021) fail to support this idea, while our Fig. 6 illustrates that a considerable portion of any interannual $\hat{\chi}^*$ residual may simply be due to mis-modeled oceanic effects. Furthermore, the possibility of intrinsic OAM signals exciting polar motion at the level of 1.0–1.5 mas even at sub-decadal periods (T = 5–10 years) adds uncertainty to attempts to detect the inner-core wobble (ICW; Rochester and Crossley, 2009) from rotation data. Recent studies (e.g., Ding et al., 2019) have made the case for an 8.7-year prograde harmonic in residual polar motion to be the signature of the ICW, at an amplitude of just 2.67 ± 0.04 mas. However, the OAM product used by Ding et al. (2019) to remove oceanic contributions from the observations did not account for intrinsic variability, implying that the error bars on the inferred ICW amplitude are likely too optimistic. Substantiating the evidence for core effects in polar motion data on sub-decadal timescales will most certainly require further improvements in modeling surface layer processes, comprising in particular oceanic and hydrological angular momentum changes (Rosat et al., 2021; Rekier et al., 2022).

4 Summary and conclusions

The 50-member OCCIPUT ensemble investigated in this
study advocates for a prominent role of oceanic intrinsic vari-
ability in equatorial OAM changes from intraseasonal to in-
terannual timescales. This marks the arrival of a new player
on the scene of Earth rotation variations, with implications
for tests of closure in the polar motion excitation budget and
as yet unmodeled geophysical processes (Table 2 and Fig. 6).
In terms of attribution, we have been able to relate most of
the mass-driven signals in $\hat{\chi}^{O,i}$ (median PVE = 63 % at pe-
riods $\ge$ 14 months) to one particular mode of intrinsic $p_b$
variations that likely manifests the free barotropic adjust-
ment of the main ocean basins to nonlinearly induced $p_b$
anomalies near the Drake Passage (cf. Fig. 4 and Zhao et al.,
2023). As in many other parts of the Southern Ocean, low-
frequency changes in eddy activity around the Drake Pas-
<figure><img src="image_9.png" /><figcaption>Figure 6. Amplitude spectra (in milliarcseconds) of $\hat{\chi} = \chi_1 + i\chi_2$ over 1995–2015 for the intrinsic oceanic excitation (50 realizations of $\hat{\chi}^{0,i}$, gray lines), in comparison to the residual geodetic excitation (dotted black line), obtained by subtracting from $\hat{\chi}^*$ estimated contributions from the atmosphere, terrestrial hydrology, ice sheets, GAL, and the forced oceanic excitation ($\hat{\chi}^{0,f}$); see the expression in the legend. The upper envelope of the intrinsic terms is represented by the blue line.</figcaption></figure>

sage are strongly random in character (Hogg et al., 2022), thus providing a plausible explanation as to why variability in $\hat{\chi}^{0,i}$ is a considerable fraction (43%–50%; see Sect. 3.1) of the total oceanic excitation of polar motion even on interannual timescales. A caveat to be acknowledged, though, is that all of these inferences are based on simulations. The credibility of OCCIPUT-based OAM series is, however, clear from Table 2, and comparisons with monthly GRACE-/FO fields have shown that the leading $p_b^i$ mode inherent to the ensemble could be present in satellite observations (Zhao et al., 2023). In addition, the Southern Ocean eddy field in the 1/4° OCCIPUT runs is slightly less intense than in altimetry data (Carret et al., 2021; Hogg et al., 2022), suggesting that our estimates for the chaotic contribution to wobble excitation are of a rather conservative nature.

Our results underscore the challenge of accurately representing oceanic effects for polar motion studies. In general, high horizontal model resolution ($\le 1/4°$) is preferred in accounting for wind-driven OAM changes (Harker et al., 2021). However, a single high-resolution baroclinic model run would also produce one specific realization of the intrinsic component that is unlikely to match the phase of the observed one. Conversely, pursuing a probabilistic approach with a large eddy-permitting ensemble such as OCCIPUT may be computationally impracticable, especially on an operational basis and at the (sub-)daily resolution now common for rotation studies (Dobslaw and Dill, 2018). Better constraining the chaotic component of $\hat{\chi}^{0,i}$ will most certainly require an eddying ocean model assimilating, amongst other observations, altimetric sea level anomalies and GRACE-/FO estimates of $p_b$ (i.e., extensions of the ocean reanalyses analyzed in Börger et al., 2023). In the absence of such frameworks, the ensemble spreads and statistics presented in the present work provide indications of the errors incurred by methods that neglect intrinsic OAM variations. We expect

these error quantifications to be useful for prediction purposes (e.g., Dobslaw and Dill, 2018; Kiani Shahvandi et al., 2022) and for the geophysical interpretation of the residual geodetic excitation of polar motion from intraseasonal to sub-decadal timescales.

## Appendix A: Effective angular momentum functions

To estimate the effect of fluid angular momentum changes on low-frequency polar motion, we use the frequently cited formalism of Gross (2007):

$$ \hat{\chi} = \hat{\chi}^m + \hat{\chi}^v = \frac{1.100\Omega\Delta\hat{I} + 1.608\hat{h}}{\Omega(C - A')} \quad (A1) $$

Here, $\Omega$ (7.292115 × 10⁻⁵ rad s⁻¹) is Earth’s nominal sidereal rotation rate, $C$ (8.0365 × 10³⁷ kg m²) and $A'$ (8.0102 × 10³⁷ kg m²) are the polar and average equatorial principal moments of inertia of the Earth (including the core), $\hat{h} = h_1 + i h_2$ represents the relative angular momentum of the considered fluid, and $\Delta\hat{I} = \Delta I_{13} + i \Delta I_{23}$ is the time-dependent equatorial perturbation to third-column elements in the tensor of inertia (Moritz and Mueller, 1987). The numerical coefficients account for the responses of the mantle, fluid core, and equilibrium ocean to rotational perturbations caused by the initial excitation and also by the deformation of the solid Earth due to surface mass loading. Implicit to Eq. (A1) is the assumption of no mechanical coupling between mantle and core during polar motion fluctuations. The validity of this assumption on timescales longer than the Chandler period is somewhat contested (Dickman, 2003), but geophysical arguments (Wahr, 2005) suggest that even complete core-mantle coupling would modulate the transfer coefficients by no more than a few percent. Similarly small errors (1%–2%) are incurred by neglecting in Eq. (A1) the
effect of mantle anelasticity on the load Love number (Wahr, 2005; Gross, 2007).

Writing $\Delta\hat{I}$ and $\hat{h}$ in spherical polar coordinates, we have, for the ocean (cf. Gross et al., 2003),

$$ \Delta\hat{I} = -\frac{a^4}{g} \int_{-\frac{\pi}{2}}^{\frac{\pi}{2}} \int_{0}^{2\pi} p_b \sin\phi \cos^2\phi \exp(i\lambda) d\lambda d\phi, \quad (A2) $$

$$ \hat{h} = -a^3 \rho_0 \int_{-\frac{\pi}{2}}^{\frac{\pi}{2}} \int_{0}^{2\pi} (\bar{u}H \sin\phi + i\bar{v}H) \cos\phi \exp(i\lambda) d\lambda d\phi, \quad (A3) $$

where a is Earth's mean radius, $\lambda$ and $\phi$ are the longitude and latitude of an arbitrary point on the Earth, g is the constant gravitational acceleration, $\rho_0$ is a mean seawater density (1025 kg m$^{-3}$), and H denotes the local water depth. The quantities ($p_b, \bar{u}, \bar{v}$) are the dynamical variables specified in Sect. 2.2, representing either forced or intrinsic variations or their sum. Expressions analogous to Eq. (A2) exist for the mass term of other Earth system components, with pressures readily computed from equivalent water heights via the hydrostatic equation. Apart from the ocean, only the atmosphere has a non-zero motion term, usually evaluated as the volume integral over horizontal currents and densities in three-dimensional space (Barnes et al., 1983).

## Appendix B: Terrestrial hydrology and cryosphere

The monthly gravity fields used for the determination of $\hat{\chi}^H$ and $\hat{\chi}^C$ over the 1995–2015 period are computed from tracking observations to up to 16 satellites from both satellite laser ranging (SLR) and the Doppler Orbitography and Radiopositioning Integrated by Satellite (DORIS) system. The problem is treated in terms of variational equations and cast as an iterative fit of forward-modeled satellite orbits to the observations. In representing the gravity field, low-degree spherical harmonics are combined with EOFs derived by principal component analysis of GRACE/-FO gravity fields. This hybrid representation, developed in Löcher and Kusche (2021), allows the same effective wavelengths as the GRACE/-FO solutions while reducing the number of parameters to an amount appropriate for techniques less sensitive to gravity field details. As the approach leaves some scope for the choice of the base functions, the monthly fields in this study are computed as weighted means of 15 solutions employing various numbers of EOFs (up to 14) and spherical harmonics (from degrees 2 to 5).

The coordinates of the SLR and DORIS stations entering the solution are part of the International Reference Frame 2014 (ITRF2014; Altamimi et al., 2016), corrected for tidal and non-tidal loading as appropriate. The force modeling for all satellites accounts for pole tides as per conventions (Pettit and Luzum, 2010) and ocean pole tides following Desai (2002). Short-term mass variations in atmosphere and ocean

are removed via the GRACE dealiasing product AOD1B RL06 (the GAC coefficients; Dobslaw et al., 2017) except for degree-1 terms. Furthermore, the force model includes the gravitational effect of glacial isostatic adjustment (GIA) based on the model of A et al. (2013); the choice of another GIA estimate has negligible effects on our polar motion results. For consistency, we also remove the A et al. (2013) model from the GRACE/-FO solutions prior to expanding them into EOFs. Neither GIA nor the AOD1B product is restored to the final gravity fields. Co- and post-seismic deformation signals are left uncorrected.

In a subsequent step, we also deduce degree-1 coefficients (i.e., geocenter displacements) by keeping the gravity field invariant and estimating station coordinates for all SLR sites using a no-net-rotation condition. The common mode of translation across the derived stations reveals the time-variable offset of the geocenter relative to the ITRF2014. Over the GRACE/-FO time span, our degree-1 time series is very similar to, but somewhat noisier than, the widely used estimate of Swenson et al. (2008). Prior to assigning this geocenter model to the monthly solutions, we subtract degree-1 contributions from the atmosphere (based on ERA-Interim with IB) and the ocean (based on forced OCCIPUT $p_b$). This correction avoids the double-counting of long-wavelength mass variations in our analysis, associated in particular with atmospheric pressure changes over land.

The recovered gravity fields inherit the directional error structures (stripes) from the GRACE/-FO EOFs and thus require filtering in post-processing. We use a standard decorrelation and anisotropic smoothing kernel (DDK2; Kusche et al., 2009) to that end. Each monthly set of Stokes coefficients (including degree-1 terms) is then mapped to gridded (1° × 1°) surface mass anomalies and iteratively corrected for signal leakage across continental boundaries into the ocean. Our approach closely follows Chen et al. (2015). In short, the algorithm (a) takes a first guess of the true gridded mass distribution, (b) derives a “predicted” mass field from item (a) by expanding and re-synthesizing it using spherical harmonics up to degree 60, and (c) successively corrects the prediction – and thus the current guess for the true mass distribution – to yield closer agreement with the particular month’s observed mass anomalies. The latter also defines the first guess but with values over the ocean replaced by their spatial mean. This delineation leads to a sharp transition between the ocean and land (or ice sheets, respectively) that is kept for all guesses in the iteration. We separately iterate land, Greenland, and Antarctica, always permitting adjustments to local mass values within 600 km from the coast and correcting the spatial mean over the ocean accordingly to conserve total mass. Acceptable convergence is typically achieved after five iterations per modeled region. The resulting mass change fields, used as inputs for Eq. (A2), are illustrated in terms of their standard deviation in Figs. B1a–c.
<figure><img src="image_10.png" /></figure>

**Figure B1.** Standard deviation of (a–c) interannual surface mass changes over continents, as deduced from long-term gravity field solutions, and (d–f) associated GAL-induced mass redistribution in the ocean, from 1995 to 2015. Mass and GAL signals are depicted separately for (a, d) land hydrology, (b, e) Greenland, and (c, f) Antarctica. Units are centimeter of equivalent water height. Trends, seasonal oscillations, and sub-annual signals have been removed from each grid point's time series and therefore do not contribute to the figure. Note that the large (~20 cm) localized mass anomaly on the northern Malay Peninsula in panel (a) is a tectonic rather than a hydrological signal; the corresponding GAL field in the area (panel b) should therefore be treated with caution.

**Appendix C: Gravitational attraction and loading**

To quantify the effect of surface mass loads on polar motion excitation, we combine mass load estimates from different sources, generally separated by Earth system compartment. For a precise net excitation estimate, allowance should be made for transports of water between continents and the ocean, leading to the residual excitation $\hat{\chi}^L$. We require our estimate of $\hat{\chi}^L$ to be consistent with the physics of GAL, such that the total load mass is conserved and the modeled sea level is an equipotential of Earth's gravity field (Nakioglu and Lambeck, 1980; Quinn et al., 2015; Adhikari and Ivins, 2016). Under the assumption that these mass fluxes and gravitational effects cause no relevant dynamic variations in the ocean (e.g., Kuhlmann et al., 2011), the resultant non-uniform water mass displacements may be estimated by solving the sea level equation (Woodward, 1888; Farrell and Clark, 1976). Here, sea level is taken relative to the deformed

crust, equivalent to a change in $p_b$ in the absence of atmospheric pressure variations. We solve the sea level equation in the spectral domain based on a global fingerprint inversion framework (Rietbroek et al., 2016; Uebbing et al., 2019) that accounts for solid Earth deformation, perturbations in the gravitational potential, sea level changes due to rotational feedback from the GAL-induced mass variations, and global mass conservation. The Earth's deformation is assumed to be elastic, as is customary for sub-decadal phenomena.

We apply four different loads, comprising (a) atmospheric surface pressure variations, mass changes over (b) Greenland and (c) Antarctica, and (d) TWS anomalies over the remaining continents. The atmospheric loads are monthly mean pressures from ERA-Interim, with values over the ocean replaced by their time-varying spatial mean (cf. Quinn et al., 2015). For terms (b)–(d), we use the satellite-based and leakage-corrected surface mass estimates described in Appendix B. In adopting these loads in the sea level equation, we assume that all continental mass changes are related to
variations in water (or ice) mass and consequently added to or removed from the ocean. Figure B1 illustrates the derived GAL fluctuations in space (bar the atmospheric contribution), while the total excitation signal ($\hat{\chi}^L$) is shown in Fig. 1. On the interannual timescales considered here, the atmospheric contribution to $\hat{\chi}^L$ is very small, featuring a standard deviation of less than 0.20 mas.

A final, more subtle issue is that GAL effects may also arise from dynamic $p_b$ variability, as driven by atmospheric forcing or intrinsic ocean processes, and not by land-ocean mass transfer (Vinogradova et al., 2015). These gravitational effects were omitted in the OCCIPUT simulations, so we account for them after the fact, in keeping with the ocean's tendency for an equilibrium response to body forces at periods longer than a few months. Restricting the loads to dynamic ocean mass changes simplifies the GAL algorithm to a convolution of $p_b$ anomalies with the proper Green's function (Vinogradova et al., 2015). We have applied a spherical harmonic version of this convolution (Schindelegger et al., 2018) to a few OCCIPUT ensemble members and the $p_b$ ensemble mean. The derived excitation series were found to differ in detail but agree in terms of their interannual peaks (~1 mas in $\chi_1$; ~1.5 mas in $\chi_2$). Thus, for the purpose of this study, we represent the ocean dynamics GAL signal using the $p_b$ ensemble mean only and impose the implied mass redistribution on each of the 50 ensemble members. The GAL effects due to ocean dynamics are therefore contained in $\hat{\chi}^O$ and not in $\hat{\chi}^L$.

**Data availability.** All geophysical angular momentum estimates analyzed in this study have been placed at https://doi.org/10.5281/zenodo.12664036 (Börger et al., 2024).

**Supplement.** The supplement related to this article is available online at: https://doi.org/10.5194/esd-16-75-2025-supplement.

**Author contributions.** LB: data curation, formal analysis, investigation, methodology, software, validation, visualization, and writing (original draft and review and editing). MS: conceptualization, formal analysis, funding acquisition, methodology, supervision, validation, and writing (original draft and review and editing). MZ: data curation, formal analysis, and writing (review and editing). RMP: conceptualization, funding acquisition, and writing (review and editing). AL: formal analysis and writing (review and editing). BU: formal analysis and writing (review and editing). JMM: data curation and writing (review and editing). TP: funding acquisition and writing (review and editing).

**Competing interests.** The contact author has declared that none of the authors has any competing interests.

**Disclaimer.** Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims made in the text, published maps, institutional affiliations, or any other geographical representation in this paper. While Copernicus Publications makes every effort to include appropriate place names, the final responsibility lies with the authors.

**Acknowledgements.** The results of this research have been achieved using the PRACE Research Infrastructure resource CURIE based in France at TGCC. This work is a contribution to the OCCIPUT and IMHOTEP projects. OCCIPUT is funded by ANR through grant no. ANR-13-BS06-0007-01. IMHOTEP is funded by CNES through the Ocean Surface Topography Science Team (OSTST).

**Financial support.** This research has been supported by the Deutsche Forschungsgemeinschaft (grant nos. 45939286 and 459420781), the National Aeronautics and Space Administration (grant nos. 80NSSC20K0728 and 80NSSC24K1154), and the Agence Nationale de la Recherche (grant no. ANR-13-BS06-0007-01).

This open-access publication was funded by the University of Bonn.

**Review statement.** This paper was edited by Claudia Pasquero and reviewed by Christian Bizouard and two anonymous referees.

## References

A, G., Wahr, J., and Zhong, S.: Computations of the Viscoelastic Response of a 3-D Compressible Earth to Surface Loading: An Application to Glacial Isostatic Adjustment in Antarctica and Canada, Geophys. J. Int., 192, 557–572, https://doi.org/10.1093/gji/ggs030, 2013.

Adhikari, S. and Ivins, E. R.: Climate-Driven Polar Motion: 2003–2015, Sci. Adv., 2, e1501693, https://doi.org/10.1126/sciadv.1501693, 2016.

Afroosa, M., Rohith, B., Paul, A., Durand, F., Bourdallé-Badie, R., Sreedevi, P. V., de Viron, O., Ballu, V., and Shenoi, S. S. C.: Madden-Julian Oscillation Winds Excite an Intraseasonal See-Saw of Ocean Mass That Affects Earth's Polar Motion, Communications Earth & Environment, 2, 139, https://doi.org/10.1038/s43247-021-00210-x, 2021.

Altamimi, Z., Rebischung, P., Métivier, L., and Collilieux, X.: ITRF2014: A new release of the International Terrestrial Reference Frame modeling nonlinear station motions, J. Geophys. Res.-Sol. Ea., 121, 6109–6131, https://doi.org/10.1002/2016JB013098, 2016.

Arbic, B. K., Müller, M., Richman, J. G., Shriver, J. F., Morten, A. J., Scott, R. B., Sérazin, G., and Penduff, T.: Geostrophic turbulence in the frequency–wavenumber domain: Eddy-driven low-frequency variability, J. Phys. Oceanogr., 44, 2050–2069, https://doi.org/10.1175/JPO-D-13-054.1, 2014.
Barnes, R., Hide, R., White, A., and Wilson, C.: Atmospheric Angular Momentum Fluctuations, Length-of-Day Changes and Polar Motion, P. Roy. Soc. Lond. A Mat., 387, 31–73, https://doi.org/10.1098/rspa.1983.0050, 1983.

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.: Development of a probabilistic ocean modelling system based on NEMO 3.5: application at eddying resolution, Geosci. Model Dev., 10, 1091–1106, https://doi.org/10.5194/gmd-10-1091-2017, 2017.

Bizouard, C. and Seoane, L.: Atmospheric and Oceanic Forcing of the Rapid Polar Motion, J. Geodesy, 84, 19–30, https://doi.org/10.1007/s00190-009-0341-2, 2010.

Börger, L., Schindelegger, M., Dobslaw, H., and Salstein, D.: Are ocean reanalyses useful for Earth rotation research?, Earth and Space Science, 10, e2022EA002700, https://doi.org/10.1029/2022EA002700, 2023.

Börger, L., Schindelegger, M., Zhao, M., Ponte, R. M., Löcher, A., Uebbing, B., Molines, J.-M., and Penduff, T.: Angular momentum estimates for global geophysical fluids, 1995–2015, Zenodo [data set], https://doi.org/10.5281/zenodo.12664036, 2024.

Brankart, J.-M., Candille, G., Garnier, F., Calone, C., Melet, A., Bouttier, P.-A., Brasseur, P., and Verron, J.: A generic approach to explicit simulation of uncertainty in the NEMO ocean model, Geosci. Model Dev., 8, 1285–1297, https://doi.org/10.5194/gmd-8-1285-2015, 2015.

Brzeziński, A.: Polar motion excitation by variations of the effective angular momentum function: Considerations concerning deconvolution problem, Manuscr. Geodaet., 17, 3–20, 1992.

Carret, A., Llovel, W., Penduff, T., and Molines, J.-M.: Atmospherically Forced and Chaotic Interannual Variability of Regional Sea Level and Its Components Over 1993–2015, J. Geophys. Res.-Oceans, 126, e2020JC017123, https://doi.org/10.1029/2020JC017123, 2021.

Chelton, D. B., Schlax, M. G., and Samelson, R. M.: Global observations of nonlinear mesoscale eddies, Prog. Oceanogr., 91, 167–216, https://doi.org/10.1016/j.poccean.2011.01.002, 2011.

Chen, J., Wilson, C. R., Kuang, W., and Chao, B. F.: Interannual Oscillations in Earth Rotation, J. Geophys. Res.-Sol. Ea., 124, 13404–13414, https://doi.org/10.1029/2019JB018541, 2019.

Chen, J. L., Wilson, C. R., Ries, J. C., and Tapley, B. D.: Rapid Ice Melting Drives Earth’s Pole to the East, Geophys. Res. Lett., 40, 2625–2630, https://doi.org/10.1002/grl.50552, 2013.

Chen, J. L., Wilson, C. R., Li, J., and Zhang, Z.: Reducing Leakage Error in GRACE-observed Long-Term Ice Mass Change: A Case Study in West Antarctica, J. Geodesy, 89, 925–940, https://doi.org/10.1007/s00190-015-0824-2, 2015.

Clarke, P. J., Lavallée, D. A., Blewitt, G., van Dam, T. M., and Wahr, J. M.: Effect of gravitational consistency and mass conservation on seasonal surface mass loading models, Geophys. Res. Lett., 32, L08306, https://doi.org/10.1029/2005GL022441, 2005.

Cravatte, S., Serazin, G., Penduff, T., and Menkes, C.: Imprint of chaotic ocean variability on transports in the southwestern Pacific at interannual timescales, Ocean Sci., 17, 487–507, https://doi.org/10.5194/os-17-487-2021, 2021.

Dee, D. P., Uppala, S. M., Simmons, A. J., Berrisford, P., Poli, P., Kobayashi, S., Andrae, U., Balmaseda, M. A., Balsamo, G., Bauer, P., Bechtold, P., Beljaars, A. C. M., van de Berg, L., Bidlot, J., Bormann, N., Delsol, C., Dragani, R., Fuentes, M., Geer,

A. J., Haimberger, L., Healy, S. B., Hersbach, H., Hólm, E. V., Isaksen, L., Källberg, P., Köhler, M., Matricardi, M., McNally, A. P., Monge-Sanz, B. M., Morcrette, J.-J., Park, B.-K., Peubey, C., de Rosnay, P., Tavolato, C., Thépaut, J.-N., and Vitart, F.: The ERA-Interim reanalysis: Configuration and performance of the data assimilation system, Q. J. Roy. Meteor. Soc., 137, 553–597, https://doi.org/10.1002/qj.828, 2011.

Deng, S., Liu, S., Mo, X., Jiang, L., and Bauer-Gottwein, P.: Polar Drift in the 1990s Explained by Terrestrial Water Storage Changes, Geophys. Res. Lett., 48, e2020GL092114, https://doi.org/10.1029/2020GL092114, 2021.

Desai, S. D.: Observing the pole tide with satellite altimetry, J. Geophys. Res.-Oceans, 107, 3186, https://doi.org/10.1029/2001JC001224, 2002.

Dickman, S. R.: Evaluation of “effective angular momentum function” formulations with respect to core-mantle coupling, J. Geophys. Res.-Sol. Ea., 108, 2150, https://doi.org/10.1029/2001JB001603, 2003.

Ding, H., Pan, Y., Xu, X. Y., Shen, W., and Li, M.: Application of the AR-z spectrum to polar motion: A possible first detection of the Inner Core Wobble and its implications for the density of Earth’s core, Geophys. Res. Lett., 46, 13765–13774, https://doi.org/10.1029/2019GL085268, 2019.

Dobslaw, H. and Dill, R.: Predicting Earth orientation changes from global forecasts of atmosphere-hydrosphere dynamics, Adv. Space Res., 61, 1047–1054, https://doi.org/10.1016/j.asr.2017.11.044, 2018.

Dobslaw, H., Bergmann-Wolf, I., Dill, R., Poropat, L., and Flechtnner, F.: Product Description Document for AOD1B Release 06, Rev. 6.1, GFZ Potsdam, Potsdam, Germany, 2017.

Dussin, R., Barnier, B., Brodeau, L., and Molines, J.-M.: Drakkar Forcing Set DFS5, Tech. rep., DRAKKAR/MyOcean Report 01-04-16, 2016.

Farrell, W. E. and Clark, J. A.: On Postglacial Sea Level, Geophys. J. Int., 46, 647–667, https://doi.org/10.1111/j.1365-246X.1976.tb01252.x, 1976.

Göttl, F., Groh, A., Schmidt, M., Schröder, L., and Seitz, F.: The Influence of Antarctic Ice Loss on Polar Motion: An Assessment Based on GRACE and Multi-Mission Satellite Altimetry, Earth Planets Space, 73, 99, https://doi.org/10.1186/s40623-021-01403-6, 2021.

Gross, R. S.: Correspondence between theory and observations of polar motion, Geophys. J. Int., 109, 162–170, https://doi.org/10.1111/j.1365-246X.1992.tb00086.x, 1992.

Gross, R. S.: The Excitation of the Chandler Wobble, Geophys. Res. Lett., 27, 2329–2332, https://doi.org/10.1029/2000GL011450, 2000.

Gross, R. S.: Earth Rotation Variations – Long Period, in: Treatise on Geophysics, edited by: Herring, T., Elsevier, Amsterdam, 3, 239–294, 2007.

Gross, R. S., Fukumori, I., and Menemenlis, D.: Atmospheric and Oceanic Excitation of the Earth’s Wobbles during 1980–2000, J. Geophys. Res.-Sol. Ea., 108, 2370, https://doi.org/10.1029/2002JB002143, 2003.

Harker, A. A., Schindelegger, M., Ponte, R. M., and Salstein, D. A.: Modeling ocean-induced rapid Earth rotation variations: An update, J. Geodesy, 95, 110, https://doi.org/10.1007/s00190-021-01555-z, 2021.
Hogg, A. M., Penduff, T., Close, S. E., Dewar, W. K., Constantinou, N. C., and Martínez-Moreno, J.: Circumpolar variations in the chaotic nature of Southern Ocean eddy dynamics, J. Geophys. Res.-Oceans, 127, e2022JC018440, https://doi.org/10.1029/2022JC018440, 2022.

Kiani Shahvandi, M., Schartner, M., and Soja, B.: Neural ODE differential learning and its application in polar motion prediction, J. Geophys. Res.-Sol. Ea., 127, e2022JB024775, https://doi.org/10.1029/2022JB024775, 2022.

Kuang, W., Chao, B. F., and Chen, J.: Reassessment of electromagnetic core-mantle coupling and its implications to the Earth's decadal polar motion, Geodesy and Geodynamics, 10, 356–362, https://doi.org/10.1016/j.geog.2019.06.003, 2019.

Kuhlmann, J., Dobslaw, H., and Thomas, M.: Improved modeling of sea level patterns by incorporating self-attraction and loading, J. Geophys. Res.-Oceans, 116, C11036, https://doi.org/10.1029/2011JC007399, 2011.

Kusche, J., Schmidt, R., Petrovic, S., and Rietbroek, R.: Decorrelated GRACE Time-Variable Gravity Solutions by GFZ, and Their Validation Using a Hydrological Model, J. Geodesy, 83, 903–913, https://doi.org/10.1007/s00190-009-0308-3, 2009.

Leroux, S., Penduff, T., Bessières, L., Molines, J.-M., Brankart, J.-M., Sérazin, G., Barnier, B., and Terray, L.: Intrinsic and Atmospherically Forced Variability of the AMOC: Insights from a Large-Ensemble Ocean Hindcast, J. Climate, 31, 1183–1203, https://doi.org/10.1175/JCLI-D-17-0168.1, 2018.

Llovel, W., Penduff, T., Meyssignac, B., Molines, J.-M., Terray, L., Bessières, L., and Barnier, B.: Contributions of Atmospheric Forcing and Chaotic Ocean Variability to Regional Sea Level Trends Over 1993–2015, Geophys. Res. Lett., 45, 13405–13413, https://doi.org/10.1029/2018GL080838, 2018.

Löcher, A. and Kusche, J.: A Hybrid Approach for Recovering High-Resolution Temporal Gravity Fields from Satellite Laser Ranging, J. Geodesy, 95, 6, https://doi.org/10.1007/s00190-020-01460-x, 2021.

Maher, N., Milinski, S., and Ludwig, R.: Large ensemble climate model simulations: introduction, overview, and future prospects for utilising multiple types of large ensemble, Earth Syst. Dynam., 12, 401–418, https://doi.org/10.5194/esd-12-401-2021, 2021.

Meyrath, T. and Van Dam, T.: A Comparison of Interannual Hydrological Polar Motion Excitation from GRACE and Geodetic Observations, J. Geodyn., 99, 1–9, https://doi.org/10.1016/j.jog.2016.03.011, 2016.

Moritz, H. and Mueller, I. I.: Earth Rotation: Theory and Observation, Ungar, New York, ISBN 978-0-8044-4671-6, 1987.

Munk, W. H. and MacDonald, G. J. F.: The Rotation of the Earth: A Geophysical Discussion, Cambridge University Press, New York, https://doi.org/10.1017/S0016756800060726, 1960.

Nakiboglu, S. M. and Lambeck, K.: Deglaciation effects on the rotation of the Earth, Geophys. J. Int., 62, 49–58, https://doi.org/10.1111/j.1365-246X.1980.tb04843.x, 1980.

Nastula, J. and Salstein, D.: Regional atmospheric angular momentum contributions to polar motion excitation, J. Geophys. Res.-Sol. Ea., 104, 7347–7358, https://doi.org/10.1029/1998JB900077, 1999.

Nastula, J., Wińska, M., Śliwińska, J., and Salstein, D.: Hydrological Signals in Polar Motion Excitation – Evidence after Fifteen Years of the GRACE Mission, J. Geodyn., 124, 119–132, https://doi.org/10.1016/j.jog.2019.01.014, 2019.

O'Rourke, A. K., Arbic, B. K., and Griffies, S. M.: Frequency-domain analysis of atmospherically forced versus intrinsic ocean surface kinetic energy variability in GFDL's CM2-O model hierarchy, J. Climate, 31, 1789–1810, https://doi.org/10.1175/JCLI-D-17-0024.1, 2018.

Penduff, T., Bernier, B., Terray, L., Bessières, L., Sérazin, G., Gregorio, S., Brankart, J.-M., Moine, M.-P., Brankart, J.-M., and Brasseur, P.: Ensembles of Eddying Ocean Simulations for Climate, CLIVAR Exchanges No. 65, 19, 26–29, 2014.

Petit, G. and Luzum, B.: IERS Conventions 2010, IERS Technical Note No. 36, Verlag des Bundesamtes für Kartographie und Geodäsie, Frankfurt am Main, ISSN: 1019-4568 (print version), 2010.

Ponte, R. M. and Piecuch, C. G.: Interannual bottom pressure signals in the Australian–Antarctic and Bellingshausen Basins, J. Phys. Oceanogr., 44, 1456–1465, https://doi.org/10.1175/JPO-D-13-0223.1, 2014.

Ponte, R. M., Stammer, D., and Marshall, J.: Oceanic Signals in Observed Motions of the Earth's Pole of Rotation, Nature, 391, 476–479, https://doi.org/10.1038/35126, 1998.

Provost, C., Renault, A., Barré, N., Sennéchael, N., Garçon, V., Sudre, J., and Huhn, O.: Two repeat crossings of Drake Passage in austral summer 2006: Short-term variations and evidence for considerable ventilation of intermediate and deep waters, Deep-Sea Res. Pt. II, 58, 2555–2571, https://doi.org/10.1016/j.dsr2.2011.06.009, 2011.

Quinn, K. J., Ponte, R. M., and Tamisiea, M. E.: Impact of self-attraction and loading on Earth rotation, J. Geophys. Res.-Sol. Ea., 120, 4510–4521, https://doi.org/10.1002/2015JB011980, 2015.

Ratcliff, T. and Gross, R. S.: Combinations of Earth Orientation Measurements: SPACE2018, COMB2018, and POLE2018, Jet Propulsion Laboratory, California Institute of Technology, Publication 19-7, https://hdl.handle.net/2014/46964 (last access: 18 December 2024), 2019.

Rekier, J., Chao, B. F., Chen, J., Dehant, V., Rosat, S., and Zhu, P.: Earth's rotation: Observations and relation to deep interior, Surv. Geophys., 43, 149–175, https://doi.org/10.1007/s10712-021-09669-x, 2022.

Rhines, P. B.: Mesoscale Eddies, in: Encyclopedia of Ocean Sciences, 3rd edn., edited by: Cochran, J. K., Bokuniewicz, H. J., and Yager, P. L., Academic Press, Oxford, 115–127, https://doi.org/10.1016/B978-0-12-409548-9.11642-2, 2019.

Rietbroek, R., Brunnabend, S.-E., Kusche, J., Schröter, J., and Dahle, C.: Revisiting the Contemporary Sea-Level Budget on Global and Regional Scales, P. Natl. Acad. Sci. USA, 113, 1504–1509, https://doi.org/10.1073/pnas.1519132113, 2016.

Rochester, M. and Crossley, D.: Earth's long-period wobbles: A Lagrangean description of the Liouville equations, Geophys. J. Int., 176, 40–62, https://doi.org/10.1111/j.1365-246X.2008.03991.x, 2009.

Rosat, S., Gillet, N., Boy, J.-P., Couhert, A., and Dumberry, M.: Interannual variations of degree 2 from geodetic observations and surface processes, Geophys. J. Int., 225, 200–221, 2021.

Schindelegger, M., Salstein, D., and Böhm, J.: Recent Estimates of Earth-atmosphere Interaction Torques and Their Use in Studying
Polar Motion Variability, J. Geophys. Res.-Sol. Ea., 118, 4586-4598, https://doi.org/10.1002/jgrb.50322, 2013.

Schindelegger, M., Green, J. A. M., Wilmes, S.-B., and Haigh, I. D.: Can We Model the Effect of Observed Sea Level Rise on Tides?, J. Geophys. Res.-Oceans, 123, 4593-4609, https://doi.org/10.1029/2018JC013959, 2018.

Sérazin, G., Penduff, T., Grégorio, S., Barnier, B., Molines, J.-M., and Terray, L.: Intrinsic variability of sea level from global ocean simulations: Spatiotemporal scales, J. Climate, 28, 4279-4292, https://doi.org/10.1175/JCLI-D-14-00554.1, 2015.

Sérazin, G., Penduff, T., Barnier, B., Molines, J.-M., Arbic, B. K., Müller, M., and Terray, L.: Inverse cascades of kinetic energy as a source of intrinsic variability: A global OGCM study, J. Phys. Oceanogr., 48, 1385-1408, https://doi.org/10.1175/JPO-D-17-0136.1, 2018.

Swenson, S., Chambers, D., and Wahr, J.: Estimating geocenter variations from a combination of GRACE and ocean model output, J. Geophys. Res.-Sol. Ea., 113, B08410, https://doi.org/10.1029/2007JB005338, 2008.

Tapley, B. D., Watkins, M. M., Flechtner, F., Reigber, C., Bettadpur, S., Rodell, M., Sasgen, I., Famiglietti, J. S., Landerer, F. W., Chambers, D. P., Reager, J. T., Gardner, A. S., Save, H., Ivins, Erik R., Swenson, S. C., Boening, C., Dahle, C., Wiese, D. N., Dobslaw, H., Tamisiea, M. E., and Velicogna, I.: Contributions of GRACE to understanding climate change, Nat. Clim. Change, 9, 358-369, https://doi.org/10.1038/s41558-019-0456-2, 2019.

Uebbing, B., Kusche, J., Rietbroek, R., and Landerer, F. W.: Processing Choices Affect Ocean Mass Estimates From GRACE, J. Geophys. Res.-Oceans, 124, 1029-1044, https://doi.org/10.1029/2018JC014341, 2019.

Venaille, A., Vallis, G. K., and Smith, K. S.: Baroclinic turbulence in the ocean: Analysis with primitive equation and quasi-geostrophic simulations, J. Phys. Oceanogr., 41, 1605-1623, https://doi.org/10.1175/JPO-D-10-05021.1, 2011.

Vinogradova, N. T., Ponte, R. M., Quinn, K. J., Tamisiea, M. E., Campin, J.-M., and Davis, J. L.: Dynamic adjustment of the ocean circulation to self-attraction and loading effects, J. Phys. Oceanogr., 45, 678-689, https://doi.org/10.1175/JPO-D-14-0150.1, 2015.

Wahr, J.: Polar motion models: Angular momentum approach, Forcing of polar motion in the Chandler frequency band: A contribution to understanding international climate change, Cahiers du Centre European de Geodynamique et de Seismologie, 24, 89-102, 2005.

Wahr, J. M.: The effects of the atmosphere and oceans on the Earth's wobble – I. Theory, Geophys. J. Int., 70, 349-372, https://doi.org/10.1111/j.1365-246X.1982.tb04972.x, 1982.

Welch, P.: The Use of Fast Fourier Transform for the Estimation of Power Spectra: A Method Based on Time Averaging over Short, Modified Periodograms, IEEE T. Audio Electroacoust., 15, 70-73, https://doi.org/10.1109/TAU.1967.1161901, 1967.

Wińska, M.: A comparative study of interannual oscillation models for determining geophysical polar motion excitations, Remote Sens., 14, 147, https://doi.org/10.3390/rs14010147, 2022.

Woodward, R. S.: On the form and position of the sea level with special references to its dependence on superficial masses symmetrically disposed about a normal to the Earth's surface, Tech. rep., Govt. Print. Off., Bulletin 48, 88 pp., https://doi.org/10.3133/b48, 1888.

Zhao, M., Ponte, R. M., Penduff, T., Close, S., Llovel, W., and Molines, J.-M.: Imprints of Ocean Chaotic Intrinsic Variability on Bottom Pressure and Implications for Data and Model Analyses, Geophys. Res. Lett., 48, e2021GL096341, https://doi.org/10.1029/2021GL096341, 2021.

Zhao, M., Ponte, R. M., and Penduff, T.: Global-Scale Random Bottom Pressure Fluctuations from Oceanic Intrinsic Variability, Sci. Adv., 9, eadg0278, https://doi.org/10.1126/sciadv.adg0278, 2023.
<figure><img src="image_11.png" /></figure>

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

Supplement of

# Chaotic oceanic excitation of low-frequency polar motion variability

Lara Börger et al.

Correspondence to: Lara Börger (lboerger@igg.uni-bonn.de)

The copyright of individual parts of the supplement might differ from the article licence.
Comparison of hydrological excitation series

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

**Figure S1.** Hydrological contribution $\hat{\chi}^H$ to polar motion excitation $\hat{\chi} = \chi_1 + i\chi_2$ (mas) deduced from the monthly GRACE/SLR/DORIS gravity field solution described in the main text (1995/01–2015/12, blue curves), the COST-G GravIS RL01 continental water storage anomalies (2002/04–2015/08 with gaps, yellow curves, Boergens et al., 2020), and the Global Land Data Assimilation System (GLDAS, 2002/01–2015/12, purple curves, Rodell et al., 2004). The COST-G and GLDAS time series have been adjusted such that their trends and mean values agree with those of GRACE/SLR/DORIS over the respective common time period, starting early 2002. Note that GRACE/SLR/DORIS was detrended over 1995–2015, as in the main text.
<figure><img src="image_14.png" /><figcaption>**Figure S2.** As in Fig. S1 but with each $\chi_1$ and $\chi_2$ time series filtered to periods longer than 14 months and cut back by 4 months at the respective end points.</figcaption></figure>
References

Boergens, E., Dobslaw, H., and Dill, R.: COST-G GravIS RL01 Continental Water Storage Anomalies. V. 0005, GFZ Data Services, https://doi.org/10.5880/COST-G.GRAVIS_01_L3_TWS, 2020.

Rodell, M., Houser, P., Jambor, U., Gottschalck, J., Mitchell, K., et al.: The Global Land Data Assimilation System, B. Am. Meteorol. Soc., 85, 381–394, https://doi.org/10.1175/BAMS-85-3-381, 2004.