The Estimated Climate Impact of the Hunga Tonga‐Hunga Ha'apai Eruption Plume

On 15 January 2022, the Hunga Tonga‐Hunga Ha'apai (HT) eruption injected SO2 and water into the middle stratosphere. The SO2 is rapidly converted to sulfate aerosols. The aerosol and water vapor anomalies have persisted in the Southern Hemisphere throughout 2022. The water vapor anomaly increases the net downward IR radiative flux whereas the aerosol layer reduces the direct solar forcing. The direct solar flux reduction is larger than the increased IR flux. Thus, the net tropospheric forcing will be negative. The changes in radiative forcing peak in July and August and diminish thereafter. Scaling to the observed cooling after the 1991 Pinatubo eruption, HT would cool the 2022 Southern Hemisphere's average surface temperatures by less than 0.037°C.


10.1029/2023GL104634
2 of 9 strengthening Brewer-Dobson circulation (Xia et al., 2019).Banerjee et al. (2019) analyzing CMIP5 models computed the stratospheric water vapor component of the climate feedback to be 0.14 W/m 2 /K for 4 × CO 2 .Li and Newman (2020) using the Goddard Earth Observing System Chemistry-Climate model computed a similar water vapor feedback value of 0.11 W/m 2 /K.A much smaller response (0.02-0.03W/m 2 /K) was found by Huang et al. (2016Huang et al. ( , 2020)).Note that these model studies are evaluating the long-term climate-system response where non-atmospheric systems (e.g., the ocean, cryosphere) have time to equilibrate.The short-term atmospheric response to sudden forcing changes may be larger because the system is out of equilibrium (Dessler & Zelinka, 2015).
Given the observed climate sensitivity to stratospheric water vapor (Dessler et al., 2013), it is logical to assume that HT might have a climate impact.Jenkins et al. (2023) used a parameterized climate-response model to investigate the impact of the HT water vapor plume.They neglected the impact of aerosols and only considered the radiative forcing due to the water vapor.Jenkins et al. (2023) computed a 0.12 W/m 2 increase in tropospheric radiative forcing.M22 arrived at a similar number, 0.15 W/m 2 .On the other hand, Sellitto et al. (2023) and Zhu et al. (2022) added the direct aerosol forcing and estimated that the plume would produce a peak forcing of −1 to −2 W/m 2 .This exceeds the estimated H 2 O IR forcing.Clearly, both the net warming due to the H 2 O and cooling due to the aerosol layer need to be considered.
In this study we extend the computation of the radiative forcing by Zhu et al. (2022) and Silletto et al. (2022) combining the H 2 O and aerosol radiative forcing.We compute the downward flux change due to stratospheric water vapor using a radiative transfer model.Because of the complexity of computing the aerosol forcing, we take a different approach to estimate the reduction in solar flux.We first use OMPS-LP measurements of stratospheric aerosol extinction to compute the stratospheric aerosol optical depth (AOD).We then convert the AOD to changes in direct radiative forcing using a parameterization based on AOD estimates and direct forcing from previous volcanic eruptions.This approach is also used in climate assessment models (e.g., Hansen et al., 2002.To check our parameterization, we also use the aerosol direct radiative forcing parameterization in Yu and Huang (2023) (hereafter YH).YH provides climatology kernels that can be used to convert AOD to aerosol direct forcing under both clear and all sky (cloudy) conditions.

Data Sets
We use Microwave Limb Sounder (MLS) V5 for temperature and H 2 O measurements.The data quality for the HT anomaly is detailed in M22 and MLS data is described in Livesey et al. (2021).We restrict our constituent analysis to below 35 km, which is roughly the maximum height of the plume a few weeks after the eruption.
We use the aerosol extinction data from OMPS-LP, level-2 V2.1.The V2.1 data (Taha et al., 2022) provides the most accurate OMPS-LP aerosol retrieval up to 36 km.Although the extinction measurements by OMPS-LP are generally consistent with those made by SAGE III/ISS, the OMPS-LP algorithm may overestimate the aerosol extinction below the aerosol peak (Bourassa et al., 2023).AOD is computed from OMPS-LP extinction at 600 nm by integrating the extinction from 36 km to the tropopause height included in the OMPS-LP files.The AOD is converted from 600 to 550 nm assuming an Ångstrom exponent of 1.0 which was derived from SAGE III/ISS HT observations (Taha et al., 2022).The daily MLS and OMPS data sets are averaged onto a 5° × 10° latitude-longitude grid.
We show the SO 2 eruption data available from the NASA Multi-Satellite Volcanic Sulfur Dioxide L4 Long-Term Global Database produced as part of the NASA MEaSUREs project (Carn et al., 2016).

Analysis
In the sections below, we describe our approach to estimating the changes in tropospheric forcing due to HT aerosols and stratospheric water vapor.

Parameterization of the Direct Solar Radiative Forcing by Aerosols
Figure 1a shows the variations in OMPS-LP AOD during 2022; both the maximum AOD and global average AOD are shown.The figure shows a rapid increase in maximum AOD following the eruption which is followed 10.1029/2023GL104634 3 of 9 by a slower growth rate until mid-April.Anomalies and gaps in the late July and early August data are due to a spacecraft anomaly.Although the peak AOD reaches 0.05, the global average AOD only reaches ∼0.02 because most of the aerosol stays within the SH (S22; Taha et al., 2022).Simulations by Zhu et al. (2022) show that SO 2 is converted to sulfate aerosols rapidly following the eruption-a process enhanced by the abundance of water vapor in the plume.After mid-April, dispersion of the aerosols combined with settling cause a slow decrease from the maximum AOD.
H2002 assumed that the aerosol driven change in the direct solar radiative forcing (∆A) can be approximated by ∆A = −RAOD, R = 21.Yu et al. (2023) used a similar relation to account for ∆A due to volcanoes and smoke and found R = ∼23.YH derived direct radiative forcing kernel maps from reanalysis to convert AOD to ∆A, but avoided large volcanic eruption periods.The YH global average clear sky conversion factor R is 29.4,and for all skies, R = 15.7.
The stratospheric AOD changes associated with volcanic eruptions can be much larger than the observed AOD perturbations in the 2000-2022 period YH analyzed.We therefore have independently derived our own parameterization from volcanic analyses.Table 1 shows a list of major observed eruptions and one simulated eruption (Aubry et al., 2021) along with their estimated SO 2 emission, the maximum globally averaged AOD, and the maximum global direct forcing (ΔΑ W/m 2 ).Using these AOD values and our estimates, we can compute the HT direct forcing.We set the background stratospheric AOD is set to 0.012 which is an offset since we want to  compute ΔA(volcaniconly).Figure 1b shows the estimated change in solar flux, ΔA, versus AOD for the data in Table 1 with each volcano listed.We also show a linear fit and log e (AOD) fit to volcanic AOD, H2002 parameterization, and both YH parameterizations.We find R = 19.5 for our linear fit, which is close to the H2002 R value of 21 and the Yu et al. ( 2023) R value of 23. Figure 1b shows that the log e fit better reproduces the estimated changes in ΔA compared to the linear fit.The log e fit is Table 1 also shows that YH clear sky parameterization produces a result very close to the results using Equation 1.

Changes in Direct Solar Forcing Due To Aerosols
For the maximum HT global average AOD value of 0.02, (1) gives a peak global decrease in solar flux (ΔA) of −0.64 W/m 2 (Figure 1b, red dot; Table 1).Using the H2002 parameterization we estimate the change is −0.42 W/m 2 .The YH parameterization gives −0.54 W/m 2 for clear skies and −0.29 W/m 2 for all skies.Zhu et al. ( 2022)'s estimate of −0.25 W/m 2 is slightly lower, because they averaged the forcing in February before the peak AOD in March-April and their simulated aerosol cloud is less extensive than that observed by OMPS-LP (their Figure 3).
Figure 2a shows the surface area weighted zonal mean AOD and aerosol solar radiative forcing due to HT (Figure 2b) from Equation 1.Although Equation 1 is derived using global averages, there is no reason to believe it would not be valid for local changes in solar flux because the approximation simply links stratospheric AOD directly to solar flux.As a check on this assumption, we also perform the calculation using YH kernel maps and these calculations are also shown in Figures 2c and 2d.Our parameterization and YH clear sky are nearly identical, suggesting that YH can be extended to volcanic events the size of HT. Figure 2d shows that if clouds are included the direct solar decreases by about a factor of ∼2.From hereafter we will use the YH parameterization.
The evolution of the direct solar forcing follows the spatial changes in AOD as might be expected, and the decrease in solar forcing is mostly confined to the SH.Our estimates of solar flux changes are in good agreement with estimates by Silletto et al. (2022).Figure 2 also shows a southward shift in the aerosol distributions near day 150 (May 30), and this is reflected in changes in forcing.Sellitto et al. (2023) also shows this southward shift (their Figure 1b).

Changes in IR Radiative Forcing Due To Water Vapor
HT produced an enhanced, mostly SH, stratospheric water vapor layer that mostly extends between 22 and 30 km.This layer generates additional downward LWIR flux and also slightly reduces the solar flux due to water vapor radiative absorption in the SWIR bands.We use the radiative transfer model (RTM) described by Mlawer et al. (1997) to compute the downward radiative flux changes produced by this layer using MLS observed trace gases and temperatures.
To quantify the flux changes, we first compute a daily climatology of MLS temperature and trace gases using 2016-2021 data.We calculate the difference between the 2022 downward tropopause fluxes and the downward fluxes computed using the climatology.We have not applied any adjustment to temperature (e.g., fixed dynamical heating) due to water vapor cooling in the radiative forcing calculations, as was done in previous work (Dessler et al., 2013;Solomon et al., 2010).Solomon et al. (2010) shows that the instantaneous forcing and adjusted forcing are very similar above 22 km.Since HT's water vapor was mostly injected above that altitude, the adjustment is unimportant.Jenkins et al. (2023) also did not perform a temperature adjustment in their estimate of HT's LWIR radiative forcing.
Figure 3a shows the 25 km zonal mean MLS water vapor versus latitude.Figure 3b shows the instantaneous tropopause downward LWIR flux change due to the observed H 2 O distribution, ΔLWIR.In Figures 3b-3d, the latitude range is restricted because the flux estimates at the poles are very noisy due to tropopause fluctuations.Figure 3c shows ΔSWIR due to the attenuation of water vapor (note the sign change in the color bar).Figure 3d shows the net change in H 2 O radiative forcing ΔLWIR + ΔSWIR.The changes in the downward radiative flux follow the evolution of the stratospheric water vapor distribution.There is a small northward shift in aerosols and water vapor shortly after the eruption, and an additional small northward shift of the water vapor distribution in April associated with the QBO (Schoeberl et al., 2023).
The question remains: how much of the LWIR forcing increase is due to H 2 O and how much is due to the temperature differences between 2022 and the climatology?To quantify this, we take the 2016-2021 temperature climatology and substitute the 2022 water vapor.We also and take the 2022 temperature field and substitute in the 2016-2021 water vapor climatology.These two experiments help isolate the impact of the HT water vapor from natural temperature fluctuations.We find that about 1/3 of the 2022 LWIR flux increase through April is due to the descending QBO thermal anomaly and 2/3 is due to the stratospheric water vapor increase.4e) and 1 December 2022 (Figures 4d and 4e).For both clear and all sky ΔA estimates, the aerosol direct forcing overwhelms the heating from stratospheric water vapor.Figure 4g shows the hemispheric and global average forcing time series.The total SH radiative forcing peaks in June/July.On the other hand, the much smaller Northern Hemisphere (NH) forcing peaks in April/May after the QBO has spread aerosols and trace gases into the NH (Schoeberl et al., 2023).the reduction in solar flux due to HT aerosols using past volcanic events (Table 1).Our results are nearly identical to Yu and Huang (2023) results for clear skies, non-volcanic conditions.We subsequently use their clear and all sky parameterizations to assess the direct solar forcing.The solar flux reduction by aerosols is larger than the net IR flux increase due to stratospheric water vapor.In other words, the direct solar radiative cooling associated with the HT aerosols overwhelms the enhanced thermal radiation from stratospheric water vapor plume.Our results are in good agreement with net radiative forcing changes estimated by Silletto et al. (2022) and Zhu et al. (2022).We find that the zonal mean peak change in net radiative forcing occurs in May 2022, but the SH average forcing peaks in June/July as the constituents spread throughout the SH.Using the observed impact on tropospheric temperatures from Pinatubo as a scale, Hunga-Tonga would produce an SH annual average surface temperature change of less than −0.038°C for clear skies and −0.021°C for all skies.

Figure 1 .
Figure 1.Part a, OMPS-LP measured AOD versus day number following the HT eruption in 2022.Crosses are the maximum AOD, blue crosses are the global average.A spacecraft anomaly resulted in missing measurements from late July to mid-August.Black lines show linear fits to the data.Part b, changes in solar forcing with AOD for volcanic eruptions shown in Table 1.The log e fit Equation 1 is the thick solid line.Linear fits from H2002 and Table 1 data are shown as dashed and thin solid lines.The YH clear and cloudy fits are shown in blue and green lines.Volcanic events are small crosses next to the names.The red dot indicates HT (Part a, blue curve).

Figure 4 a
Figure 4 a,b show the time series of the net change in radiative forcing (ΔLWIR + ΔSWIR + ΔA) using the YH clear and all sky parameterizations for ΔA.We also show two latitude cross sections on April 15 (Figures 4cand 4e) and 1 December 2022 (Figures4d and 4e).For both clear and all sky ΔA estimates, the aerosol direct forcing overwhelms the heating from stratospheric water vapor.Figure4gshows the hemispheric and global average forcing time series.The total SH radiative forcing peaks in June/July.On the other hand, the much smaller Northern Hemisphere (NH) forcing peaks in April/May after the QBO has spread aerosols and trace gases into the NH(Schoeberl et al., 2023).

Figure 2 .
Figure 2. Part a, changes in stratospheric OMPS-LP AOD versus time during 2022.Part b, change in solar forcing due to AOD using Equation 1. Part (c, d) shows the forcing using the YH kernels for clear sky and all skies.Vertical lines divide months with initial shown, dotted line is the latitude of HT, dashed line locates the equator.

Figure 4 .
Figure 4. Part a, net radiative forcing using YH clear sky solar forcing (Figure 2c) and IR forcing change (Figure 3d).Part b, same as part a using YH all sky solar forcing.Black contours, 0, 0.5, 1.0, 1.5.Parts c, e show the components of the forcing versus latitude on 4/15/2022 and parts d, f for 12/1/2022, clear and all sky forcing.Part g shows hemispheric average and global average forcing, thick black line is global, red line is NH and thin line is SH, dashed lines for clear, and solid for all sky.

Table 1 .
The log e fit Equation 1 is the thick solid line.Linear fits from H2002 and Table1data are shown as dashed and thin solid lines.The YH clear and cloudy fits are shown in blue and green lines.Volcanic events are small crosses next to the names.The red dot indicates HT (Part a, blue curve).Estimated Emission SO 2 , the Maximum Globally Averaged AOD (550 nm) and Decrease in Global Solar Flux for the Indicated Large Volcanic Events Kloss et al. (2021)he change in direct forcing Equation 1 and the YH parameterization.We use the Raikoke AOD maps shown inKloss et al. (2021)to estimate the Raikoke AOD.SO 2 amounts used by the authors shown; amounts in parenthesis are from the NASA database, Carn et al. (2022), and Sellitto et al. (2023) for HT.Table 1