Evidence for Hunga-Tonga eruption influence on the anomalously strong Northern-Hemisphere stratospheric polar vortex in 2024/2025

The 2022 Hunga Tonga (HT) eruption injected an unprecedented amount of water vapour (WV) into the stratosphere and mesosphere, emerging as a potential multi-year driver of variability in those layers and associated climate feedbacks. Using satellite and reanalysis datasets, and ensemble simulations with the SOCOLv4 model, with and without the HT forcing, we diagnose the chain of processes linking the eruption to the exceptionally strong Northern Hemisphere stratospheric polar vortex (SPV) observed in winter 2024/2025. Satellite data show “tongues” of enhanced WV descending from the mesosphere into the polar stratosphere, collocated with ozone reductions and persistent cold anomalies, alongside a record-strong SPV. Our model with the HT forcing can reproduce the main structure of the descending plume, its effects on chemistry, and the SPV strengthening, albeit with underestimated magnitudes and an earlier onset than observed. Offline radiative transfer calculations indicate that the WV and ozone anomalies drive net radiative cooling over the polar stratosphere, which steepens meridional temperature gradients and thereby intensifies the Arctic polar vortex. The observed 2024/2025 winter thus represents a plausible manifestation of HT-induced vortex variability. Our simulations also predicted the observed shift towards more favorable conditions for sudden stratospheric warmings in the winter 2025/2026 as the WV forcing declines.

1 Introduction

The 2022 Hunga Tonga (HT) eruption was unprecedented in terms of (i) injection altitude up to 57 km (Proud et al., 2022), (ii) the observational coverage of a direct water vapor (WV) injection of this magnitude (Vömel et al., 2022), and (iii) the radiatively-driven plume descent shortly after the eruption due to the highly localized moist volcanic cloud (Kroll and Schmidt, 2024). Many studies have been dedicated to the impacts of the HT eruption, including changes in stratospheric and mesospheric composition (Evan et al., 2023), radiation (Stocker et al., 2024; Randel et al., 2024b), and direct (Jenkins et al., 2023; Schoeberl et al., 2023; Stenchikov et al., 2025) and indirect (Bednarz et al., 2026; Kuchar et al., 2025) climate impacts.

The WV inputs of HT were distributed throughout the stratosphere, and even to the mesosphere of the Southern Hemisphere subtropics in the days/weeks following the eruption (APARC, 2025). Randel et al. (2024a) documented positive residuals of WV over high latitudes in both hemispheres in late 2023, tied to the downward transport of WV from the HT eruption. Modelling (Fleming et al., 2024; Zhou et al., 2024; Kuchar et al., 2025; Zhuo et al., 2025) and observational (Millán et al., 2022; Schoeberl et al., 2025) studies suggest that the anomalous WV induced by the HT eruption in 2022 should diminish within 5–7 years. The Antarctic dehydration likely remains the dominant WV removal pathway (Zhou et al., 2024), together with transport of stratospheric air to the troposphere at mid-high latitudes by the Brewer-Dobson circulation (BDC). Therefore, we expect that BDC is responsible for the downward transport of WV over middle and high latitudes during the wintertime, in analogy to the downward transport of NO_x, produced by energetic particle precipitation, inside the polar vortex (e.g. Funke et al., 2014).

Kuchar et al. (2025) explored a novel pathway by which water-rich volcanic eruptions, such as HT may modulate stratospheric and mesospheric conditions and influence regional surface climate. Their targeted model simulations illustrated that the excess WV can reduce meridional temperature gradients near the stratopause, leading to a negative Northern Annular Mode anomaly that propagates to the surface in late winter. Several features discussed in Kuchar et al. (2025) were reproduced by an independent atmosphere-only simulation using the NCAR model (Bednarz et al., 2026, and accompanying discussion), including mesospheric cooling at lower latitudes, which, via the thermal wind relation, has the potential to weaken the polar-night jet. Taken together, these studies suggest that HT may have provided favorable conditions for Sudden Stratospheric Warmings (SSWs) observed in NH winters following the HT eruption (Lee et al., 2025).

Figure 1 (A) The average (east-west) zonal wind for 60° N at 10 hPa for the extended winter 2024/2025 (OW3) in MERRA2. This is near the peak of the polar jet maximum. (B) The eddy heat flux (EHF) averaged between 45 and 75° N at 100 hPa for the 45 d before the date indicated.

Various atmospheric (and surface) signals have also been linked to a stratospheric WV increase. Joshi et al. (2006) showed a strengthening of the stratospheric jet in boreal winter associated with the cooler lower stratospheric midlatitudes (increase of meridional temperature gradient) in response to an increased stratospheric WV. Maycock et al. (2013) and Seabrook et al. (2023) reported a similar circulation response. This feature was indeed observed in the winter 2024/2025 when the stratospheric polar vortex was anomalously strong, reaching the climatological maximum (∼ 61.19 m s⁻¹) on 24 December 2024 (see Fig. 1A and Sect. 2.1.1). There were only a few recent winters with a similarly strong and stable polar vortex in the Northern Hemisphere, i.e. 2010/2011, 2015/2016, and 2019/2020 (Manney et al., 2011; Matthias et al., 2016; Lawrence et al., 2020; Rieder et al., 2025).

Motivated by the diverse findings of HT impacts summarized in the recent dedicated report (APARC, 2025), the aim of this study is to document the portfolio of pathways by which the HT eruption may have modulated stratospheric and mesospheric conditions and consequently impacted the surface climate. Specifically, we aim here to answer the following questions: (i) How does the long-term response of the stratospheric polar vortex (SPV) to atmospheric composition changes following the HT eruption differ from the short-term response? (ii) Did atmospheric composition changes precondition SPV states, particularly the strong SPV observed in 2024/2025? and (iii) if so, what predictive power can be derived from comprehensive chemistry-climate model (CCM) simulations for SPV states in the coming years?

To this end, we provide a combined observational and modelling assessment of the 2024/2025 event and its linkages to the lingering Hunga-induced moist plume. We use a reanalysis dataset and observations from MLS and SWOOSH (Sect. 2.1), together with large-ensemble SOCOLv4 simulations with and without the HT forcing (Sect. 3.1), to investigate the sensitivity of the circulation to excess WV from HT and assess the chain of processes leading to SPV strengthening (Sect. 4).

2 Data

2.1 Observational datasets

2.1.1 MERRA2

Here we analyze a widely used reanalysis data set, the National Aeronautics and Space Administration (NASA) Modern Era Retrospective Analysis for Research and Applications Version 2 (MERRA2, resolution 0.625°×0.5°, model top at 0.01 hPa; Gelaro et al., 2017). We extract zonal wind (ZW) for 60° N at 10 hPa and eddy heat flux (EHF) between 45 and 75° N at 100 hPa postprocessed by the NASA Ozone Watch (Kramarova and Lait, 2025) for the extended winter 2024/2025.

2.1.2 MLS

The Microwave Limb Sounder (MLS; Waters et al., 2006) instrument onboard NASA's Aura Satellite has been operational since 2004 and provides near-global (82° S–82° N) measurements of key atmospheric constituents and temperature based on thermal microwave emissions along the Earth's limb. The MLS provides about 3500 vertical profiles per day for each species. The profiles have a vertical resolution of approximately 3 to 5 km in the stratosphere and 6 to 8 km in the mesosphere. In this study, we use Version 6 of the retrieval products for water vapor, ozone, and temperature (Livesey et al., 2026).

Since May 2024, the 190 GHz subsystem in MLS, which is responsible for observing water vapor, N₂O, HCN, and upper stratospheric HNO₃, has been operating on a duty cycle to extend the instrument's lifetime (see https://mls.jpl.nasa.gov/, last access: 22 October 2026). The measurements are reactivated for about six days in the middle of each calendar month. Other MLS channels, including those used for ozone and temperature, remain unaffected.

2.1.3 SWOOSH

The Stratospheric Water and Ozone Satellite Homogenized (SWOOSH) version 2.72 dataset provides zonally averaged monthly profiles for water vapor and ozone. SWOOSH combines measurements from several satellite sensors, including the Stratospheric Aerosol and Gas Experiment (SAGE II & III), the Upper Atmosphere Research Satellite Halogen Occultation Experiment (UARS HALOE), UARS MLS, and Aura MLS version 5.0. The retrieval artifact in the 10–8 hPa layer apparent in SWOOSH (see Fig. 2 and Niemeier et al., 2023; Millán et al., 2024) has been corrected in the MLS data version 6 we use in this study. To account for measurement offsets and biases, the data records were homogenized using “reference” satellites during times in which observations of different satellite systems overlapped. The satellite that showed the best agreement with balloon-based measurements was chosen as a reference (Davis et al., 2016b).

Figure 2 Globally averaged water vapour (WV) anomalies (shading and labeled contours for 0.8 ppmv) in SOCOLv4 (A), MLSv6 (B) and SWOOSH (C) for the period 202201–202511. The 10–8 hPa artifact (Millán et al., 2024) in WV remains apparent in SWOOSH, as it contains data from MLSv5. 0 (Davis et al., 2016a).

3 Methods

3.1 SOCOLv4 simulations

Our set of SOCOLv4 (Sukhodolov et al., 2021) simulations comprises an ensemble of transient simulations with and without the HT forcing, i.e. 150 Tg of WV and 0.5 Tg of SO₂ have been distributed over the emission region (22–14° S, 182–186° E and 25–30 km) within 15 January similarly to the community recommendation (Zhu et al., 2025). We perform a 5-year spin-up before the HT eruption. Both ensembles comprise 10 ensemble members. We stress that the Quasi-Biennial Oscillation (QBO) is nudged in both simulations with and without the HT forcing. We use equatorial zonal wind at 10 hPa as a proxy for Quasi-biennial Oscillation (QBO) to demonstrate a link between QBO and the polar vortex in SOCOLv4, i.e. the so-called Holton-Tan mechanism (Holton and Tan, 1980). For a detailed description of our SOCOLv4 simulations, we refer to Appendix A1 in Kuchar et al. (2025).

As indicated in Kuchar et al. (2025), the anomalous WV following HT decays and propagates vertically too rapidly in SOCOLv4 compared with MLS/SWOOSH (cf. globally averaged WV anomaly in SOCOLv4, MLSv6 and SWOOSH in Fig. 2) or some of the other models (Randel et al., 2024b). Note that the contour at 0.5 ppmv does not reach the year 2024 in SOCOLv4, unlike in SWOOSH. This limitation is well documented for the SOCOL model family: SOCOLv4 retains, as in previous SOCOL versions, a too-fast residual circulation and excessive mixing in the surf zone, and its tropical water-vapour tape-recorder signal propagates upward too rapidly (Stenke et al., 2013; Arsenovic et al., 2019; Sukhodolov et al., 2021). Related stratospheric transport biases are also common in many CCMs, as indicated by too-young age of air (AoA) and overly rapid tropical ascent in CCM intercomparisons (Dietmüller et al., 2018; Abalos et al., 2026), although their magnitude and vertical/latitudinal structure vary across models.

Overall, SOCOLv4 is adequate for broad NO_x/HO_x/ClO_x driven ozone depletion/recovery studies (Sukhodolov et al., 2021). The main remaining weaknesses are transport and mixing, and photolysis at large solar zenith angles over the high latitudes of the Southern Hemisphere (e.g. Imanova et al., 2025).

3.2 Offline radiative transfer simulations

We prescribed WV and ozone anomalies from MLS for the extended winter 2024/2025 (OW3) and calculated longwave and shortwave radiative heating rates using konrad's (Kluft et al., 2019; Dacie et al., 2019) interface to the Rapid Radiative Transfer (RT) Model for General circulation models (RRTMG; Kroll et al., 2021). We state that these radiative calculations are robust for various solar zenith angles. We accompany these calculations using a precise line-by-line RT simulations (see details in Bell et al., 2025)

3.3 Calculation of anomalies

Throughout our analysis, we evaluate significance fields using the minimum local p-values from Student's t-test with global test statistics. Most of the illustrations in Sect. 4 show differences between simulations with and without the HT forcing. For the MLSv6 observation, we compute local p-values as empirical two-sided rank p-values for each target month-pressure bin. We confirm this approach using local Monte Carlo p-values with 1000 samples.

3.4 Eliassen-Palm flux diagnostics

We investigate the response of resolved planetary waves (PWs) using the Eliassen-Palm flux diagnostics (EPF; Andrews and McIntyre, 1987) applied to the SOCOLv4 6-hourly outputs. EP fluxes are computed and scaled following Jucker (2021). The EPF convergence serves as an indicator of wave dissipation, and the EPF divergence (EPFD) indicates sourcing. We illustrate all EPF diagnostics from 1 October 2023 to 30 March 2024 in the Supplement.

3.5 Calculation of Southern Oscillation Index (SOI)

We calculate the Southern Oscillation Index (SOI) using sea level pressure (SLP) standardised difference between Tahiti (17.65° S, 149.45° W) and Darwin (12.46° S, 130.84° E) as a proxy for the El Niño Southern Oscillation (ENSO). Negative and positive SOI values coincide with abnormally warm and cold ocean waters across the eastern tropical Pacific typical of El Niño and La Niña episodes, respectively.

4 Results

4.1 Observed changes in water vapor, ozone, and temperature

We set the scene by illustrating the daily evolution of stratospheric dynamics in the extended winter 2024/2025 (the third observed winter after the HT eruption; thus, marked as OW3) using the MERRA2 reanalysis in Fig. 1. Here, the zonal-mean winds diagnosed at 60° N at 10 hPa illustrate that the SPV strength exceeds the average strength during DJF 2024/2025. Next, we turn the focus to anomalies of WV, temperature, and ozone (derived from the MLS observations in Fig. 3), to put the exceptionally strong SPV in the context of the lingering HT perturbation in polar latitudes. The WV anomalies (see Fig. 3A) are characterized by “tongues” with concentrations of up to 2 ppmv. The first tongue appears in the mesosphere during September 2024 and continuously propagates downward with time to ∼ 10 hPa, where it gradually weakens and dissipates by March 2025, when it merges with the stratospheric background. A second tongue forms in the mesosphere in January 2025 and also propagates downward (not shown for the chosen time frame of Fig. 3). The WV anomaly is statistically significant (at respective p-values for given figure) over much of the descending tongue. In contrast, the negative ozone and temperature anomalies are statistically significant only locally, mainly along the upper flank of the WV tongue starting around 0.1 hPa in September and ending below 10 hPa in February (see Fig. 3B, C), reflecting the large daily and interannual variability of polar wintertime conditions. We therefore do not interpret the 2024/2025 anomalies as uniquely attributable to HT; rather, the observations indicate a plausible HT contribution superimposed on internal dynamical variability. The temperature anomalies are followed by their positive counterparts. By October, an enhanced ozone anomaly (∼0.5 ppmv) emerged around 10 hPa, leading to effective warming due to ozone absorption of solar radiation and increased mixing and thereby masking the WV-induced cooling propagating downward. In March, we observe a significant increase and decrease in stratospheric and mesospheric polar temperatures, respectively. However, the stratospheric ozone anomalies have both alternating positive and negative bands. This corresponds to the stratospheric final warming during which the SPV does not recover and transitions to summer conditions, with easterly winds across the stratosphere (see Fig. 1A). The observed stratospheric and lower-mesospheric longwave cooling rates are dominated by CO₂ together with ozone and WV (Clough and Iacono, 1995; Brasseur and Solomon, 2005). The close correspondence between the WV and ozone anomalies and the diagnosed temperature response suggests the radiative cooling observed in OW3. However, a comparison with OW2 (see Fig. 4) implies that slowly descending WV, i.e. dynamical preconditioning (see EHF in Fig. 1B), positively feeds back on radiative cooling and strengthens the vortex. We examine this mechanism explicitly in the model simulations below.

Figure 3 Zonal mean weighted average for latitudes between 60–90° N of water vapour volume mixing ratio (VMR; shading in A and contour of 0.8 ppmv in B and C; in ppmv; monthly), ozone VMR (shading in B, in ppmv; daily) and temperature (shading in C; in K; daily) and anomalies for the extended winter 2024/2025 (OW3) from MLS. Anomalies are expressed as differences with respect to the climatology for the period 2006–2020. Hatching $\backslash \backslash \backslash \backslash$ and $////$ show where the empirical two-sided p-values are <0.1 and <0.05, respectively. Due to missing data, water vapour is visualized as monthly anomalies with interpolated gaps.

Figure 4 As Fig. 3 but for the extended winter 2023/2024 (OW2).

4.2 Modeled changes in water vapor, ozone, and temperature

The downward propagating tongue of WV during winter is captured in our SOCOLv4 ensemble. The feature emerges in the extended winter 2023/2024 (marked as the second model winter after the HT eruption, MW2; see Fig. 5A) and winter of 2024/2025 (MW3; see Fig.  6A). The WV enhancement is, however, underestimated compared to observations in Fig. 3A, especially in MW3. Given the well-documented model deficiency of accelerated transport of SOCOLv4 (discussed in detail in Sect. 3.1), we argue that the winter MW2 in our simulations is comparable to OW3 in observations (see Sect. 4.1 and Fig. 3). For the MW2 SOCOLv4 yields a similar strengthening of the vortex with a less pronounced weakening than observed at the end of OW3 (cf. Fig. 5D at 10 hPa and Fig. 1A). In analogy, we assume MW3 matches with the upcoming OW4, as we detail further below.

Figure 5 Weighted zonal average over 60–90° N of water vapour volume mixing ratio (VMR; shading in A and contours of 0.5 and 1 ppmv in A, B, C and D; in ppmv), ozone VMR (shading in B, in  ppmv), temperature (shading in C; in K), and zonal wind at 60° N (shading in D) and Eddy Heat Flux at 100 hPa averaged over 45–75° N (EHF in m s⁻¹ K; red line in C) daily anomalies for the extended winter 2023/2024 (MW2) from SOCOLv4. Anomalies are expressed as the difference between the SOCOLv4 simulation with and without the HT forcing. Hatching $\backslash \backslash \backslash \backslash$ and $////$ show where the p-values of the t-test are <0.1 and <0.05, respectively. Gray contours in (D) represent absolute values of zonal wind in the SOCOLv4 simulation with the HT forcing.

Figure 6 As in Fig. 5 but for the extended winter 2024/2025 (MW3).

For MW2 (equaling OW3) we attribute the initial ozone destruction in the upper stratosphere and lower mesosphere to increased OH concentrations induced by the excess WV from the HT eruption (HO_x cycle) and the ClO_x activation from chlorine reservoirs (see Fig.  7). These coincide with a negative anomaly in NO_x mainly attributable to cooling at these levels (Fleming et al., 2024). In contrast, as reported in Fleming et al. (2024); Kuchar et al. (2025), the positive ozone anomaly limited to September–November and rather in the lower to mid stratosphere, similar to MLS, emerges due to increased conversion of NO_x to the HNO₃ reservoir given an overall higher abundance of OH. While negative ozone and negative temperature anomalies are not statistically significant (at respective p-values for given figure) along the entire WV tongue similarly as in observations (see Fig. 3), they are associated with a stronger SPV and correspondingly weaker poleward transport of ozone from lower latitudes reflected by the EHF (red line in Fig. 5D) and AoA (Fig. 8). Thus, the paired SOCOLv4 simulations indicate that HT forcing shifts the ensemble-mean circulation response toward stronger polar westerlies. However, this ensemble-mean response does not imply that individual ozone, temperature, or wind anomalies, or the full observed evolution of OW3, can be uniquely attributed to the HT forcing.

Figure 7 Monthly zonal-mean anomalies of polar chemistry for the period 2023/2024 (MW2). Anomalies are expressed as the difference between the SOCOLv4 simulation with and without the HT forcing. Hatching $\backslash \backslash \backslash \backslash$ and $////$ show where the p-values of the t-test are <0.1 and <0.05, respectively.

Figure 8 Weighted zonal average over 60–90° N of age of air estimated as idealized age tracer volume mixing ratio monthly anomalies for the extended winter 2023/2024 (MW2). Anomalies are expressed as the difference between the SOCOLv4 simulation with and without the HT forcing and normalized to its climatological standard deviation (SD). Dots … and hatching $\backslash \backslash \backslash \backslash$ and $////$ show where the p-values of the t-test are <0.3 and <0.1 and <0.05, respectively. The black contour line represents the contour of 1 ppmv of water vapour anomaly.

Next, we investigate the radiative effect of the observed changes in WV and O₃ on temperature via offline radiative transfer (RT) calculations (see Sect. 3.2). These RT calculations presented in Fig. 9 reveal that perturbations in both WV and ozone abundances are capable of inducing cooling rates in the order of −0.5 to −1 K d⁻¹, with WV and ozone dominating the change in the longwave and shortwave heating rates throughout the stratosphere and in the upper layers, respectively. These heating rate anomalies are collocated with those observed (see Fig. 3C) and simulated (see Fig. 5C). SOCOLv4, however, does yield a weaker temperature response initially, attributable to weaker changes in shortwave heating rates following the weaker ozone anomaly than seen in the MLS observations.

Figure 9 RRTMG radiative heating rates (longwave, shortwave and their sum) for prescribed weighted zonal average over 60–90° N of water vapour and ozone monthly anomalies from MLS for the extended winter 2024/2025 (OW3).

Figure 10 Weighted zonal average over 60–90° N of clear-sky (offline) longwave heating rate (in K d⁻¹) anomalies calculated using solely WV (shading in A); WV and O₃ (shading in B). Contour in all panels represents a WV anomaly (0.5 ppmv) as in Fig. 5. Anomalies are expressed as the difference between the SOCOLv4 simulation with and without the HT forcing. Hatching $\backslash \backslash \backslash \backslash$ and $////$ show where the p-values of the t-test are <0.1 and <0.05, respectively.

We confirm these calculations using the line-by-line RT simulations (see Sect. 3.2). The resulting longwave heating rates from prescribed time-evolving WV and ozone from SOCOLv4 are shown in Fig. 10.

As detailed in Bloxam and Huang (2021), the radiative-cooling recovery of SSW-induced temperature anomalies in the absence of dynamical heating is dominated by longwave radiative cooling, and is characterized by increasing e-folding times from the upper to middle stratosphere. If we invert the e-folding times (4–7.5 d in the upper stratosphere) of Bloxam and Huang (2021), we obtain an equivalent heating-rate effect of roughly 0.13–0.25 K d⁻¹ K⁻¹. Given that the e-folding time increases in the lower stratosphere and decreases towards the mesosphere, we expect that temperature perturbations translate into weaker and stronger heating-rate effects at lower and upper levels, respectively. In the lower stratosphere, the air cools more slowly because there is a deeper layer of atmosphere above it that absorbs the radiation, re-emits it upwards but also downwards, and thereby effectively reduces the amount of energy that can escape into space. Closer to the tropopause, the cooling effect reduces as the radiation emitted depends strongly on the background temperature itself (Bloxam and Huang, 2021). This asymmetry may explain the different persistence of the heating rate anomalies at different heights that appear in Fig. 10B. Furthermore, the pronounced positive heating rate anomaly around 1 hPa spanning from December 2024 to March 2025, consistent with enhanced longwave heating related to decreased ozone, aids in explaining the rapid change from negative to positive temperature anomalies (see Figs. 3C and 5C).

The negative temperature anomalies throughout the upper stratosphere and lower mesosphere result in an enhanced meridional temperature gradient (see Fig. 14), which strengthens SPV via the thermal wind balance. The SPV strengthening may initially enhance vertical wave propagation during October, a time when ZW are usually weak and thus allow more PWs to propagate upwards into the stratosphere (Charney and Drazin, 1961). In this respect, we understand the simultaneously observed anomalously high EHF at 100 hPa (see Fig. 1B) as a proxy for upward propagation of PWs (e.g. Newman et al., 2001). This way we associate the anomalously weak ZW below in the stratosphere to anomalously breaking PWs (see visualization of EPF in the Supplement, Sect. 3.4). As the cooling progresses downward, it aids the strengthening of the winds, which tend to be too strong for wave propagation (Plumb, 1989). This may confirm the nonlinearity of SPV for generally large wave forcing during the NH winter conditions (Scott and Haynes, 2002). The stronger SPV is associated with pronounced negative and positive temperature anomalies below and above, respectively. A similar dipole structure has been discussed in studies diagnosing SH high-latitude temperature trends as a consequence of ozone depletion (e.g. Young et al., 2013; Keeble et al., 2014). Manzini et al. (2003) attributed the positive temperature anomaly at upper levels to increased downwelling, and thus, consequently, adiabatic warming due to enhanced gravity wave propagation. Correspondingly, stronger westerlies in the stratosphere filter eastward gravity waves while westward gravity waves propagate freely upward and deposit their momentum above.

4.2.1 QBO and ENSO preconditioning

As for MW2, we find similar “tongues” of WV and downward propagation of ozone in the SOCOLv4 simulations for MW3 (see Fig. 6). Here, however, the excess of WV only yields up to 0.25 ppmv, and thus a much smaller magnitude compared to OW3 by MLS. This reduced WV forcing results only in a subtle cooling, which slightly strengthens the ZW, before they significantly weaken as a consequence of anomalously high EHF in December 2024 and the beginning of January 2025. This episode is followed by a recovery of SPV and collocated negative temperature anomaly propagating from the lower mesosphere to the lower stratosphere.

Figure 11 Daily zonal-mean zonal wind at 10 hPa and 60° N based on the SOCOLv4 simulation with (ensemble mean (EM); red solid line) and without (black solid line) the HT forcing for the period 2024/2025 (MW3). The SOCOLv4 simulation with the HT forcing reveals two clusters: weak SPV (red dotted line) and strong SPV (red dashed line).

Figure 12 Scatter plot of equatorial zonal wind at 10 hPa as a proxy for Quasi-biennial Oscillation (QBO) in SOCOLv4 and Southern Oscillation Index (SOI) in the simulation with the HT forcing averaged over SON. Individual dots correspond to each of the 10 ensemble members in the corresponding year and zonal wind (ZW) anomalies at 60° N averaged over DJF.

Figure 13 Monthly anomaly of sea level pressure (SLP in hPa) in February 2025 (MW3). Anomalies are expressed as the difference between the SOCOLv4 simulation with and without the HT forcing. Dots … and hatching $\backslash \backslash \backslash \backslash$ and $////$ show where the p-values of the t-test are <0.3 and <0.1 and <0.05, respectively.

The episode with the weakened ZW as a consequence of anomalously high EHF in December 2024 (MW3 in Fig. 6) and the beginning of January 2025 can be traced back to two main clusters of ensemble members in the simulations with the HT forcing (see Fig. 11). While one major cluster of 4 ensemble members (with mixed Southern Oscillation Index (SOI) preconditioning) strengthens SPV, the other cluster of 6 ensemble members (with mostly positive SOI (see Fig. 12) corresponding to La Niña preconditioning) weakens the stratospheric winds leading to the emergence of SSWs between December 2024 and January 2025. Considering the easterly phase of QBO at 10 hPa for MW3, we hypothesize that this may have preconditioned the latter cluster to a higher probability of weaker SPV or SSWs, respectively. Such interpretation is consistent with the observation that SSWs are more likely to occur during easterly QBO in tandem with either El Niño or La Niña conditions (Lee et al., 2025).

That said, the simulated equatorial winds in the westerly QBO phase during the early stage of MW2 (see Fig. 12) create conditions that favor a stronger SPV, consistent with the Holton–Tan mechanism (Holton and Tan, 1980). A similar configuration occurred in OW3, when the QBO was also in its westerly phase (https://www.climate.gov/news-features/blogs/polar-vortex/predicting-chances-polar-vortex-disruption-winter, last access: 22 October 2026). Since the QBO preconditioning is identical in both simulations (with and without the HT forcing) and ENSO is uniformly sampled in SOCOLv4, encompassing both El Niño (SOI<0) and La Niña (SOI>0; as ultimately observed in OW3), we infer that the forced SPV response in MW2 arises from the HT-induced WV compared to OW3 where perhaps both QBO and the HT forcing contributed to the stronger SPV.

From the SOCOLv4 ensemble, yielding one cluster with easterly QBO and La Nina conditions in MW3, we hypothesize that the winter 2025/2026 (OW4) bears an enhanced probability of SSW occurrence. As the associated wind pattern propagates downward, we consequently identify a SLP pattern characteristic of a weaker SPV, i.e., a positive SLP anomaly at the poles and a negative SLP anomaly at the mid-latitudes (see Fig. 13).

5 Discussion and summary

Our study addresses the question of whether the anomalously strong Northern Hemisphere stratospheric polar vortex of the boreal winter 2024/2025 has been linked to lingering dynamical and radiative impacts of the HT eruption. For this time period, satellite data show a coherent descent of HT-related water-vapour anomalies into the high-latitude stratosphere, collocated with a marked reduction in ozone and pronounced cold anomalies in the upper stratosphere and lower mesosphere. We detail that an ensemble of SOCOLv4 simulations, with the HT forcing, is able to reproduce the main spatial patterns of this WV plume and associated ozone response. Furthermore, through offline radiative transfer calculations, we show that the WV and ozone anomalies are, in combination, sufficient to generate substantial radiative cooling, which sharpens the meridional temperature gradient (see Fig. 14) and strengthens polar stratospheric westerly winds through thermal-wind balance.

Figure 14 As Fig. 5, but daily anomalies of the weighted meridional temperature gradient (K) between lower (0–20° N) and polar (60–90° N) latitude bands for the extended winter 2023/2024 (MW2) from SOCOLv4. Anomalies are calculated as the difference between the SOCOLv4 simulation with and without the HT forcing. Hatching $\backslash \backslash \backslash \backslash$ and $////$ indicates where the p-values of the t-test are <0.1 and <0.05, respectively. Black contour line represents the zero anomaly.

Figure 15 As Fig. 3 but for the extended winter 2025/2026 (OW4).

The proposed mechanism, in which excess WV from HT and associated ozone changes drive radiative cooling that preconditions the polar vortex, is partly analogous to the well-documented Southern Hemisphere response to Antarctic ozone depletion in the lower stratosphere (Haase et al., 2020) and effects driven by energetic particle precipitation from mesospheric levels downward (Seppälä et al., 2025). Here, however, excess stratospheric WV from HT is considered as a joint driver together with slow Brewer–Dobson descent of the HT WV anomaly, enhanced HO_x and ClO_x chemistry leading to upper-stratospheric/mesospheric ozone loss, and the resulting longwave and shortwave cooling where WV and ozone anomalies overlap. This cooling at the high latitudes steepens the meridional temperature gradient and, via the thermal-wind balance, favours a strong-vortex regime. At the same time, the response is strongly nonlinear: early-autumn wave activity can still transiently weaken the vortex, but once the radiatively strengthened westerlies have reached critical strength, additional wave forcing shows limited impact.

Taken together, the observations and paired SOCOLv4 simulations support a probabilistic interpretation: HT-related WV and ozone anomalies increased the likelihood of a stronger polar vortex while internal variability remains too large to attribute the full evolution of the observed and simulated winters to the HT forcing. Negative ozone and temperature anomalies can also arise dynamically in other strong-vortex winters; the distinctive feature of OW3 is the statistically significant HT-related WV tongue, partly overlapping with these anomalies. The weaker and more variability-dominated MW3 response further indicates that both composition forcing and dynamical preconditioning matter.

Several limitations temper the interpretation and generality of these findings. For one thing, the too-fast transport in SOCOLv4 causes the model to simulate the analogous 2024/2025 winter approximately one year earlier. SOCOLv4 also underestimates the amplitude of the HT-induced WV anomaly when compared with satellite records, implying that our simulated radiative cooling and dynamical response are likely conservative. However, SOCOLv4 agrees well with other models (Zhuo et al., 2025) in terms of the water vapor e-folding time (cf. Kuchar et al., 2025, Appendix A1). Further, the radiative cooling and heating estimates rely on idealized offline calculations and thus cannot fully capture the coupling between radiation, atmospheric composition, and dynamics. In addition, our analysis is also based on a single model with a specific configuration (Sukhodolov et al., 2021, and references therein). Also, internal variability associated with ENSO and other tropospheric sources of internal variability may be only partially sampled by 10 ensemble members (e.g. Bittner et al., 2016a; Bednarz et al., 2026). On the other hand, a sufficiently strong (volcanic) forcing can produce a robust, stronger-than-average NH polar vortex (Bittner et al., 2016b).

Despite these caveats, our findings support the view that HT-induced WV and ozone anomalies have become an important multi-year driver of high-latitude winter variability. In fact, our current study and previous research detail a clear dichotomy: manifesting on the early-time scale as SPV weakening (Kuchar et al., 2025) and on the long-time scale as SPV strengthening. Furthermore, our ensemble simulations indicated that as the WV plume continues to descend and QBO/ENSO conditions evolve, the balance between strong-vortex winters and SSW-prone winters shifted once more in 2025/2026 because of the weaker WV forcing (see recent observed WV (less than 1 ppmv), ozone and temperature anomalies in Fig. 15), with corresponding impacts on surface patterns. More generally, unprecedented eruptions like HT may act as intermittent “sources of predictability” also called “windows of forecast opportunity” for polar winter circulation (Rupp et al., 2026, and references therein), albeit in a fundamentally probabilistic sense.

Future work should therefore (i) extend the analysis to coordinated multi-model HT experiments as in APARC (2025), (ii) investigate how volcanic moistening interacts with ongoing anthropogenic climate change to shape the statistics of SPV strength and SSW occurrence over the coming decades. Better constraining the structure and transport of the descending WV plume using MLS on NASA's Aura satellite may be problematic due to the “duty cycle” of the 190 GHz observations to conserve remaining life for these measurements. Its decommissioning will cause a substantial gap in the measurement of critical atmospheric components such as WV and poses significant challenges for continuous monitoring and assessment of the effects of future extreme events including volcanic eruptions (Salawitch et al., 2025). However, other satellite (Knowland et al., 2025; Himes et al., 2026) and ground-based (De Mazière et al., 2018) instruments can continue to monitor WV after Aura MLS, respectively.