Sea level extremes and compounding marine heatwaves in coastal Indonesia

Extreme sea level events are among the most consequential manifestations of climate change, with anthropogenic global sea level rise magnifying flooding and clear-sky flood frequency. While daily-scale extremes driven by storms and tides have been well studied, climate-variability-driven sea level extremes and their evolution under anthropogenic change remain less explored. ENSO exerts global climate impacts and often induces strong marine heatwaves in the Indonesian–Australian basin, as occurred during the 2015–2016 El Niño and subsequent 2016 negative IOD. Sea level Height EXtreme (HEX) events and marine heatwaves each pose substantial ecological, economic, and social risks; combined as compound extremes (CHHEX) they can be even more damaging. The Indian Ocean rim, hosting one-third of the world’s population, includes vulnerable low-lying coastal regions. Indonesia, at the Indo-Pacific warm pool confluence and influenced by monsoons, the IOD, and ENSO, with rapid coastal urbanization and land subsidence, is an ideal testbed for studying HEX and CHHEX events. This study integrates in situ and satellite observations with reanalysis and model experiments to detect HEX and CHHEX along Indonesia’s Indian Ocean coasts since 1993, extend analyses to the 1960s, and elucidate the roles of anthropogenic sea level rise and decadal variability and the physical mechanisms and forcings involved.

Background work highlights: anthropogenic sea level rise has accelerated in recent decades and exacerbates coastal flooding. ENSO is the dominant interannual mode with global impacts, and El Niño events can trigger strong marine heatwaves in the tropical Indian Ocean. The IOD is a key coupled mode in the Indian Ocean, with its negative phase producing warm SST and deeper thermocline in the eastern basin that affects Indonesia. Compound extremes (e.g., heatwaves with drought on land) are becoming more common in a warming climate, yet integrated studies of sea level extremes with coincident marine heatwaves (CHHEX) remain limited. Prior reanalyses and observational studies have characterized storm surges and extreme sea levels, but climate-driven, monthly-scale extremes and their compounding with marine heatwaves along coastal Indonesia, and their links to monsoons, ENSO, and IOD under anthropogenic change, have received comparatively less attention.

Data and detection: Monthly multi-mission satellite altimetry on 0.25° grids (1993–2018), tide gauge data at Cilacap B (2007–2016, GIA and IB corrected), CCMP surface winds, NOAA OISST monthly and daily SST, and ECMWF ORAS4 ocean reanalysis (1958–2017; SLA, subsurface T used to estimate thermocline depth via D20). HEX events are defined as monthly mean SLAs exceeding the 90th percentile relative to 1993–2018. Marine heatwaves (MHWs) are defined as monthly SST anomalies (relative to 1989–2018 climatology) exceeding the 90th percentile. CHHEX events occur when MHWs co-occur with HEXs (allowing a short lead for December 2010 as described). A minimum one-month gap separates consecutive events. Seasonal cycles are retained for assessing total magnitudes relevant to flooding and ecosystems. Attribution and filtering: Anthropogenic global mean sea level rise (GMSLR) was estimated from CSIRO and NASA GMSL records (1880–2019) using two methods: 90% of quadratic fits for 1960–1992 and 1993–2019 periods; and a climate-change-induced acceleration (0.084 mm/yr²) for 1993–2019 combined with the earlier fit, yielding nearly identical curves. To isolate variability, anthropogenic GMSLR and an 8-year low-pass filtered decadal SLA component were removed from model SLAs to analyze seasonal-to-interannual signals. Model experiments: Two OGCMs were used. ROMS configured for 25°S–25°N at 1/3° with 40 sigma layers was forced by 3-hourly JRA55-do fields (1958–2017) with open boundaries relaxed to monthly ORAS4, allowing inclusion of global SLR and remote signals. Two ROMS runs: a main run (full forcing) and WSTRESS (only wind-stress variability; surface heat/freshwater fluxes climatological) to isolate wind-driven variability. HYCOM global configuration at 1/2° with 50 hybrid layers, forced by JRA55 (1958–2017); land ice melt contribution to global SLR is not included in HYCOM. Validation included correlations and trend comparisons to altimetry and tide gauge. Additional analyses: A Bayesian Dynamic Linear Model (DLM) quantified contributions of two predictors—equatorial Indian Ocean zonal wind stress anomalies (65°E–95°E, 5°S–5°N; one-month lead) and local longshore wind stress along Sumatra–Java—to Indonesian coastal SLA, with time-varying coefficients estimated via Kalman filtering and smoothing. CESM1 Pacific pacemaker 10-member ensembles (1920–2019) restored observed SST in the central/eastern equatorial Pacific to assess ENSO-driven impacts via atmospheric bridge and Indonesian Throughflow on Indian Ocean sea level. CMIP6 large ensembles were analyzed to estimate external forcing impacts on regional dynamical sea level near Indonesia. A ROMS mixed-layer heat budget was performed to attribute SSTA tendencies to net surface heat flux, horizontal advection and mixing, and subsurface processes (upwelling, vertical mixing, entrainment).

• Rapid coastal sea level rise: Satellite altimetry shows a 5.12 ± 0.17 mm/yr rise near the Java coast (1993–2018), exceeding the global mean 3.1 ± 0.3 mm/yr.
• Event detection: Fifteen HEX events were identified (1993–2018), with 10 clustered in 2010–2017. The strongest monthly HEX occurred in June 2016 with SLA ~0.44–0.45 m (satellite 0.438 m; tide gauge 0.446 m). Six events were CHHEXs (HEX with concurrent MHW), four of which occurred during 2010–2017.
• Spatial characteristics: CHHEX SLA signals extend broadly along Southeast Asian coasts, while associated marine heatwaves are confined to Indonesian coasts and several hundred kilometers offshore. HEX-alone SLAs are weaker and more confined.
• Forcing mechanisms of HEXs: Seasonal-to-interannual wind stress forcing is the deterministic driver of individual HEXs. Equatorial westerly wind anomalies induce eastward-propagating equatorial Kelvin waves that raise sea level; upon reaching the eastern boundary, signals propagate poleward along coasts as trapped waves. Concurrently, local northwesterly longshore winds drive Ekman convergence toward the coast, amplifying coastal SLAs.
• CHHEX versus HEX-alone: CHHEXs are associated with negative IOD years, five of six co-occurring with La Niña. Negative IOD and La Niña produce equatorial westerly and northwesterly longshore wind anomalies that weaken or reverse seasonal southeasterly monsoon winds, reducing coastal upwelling and cooling, deepening the thermocline, and enhancing horizontal advection of warm waters. This yields large interannual warm SST anomalies that, when superimposed on seasonal cycles, create marine heatwaves in May–June and November–December. In contrast, HEX-alone events mainly occur in December–March with strengthened northwest monsoon downwelling but enhanced turbulent heat loss and mixing limiting surface warming.
• Decadal modulation and attribution: The concentration and increased magnitude of HEXs during 2010–2017 arise from the combination of anthropogenic GMSLR and decadal SLA increase due to climate variability; removing anthropogenic SLR and 8-year low-pass decadal variability in models eliminates the enhanced HEX frequency/magnitude. Anthropogenic SLR and decadal variability contribute roughly equally to the enhanced HEX activity. CMIP6 large ensembles indicate the external-forcing effect on regional dynamical sea level near Indonesia is weak (<2 cm) with large uncertainties.
• Source of decadal variability: Prior to 2012, decadal SLAs are largely explained by ENSO decadal variability (La Niña–like), as shown by CESM1 pacemaker ensembles and consistent with intensified Pacific easterlies and enhanced ITF transports. From 2013–2017, a transition toward negative IOD decadal conditions sustains high SLAs via equatorial westerly and longshore northwesterly wind anomalies.
• Heat budget diagnostics: For CHHEXs (1998, 2010, 2016), reduced upwelling and vertical mixing, augmented by horizontal advection, drive interannual warm SST anomalies; in May 2013, increased surface heat flux plus reduced upwelling/mixing contribute.
• Model–data agreement: ROMS, HYCOM, and ORAS4 reproduce observed coastal SLAs with high correlations (ORAS4 0.98, ROMS 0.95, HYCOM 0.90). ROMS trend (6.50 ± 1.16 mm/yr) is within satellite uncertainty (5.59 ± 0.99 mm/yr); ORAS4 and HYCOM underestimate trends.
• Seasonal context: Extremes are absent during July–October when seasonal upwelling lowers sea level and cools SST; even with interannual warming, total SST anomalies remain modest during peak upwelling but can exceed thresholds during IOD onset (June) and decay (Nov–Dec).

The study addressed how anthropogenic sea level rise interacts with climate variability to generate and modulate sea level extremes and their compounding with marine heatwaves along coastal Indonesia. Results show that wind-driven processes tied to ENSO and IOD, operating across seasonal to decadal timescales, underpin both HEX and CHHEX events, while anthropogenic SLR elevates the baseline, increasing the frequency with which interannual anomalies cross extreme thresholds. The clustering of HEXs in 2010–2017 is not random but driven by the superposition of anthropogenic SLR and a positive decadal sea level state linked first to La Niña–like Pacific decadal variability and later to a negative IOD phase. Mechanistically, equatorial westerlies and local northwesterly longshore winds raise sea level through Kelvin wave propagation and coastal Ekman convergence; CHHEXs occur when these winds also weaken the seasonal southeasterly monsoon, suppressing upwelling and enhancing warm advection to produce marine heatwaves coincident with high sea levels. The findings highlight the need to consider compounding hazards at monthly timescales, the value of sustained altimetry and tide gauges, and the importance of decadal climate predictions for risk assessment and adaptation planning in vulnerable coastal regions.

This multi-dataset, multi-model analysis establishes that compound sea level height and marine heat extremes (CHHEXs) along coastal Indonesia arise from the interplay between anthropogenic sea level rise and climate variability spanning seasonal, interannual (ENSO, IOD), and decadal timescales. Anthropogenic SLR and decadal variability contributed approximately equally to the elevated HEX frequency and magnitude during 2010–2017, while negative IOD conditions—often co-occurring with La Niña—set the wind patterns that suppress upwelling and promote marine heatwaves during key seasons (May–June, November–December). Given projections of mean-state changes toward a shallower thermocline in the eastern Indian Ocean and increased IOD amplitude under continued warming, the region may be preconditioned for stronger HEXs and CHHEXs, amplifying ecological and societal risks. Future work should prioritize: improved decadal prediction of ENSO and IOD; refined regional sea level projections including vertical land motion; higher-resolution coupled models resolving coastal processes; and enhanced observing systems (altimetry, tide gauges, coastal SST) to monitor and anticipate compound marine hazards.

• Data and model uncertainties: Tide gauges may include long-period tides and storm surges; altimetry filters some of these signals. Land motion corrections were not applied at the tide gauge due to lack of nearby GPS, introducing potential local bias.
• Trend biases: ORAS4 (coarse 1°) and HYCOM underestimate coastal sea level rise trends; HYCOM omits land ice melt contribution to GMSLR. Spatial resolution and boundary conditions affect coastal variability representation.
• Attribution uncertainty: CMIP6 large ensembles indicate small, uncertain external-forcing effects on regional dynamical sea level near Indonesia. Separating anthropogenic and natural decadal components relies on filtering and model-based estimates.
• Event definitions: MHWs defined using monthly data and the 90th percentile threshold may differ from daily-based definitions; some events have timing nuances (e.g., Dec 2010 CHHEX counted with near-threshold SSTA).
• Model biases: Despite high correlations, models contain structural and forcing uncertainties; CESM1 pacemaker experiments and ROMS/HYCOM configurations have known biases that may affect amplitude and trend fidelity.