Meridional Asymmetry in Recent Decadal Sea‐Level Trends in the Subtropical Pacific Ocean

Recent sea surface height (SSH) trends in the South Pacific are substantially greater than trends in the North Pacific. Here, we use the Estimating the Climate and Circulation of the Ocean Version 4 Release 4 ocean state estimate and the Ocean Reanalysis System 5 to identify the forcing and mechanisms underlying that meridional asymmetry during 2005–2015. Thermosteric contributions dominate the spatial structure in Pacific SSH trends, but contributions from local surface heat fluxes are small. Wind stress trends drive a spin‐up of the South Pacific subtropical gyre and a northward shift of the North Pacific subtropical gyre. A reduced gravity model forced with reanalysis winds qualitatively reproduces the meridional seesaw in sea level, suggesting that asymmetric trends in subtropical wind stress drive a cross‐equatorial heat transport. A reversal in forcing associated with this process could impact near‐term rates of coastal sea‐level change, particularly in Pacific Island communities.

subtropical North Pacific are not discussed in the context of hemispheric imbalances in warming, perhaps due to an assumed association with zonal structure in sea level related to canonical modes of North Pacific decadal climate variability (Hamlington et al., 2012;Han et al., 2017;Merrifield et al., 2012). However, when the dominant modes of Indo-Pacific sea-level variability are subtracted from the altimeter observations, the band of negative sea-level trends in the subtropical North Pacific-and the associated hemispheric asymmetry in the Pacific-remains .
Here we focus on the 2005-2015 period discussed above and investigate the dynamical mechanism underlying hemispheric asymmetry in sea-level and upper-ocean warming trends in the Pacific basin. Specifically, we use a reduced gravity model to quantify the impact of changes in wind-driven subtropical gyre circulation on asymmetric sea-level trends.

Data and Model Products
We use the Version 4 Release 4 data-constrained ocean state estimate from the Estimating the Circulation and Climate of the Ocean consortium (ECCOv4r4). The estimate covers the period 1992-2017. This product updates and extends the Release 1 solution described by Forget et al. (2015), which covers 1992-2011. The ocean model resolution is nominally 1° in the horizontal and there are 50 levels in the vertical surface forcing is derived from ERA-Interim. A particularly useful property of ECCOv4r4 is that it is an exact solution to a freely running general circulation model, and observational constraints are imposed in such a way that conservation laws are obeyed and physical consistency is preserved, allowing for meaningful closed budget diagnostics (e.g., Piecuch et al., 2019).
We use Archiving, Validation and Interpretation of Satellite Oceanographic data (AVISO) Level 4 Absolute Dynamic Topography. These multimission altimeter monthly sea-level anomalies are provided on a regular global grid with 0.25° horizontal resolution starting from year 1993.
We use the European Centre for Medium-Range Weather Forecasts (ECMWF) Ocean Reanalysis System 5 (ORAS5) for 1979-2017 (Zuo et al., 2017). The ocean model resolution is 0.25° in the horizontal and there are 75 levels in the vertical. Results presented in this manuscript are based on ensemble member no. 0.

Reduced Gravity Model
Reduced gravity models (RGMs) have been shown to be useful and able to simulate the first-order, largescale response to wind forcing in many studies (e.g., Long et al., 2020;Qiu et al., 1997). Here, we consider solutions to equations where f is the Coriolis parameter; U, V are zonal and meridional transports per unit width; g′ = 0.02 ms −2 is the reduced gravity; τ = (τ x , τ y ) is the wind stress, and h(x, y, t) is the spatially and temporally varying layer thickness. The constant background layer thickness H = 330 m is chosen such that the baroclinic gravity wave speed    c g H 2.5 ms −1 , generally consistent with observations (Chelton & Schlax, 1996). Frictional terms in brackets are retained formally in the limit ν → 0 s −1 , thereby allowing for infinitely narrow western boundary layers to close the circulation. Sea surface height (SSH) η is linearly related to layer depth in the reduced gravity model (RGM), with ρ 1 = 1025.9 kg m −3 and ρ 2 = 1028 kg m −3 being the density in the surface and deep layer, respectively, and η o some constant reference level.
A condition of no flow is imposed on all lateral boundaries of a Pacific Ocean basin extending meridionally from 40°S to 60°N. The eastern (x e ) and western boundary (x w ) is set as in ORAS5. Note that Equations 1-3 assume that adjustments associated with clockwise and counterclockwise propagation of boundary waves in the Southern and NH respectively (e.g., Kawase, 1987;Johnson & Marshall, 2002) occur instantaneously. This assumption is reasonable: For example, the transit time for a boundary wave propagating at gravity wave speed from Sydney to San Diego (propagating first equatorward along the western boundary, then crossing the Pacific at the equator and propagating poleward along the eastern boundary) is a few months, which is fast relative to the decadal timescales of interest here.
Cross differentiation of the inviscid (1) and (2) and substitution into (3) yields which governs the interior ocean solution. Solutions are obtained by integrating westward from the eastern boundary along Rossby wave characteristics, i.e., h t = −c r h x − w ek , with c r = −βc 2 f −2 being the Rossby wave speed and x y w f f the Ekman pumping velocity. We tested the impact of altering the Rossby wave speed c r to match more closely the observed wave speed (Chelton & Schlax, 1996) and found that reasonable alternative selections for the Rossby wave speed did not qualitatively change the outcome or conclusions reported here. In the equatorial region (5°S < y < 5°N), where f is small and Rossby waves cross the Pacific ocean basin within a few months (Chelton & Schlax, 1996), we neglect the h t -term in Equation 5, thereby obtaining an equilibrium solution at every time step.
At all latitudes, solutions require an eastern boundary condition h e = h(x e , y, t) to be determined. Since the boundary wave response is infinitely fast in our model, h e is directly obtained from the preservation of domain-integrated surface-layer volume in time, M o = ∫h dA, which follows from continuity and the noflow boundary condition (In practice, h e is iterated at every time step such that the solution satisfies mass conservation). We furthermore neglect the effect of alongshore winds, such that (2) implies h e (t) is constant along the eastern boundary.
The model is forced by ERA-Interim kinematic wind stress, which is also the basis of the boundary conditions for ECCOv4r4 and ORAS5, and RGM solutions are integrated from 1979 to 2017. We also tested the impact of using monthly NCEP/NCAR reanalysis wind stress (years 1948-2016) to confirm the robustness of our results. Initially, h is set to the equilibrium response to the wind stress averaged over the first 15 years. We tested the impact of using different initial conditions (e.g., average over first 5 years), and the impact is negligible for the period of primary interest (2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015), presumably because the time of integration before the start of the study period is sufficiently long. Since there is no explicit diffusion in our model, the memory of earlier states is erased primarily as anomalies propagate into the western boundary, where dynamics average over all incoming anomalies.
Meridional heat transports across a latitude are computed from the RGM, and compared to meridional heat transports from ECCOv4r4. Temperature is uniform throughout the RGM surface layer, hence heat transports are proportional to the zonally integrated volume transports in the surface layer, which are (due to the closed boundaries) compensated by opposite transports at depth associated with colder temperatures.
Since solutions depend only on density stratification (g′) and not directly on temperature, we refrain from defining layer temperatures and scale heat transports such that they have similar amplitudes as those in ECCOv4r4. Specifically, upper-layer volume fluxes are multiplied with a factor 0.45 × 10 −7 PW/(m 3 /s) to calculate heat fluxes, which with a specific heat of 3994 J/(kg C) and density of 1030 kg/m 3 corresponds to a temperature difference of about 12 °C between layers one and two. Figure 1a shows the first empirical orthogonal functions (EOFs) of 10 year averages of ECCOv4r4 Pacific sea-level and thermosteric sea-level anomalies, calculated relative to the Pacific Ocean mean, and then averaged over boxes extending 10° in latitude and zonally across the Pacific Ocean basin. Amplitudes are largest, and of opposite sign, in the subtropical/mid-latitude regions in both hemispheres, indicated by gray shading (For brevity, we refer to both of these regions as subtropical in the discussion below, even though the center of gravity occurs at somewhat higher latitudes in the SH.) The corresponding principle components in Figure 1b indicate a swing from anomalously high sea levels in the NH to anomalously low sea levels beginning around 2005, and an opposite swing in the SH (The same pattern is revealed by decadal sea-level anomalies shown in Figure S2.)

Decadal Averages and Trends
These sea-level changes also manifest themselves in trends ( Figure 1c) during 2005-2015, which was previously identified as a period of hemispheric imbalance in ocean warming (Rathore et al., 2020). Specifically, trends are positive in the SH and negative in the NH subtropics, with an average amplitude of 2-3 mm yr −1 , which is about the same magnitude as the coeval global-mean trend. These trends are consistent among observations (AVISO) and model products (ECCOv4r4 and ORAS5).
The similarity between dynamic and thermosteric sea-level anomalies (Figure 1a and 1b) suggests that these meridionally asymmetric sea-level trends are associated with changes in ocean heat content. The ECCOv4r4 solution provides runtime-averaged horizontal and surface heat fluxes, which allows for offline diagnosis of a closed heat budget (e.g., Piecuch et al., 2019). As expected given the dominance of thermosteric sea-level changes, and consistent with the findings by Rathore et al. (2020), Pacific heat content trends in Figure 1d show a meridional asymmetry similar to that in sea level, with increased accumulation of heat in the SH. Furthermore, the changes in heat content are dominated by the convergence of horizontal heat transports, suggesting meridional redistribution of heat from NH to SH latitudes. The contribution of surface heat fluxes, on the other hand, is generally smaller in magnitude and of opposite sign to the overall heat content trends, thereby suggesting a damping process. Figure 2 shows global maps for sea-level trends relative to the Pacific mean for ECCOv4r4, AVISO, and ORAS5. The maps indicate that the different products agree well in terms of the large-scale trend patterns. Differences occur on smaller scales associated with mesoscale eddies, which is unsurprising, since mesoscale features are either unresolved or only partially resolved by these models, and so we do not necessarily expect consistency between data and model products on these spatial scales.

Forcing
Since the meridional asymmetry in sea-level trends is mostly associated with meridional redistribution of heat (Figure 1), wind-forcing is a reasonable candidate for the underlying cause. Figure 3a shows trends in ERA-Interim sea-level pressure and wind stress during 2004-2014. This trend period is chosen assuming a 1-year lag in the sea-level response to wind forcing via baroclinic adjustment (This assumption is based on lag-correlation coefficients computed between sea-level anomalies over 10°-30°N from ECCOv4r4, and a steady state version of the RGM driven by annual-mean winds, which peaks at 1-year lag; see Figure S3. Note that large-scale trends in sea-level pressure and wind stress are not strongly sensitive to precise boundaries of that period.) The most coherent change in the Pacific is a strengthening of the South Pacific High, which couples to intensification of SH trade winds and westerlies. A comparatively small reduction in trade-wind strength occurs over the NH latitude range 10°-30°N, which is surprising given the reduction in ocean heat content and negative sea-level trends relative to the Pacific mean. Figure 3b shows the streamfunction for the mean Sverdrup gyres in the Pacific over the full record length of ERA-Interim (1979contours) and trends in Sverdrup transport streamfunction (shading) based on ERA-Interim wind-stress trends during 2004-2014. In the SH, the strengthening of the large-scale winds leads to a spin-up of the gyre circulation. In contrast, the NH subtropical gyre weakens south of 30°N and strengthens north of this latitude, indicating a northward shift of the gyre circulation. As the strength and position of the subtropical gyres strongly impact the sea-level topography; these trends in gyre circulation are generally consistent with the asymmetric sea-level and heat-content trends in Figure 1. Although showing some regional differences, the general trends in wind forcing and Sverdrup transports are similar for NCEP ( Figure S4), indicating robustness of the results. Next, we utilize the RGM (Section 2.2) to better quantify their impact. trends, in particular the hemispheric asymmetry resulting from negative trends in the NH subtropics and strong positive trends in the western parts of the subtropical and midlatitude South Pacific. The meridional asymmetry is reproduced for RGM solutions driven by ERA-Interim as well as NCEP winds (Figure 1c). The greatest differences in spatial structure of the trends from the RGM compared to the observations and reanalysis are along the equator (Figure 2). Due to its simplicity, the RGM is not expected to reproduce observations in the equatorial region due to not resolving equatorial waves, the relatively stronger impact of interannual variability, surface heat fluxes playing a leading role in decadal SSH variability (e.g., Piecuch et al., 2019), and the effect of the Indonesian Throughflow (see Section 3.4). Discrepancies in the equatorial region do not, however, reduce the utility of the RGM for understanding the basic processes leading to hemispheric asymmetry in sea-level trends away from the equator.  Figures 4a and 4b).

Reduced Gravity Model
To isolate the impact of SH winds on these trends, we conducted a RGM sensitivity experiment with the Ekman pumping south of 10°S being replaced by its climatological value. We refer to this sensitivity experiment as RGM-S10. In RGM-S10, the sea-level trend during 2005-2015 in the NH subtropics almost vanishes. This supports the hypothesis that the mechanism most responsible for the asymmetry in sea-level trends during 2005-2015 in the RGM is the wind-driven spin-up of the SH subtropical gyre. As warm surface waters converge into the SH subtropical gyre region and deepen the thermocline there, the basin-wide adjustment (i.e., thinning of the thermocline elsewhere) occurs via boundary waves propagating equatorward along the western boundary, then eastward across the basin at the equator, and finally poleward along the eastern boundary. In the absence of strong surface forcing over subtropical latitudes of the NH, Rossby waves radiate into the basin from the eastern boundary, shoaling the thermocline and reducing sea level (Fu & Qiu, 2002) (note that our model does not explicitly resolve SCHLOESSER ET AL. 10.1029/2020GL091959 6 of 10 boundary waves, and these implicit boundary adjustments occur instantaneously, Section 2.2). A similar pattern occurs with NCEP winds (Figures S4 and S5), however, with the relative impact of NH winds on the NH sea-level trends being larger. To further test the hypothesis that changes in gyre circulation drive a hemispheric redistribution of upper-ocean waters, we analyze Pacific meridional heat fluxes.

Meridional Heat Fluxes
Figure 4b and 4c show ECCOv4r4 Pacific annual mean meridional heat flux anomalies (i.e., with the timemean removed) across 10°N and 10°S. While both time series exhibit strong interannual variability, they reveal a tendency for fluxes to be negative relative to their time means (i.e., anomalous redistribution of heat from north to south) during the period 2005-2015. These tendencies are emphasized in the 9 year running means shown in Figure 4d, and are consistent with our previous results as well as those of Rathore et al. (2020). In agreement with the earlier analysis, the RGM captures the decadal trends, and some of the interannual variability (Figure 4b and 4c; 9 year running means are shown in Figure S6; see Section 2.2 for how meridional heat fluxes are calculated in the RGM).
The assumption of closed boundaries in the RGM implies that heat is redistributed within the Pacific on the decadal timescales, which may not hold strictly in reality, since the Pacific Ocean is not a closed basin. The goodness of the closed-boundary assumption used in the RGM is tested by calculating ECCOv4r4 heat fluxes across the RGM domain boundaries (Figure 4d). The flux across 60°N and 40°S are at least an order of magnitude smaller than those across 10°S and 10°N, suggesting that the assumption is reasonable, and SCHLOESSER ET AL.
10.1029/2020GL091959 7 of 10 hemispheric exchange of heat between the subtropical basins in the Pacific is much larger than exchange with other ocean basins on decadal timescales.
Anomalies of Indonesian Throughflow heat fluxes (calculated as the sum of the fluxes through the four open sections between Malaysia, Sumatra, Java, Lombok, and Australia), on the other hand, are of similar magnitudes and-on decadal scales-correlate well with the discrepancy between heat fluxes at 10°S and 10°N. This is consistent with previous findings (Mayer et al., 2018)   Marshall, 2002) suggest that fluxes through the Indonesian Throughflow do not substantially impact the Pacific sea-level trend difference between NH and SH subtropics of primary interest here.

Summary
We have assessed the impact of changes in wind-driven gyre circulation on sea-level trends and redistribution of heat in the Pacific over the period 2005-2015. Our results suggest that this mechanism substantially contributes to previously identified meridionally asymmetric trends in sea level and heat content (Rathore et al., 2020). Specifically, an intensification of the South Pacific High leads to a strengthening of SH winds, and a spin up of the subtropical gyre. In the NH, the subtropical gyre shifts northward, but does not show an increase in strength. This difference leads to a deepening of the SH thermocline relative to the NH, which is associated with an anomalous, interhemispheric heat transport and thermosteric sea-level trend in the SH subtropics at the expense of NH heat content. On the large scale, the sea-level trends associated with this seesaw pattern are of the order of several mm/year, which is comparable to the magnitude of recent global-mean sea-level trends (Beckley et al., 2017;Church & White, 2006).
Finally, while our analysis indicates a coherent shift in the meridional structure of Pacific sea-level trends, dominated by changes in the subtropics, the results do not necessarily suggest that hemispheric asymmetry is the result of an internal ocean or climate mode. The results do, however, provide a mechanism whereby interhemsipheric heat exchange and hemispheric asymmetry in sea-level trends relate to subtropical wind forcing that could result from forced or unforced climate variability. Ultimately, the satellite-altimetry record may be too short to unambiguously distinguish decadal and longer period variability from forced secular trends, and the trends identified here could result from integration of shorter scale variability. Since the Rossby wave speed decreases with latitude, zonal thermocline slopes and depth anomalies tend to increase (Equation 5), and so does the potential for larger sea-level trends (relative to the basin mean) on decadal and longer timescales. In other words, while the ocean generally integrates and acts as a low-pass filter for atmospheric forcing (Hasselmann, 1976), the timescales of the ocean response increase with latitude. This view is consistent with our finding of subtropical dynamics to be increasingly important on decadal timescales, whereas the interhemispheric transport has been shown to be dominated by tropical processes on shorter time scales (McGregor et al., 2014).