Revisiting the observed surface climate response to large volcanic eruptions.

. In light of the range in presently available observational, reanalysis and model data, we revisit the surface climate response to large tropical volcanic eruptions from the end of the 19th century until present. We focus on the dynamically driven response of the North Atlantic Oscillation (NAO) and the radiative-driven tropical temperature response. Using 10 different reanalysis products and the Hadley Centre Sea Level Pressure observational dataset (HadSLP2) we conﬁrm a positive tendency in the phase of the NAO during boreal winters following large volcanic eruptions, although we conclude that it is not as clear cut as the current literature suggests. While different reanalyses agree well on the sign of the surface volcanic NAO response for individual volcanoes, the spread in the response is often large ( ∼ 1/2 standard deviation). This inter-reanalysis spread is actually larger for the more recent volcanic eruptions, and in one case does not encompass observations (El Chichón). These are all in the satellite era and therefore assimilate more atmospheric data that may lead to a more complex interaction for the surface response. The phase of the NAO leads to a dynamically driven warm anomaly over northern Europe in winter, which is present in all datasets considered. The general cooling of the surface temperature due

Abstract. In light of the range in presently available observational, reanalysis and model data, we revisit the surface climate response to large tropical volcanic eruptions from the end of the 19th century until present. We focus on the dynamically driven response of the North Atlantic Oscillation (NAO) and the radiative-driven tropical temperature response. Using 10 different reanalysis products and the Hadley Centre Sea Level Pressure observational dataset (HadSLP2) we confirm a positive tendency in the phase of the NAO during boreal winters following large volcanic eruptions, although we conclude that it is not as clear cut as the current literature suggests. While different reanalyses agree well on the sign of the surface volcanic NAO response for individual volcanoes, the spread in the response is often large ( ∼ 1/2 standard deviation). This inter-reanalysis spread is actually larger for the more recent volcanic eruptions, and in one case does not encompass observations (El Chichón). These are all in the satellite era and therefore assimilate more atmospheric data that may lead to a more complex interaction for the surface response. The phase of the NAO leads to a dynamically driven warm anomaly over northern Europe in winter, which is present in all datasets considered. The general cooling of the surface temperature due to reduced incoming shortwave radiation is therefore disturbed by dynamical impacts. In the tropics, where less dynamically driven influences are present, we confirm a predominant cooling after most but not all eruptions. All datasets agree well on the strength of the tropical response, with the observed and reanalysis response being statistically significant but the modelled response not being significant due to the high variability across models.

Introduction
Understanding of the atmospheric naturally forced variability is a key issue in estimating the human-induced contribution to the recent climate change. Large volcanic eruptions can have an impact on the climate over several years (Robock, 2000). The injected aerosols into the lower stratosphere influences the radiation balance, resulting in a cooling of the tropical surface (Humphreys, 1913(Humphreys, , 1940 and a heating in the tropical lower stratosphere (Labitzke and McCormick, 1982;Parker and Brownscombe, 1983). For the eruption of Mt Pinatubo in June 1991, which was the strongest tropical eruption in the satellite era, the lower tropical stratosphere was warmed by up to 3 K (Mitchell et al., 2014;Fujiwara et al., 2015). The cooling signal at the surface is less pronounced and therefore more difficult to separate from other internal and external climate variability. To perform this separation Mitchell et al. (2014) and Fujiwara et al. (2015) used multiple linear regression with nine reanalysis datasets (ERA-Interim, ERA-40, JRA-25, JRA-55, MERRA, NCEP-R1, NCEP-CFSR, NCEP-R2 and NOAA-20CR). Good agreement between most reanalysis datasets for the temperature response in the stratosphere was found. Uncertainties which arise from variability of the differences between the reanalysis models are currently being assessed in the Stratosphere-troposphere Processes And their Role in Climate (SPARC) Reanalysis Intercomparison Project (S-RIP) (Fujiwara and Jackson, 2013;Fujiwara et al., 2016). The agreement between the reanalysis datasets in the troposphere is strong but no clear tropospheric cooling could be found when all large tropical eruptions after 1960 were taken into account (Fujiwara et al., 2015). A factor which contributes to the uncertainty of the tropospheric cooling is a winter warming signal over Europe. The stratospheric polar vortex in boreal winter is strengthened due to weaker wave flux from the troposphere into the stratosphere and ozone losses in the high latitudes of the stratosphere (Stenchikov et al., 2002). This strengthening of the polar vortex can lead to descending wind anomalies into the troposphere and intensified westerly winds near the surface (Baldwin and Dunkerton, 1999). With a strengthening of the polar vortex comes a shift of the North Atlantic Oscillation (NAO) towards a positive state, leading to anomalous advection of warm air from the Atlantic to Europe (Thompson and Wallace, 2001). Christiansen (2008) found a significant shift of the NAO towards a positive state, using 13 major volcanic eruptions. The positive NAO causes a winter warming in northern Europe of up to 2 K (Fischer et al., 2007). On a global scale the dynamically driven warming weakens the radiative-driven cooling after large volcanic eruptions; hence the interplay between the two is important to consider.
Unfortunately many CMIP5 models fail to reproduce the dynamical mechanism resulting in a positive NAO in winter (Driscoll et al., 2012); it is currently debated as to why this may be, but possible explanations could be poorly represented spatial variability of volcanic aerosol in the models, in terms of vertical extend or latitudinal banding, a misrepresentation of the dynamics of volcanic aerosol in the stratosphere or sampling issues (only Mt Pinatubo has been well observed at those levels). The tropical surface response to volcanic eruptions comes from reduced shortwave solar radiation reaching the surface, which is instead absorbed by volcanic aerosol in the lower stratosphere (Robock, 2000). A recent detection and attribution study showed that this lower stratospheric response to volcanic eruptions is magnified in models compared with observations (Mitchell, 2016), which agrees well with the larger surface cooling also seen in models (e.g. Stenchikov et al., 2006). Other detection and attribution studies show that the temperature response to natural forcings (combined solar and volcanic forcing) can be detected at the surface (e.g. Jones et al., 2013;Bindoff et al., 2013) in models, although the global cooling signal is overestimated (Ribes and Terray, 2013). A more recent study has shown that El Niño correlations are partially to blame for this signal (Lehner et al., 2016), but there are still uncertainties in the dynamically and radiative-driven response to volcanic eruptions.
While studies have looked at the NAO and tropical surface responses individually to volcanic eruptions before, here we investigate both in conjunction, assessing the variability of the NAO and the general temperature response following volcanic eruptions. We test the robustness of previous studies by including observations and all available reanalysis products to examine the uncertainty of the radiative and dynamical response to volcanic eruptions. Our results are also informed by the latest coupled climate models. The analysis of the uncertainty in the reanalysis products will contribute to the S-RIP report, comparing the representation of various atmospheric modes in reanalysis.
The remainder of this paper is as follows. Datasets and methods used in this study are described in Sect. 2, followed by the results in Sect. 3 and the discussion and summary in Sect. 4.

Data and analysis method
We use near-surface monthly mean temperature of air at the surface (TAS) and sea level pressure (SLP) data from 10 available reanalyses (Table 1). Furthermore we analyse the TAS fields of the CMIP5 historical model experiments provided by the World Climate Research Programme (Taylor et al., 2012).
For comparison of the temperature we use the Met Office Hadley Centre HadCRUT4 dataset, available with 100 ensemble members covering the period 1850 until present (Morice et al., 2012). This product has a global coverage with a resolution of 5 • × 5 • , but it includes missing data depending on the raw observational data used. The provided data are anomalies referring to the mean climatology between 1961 and 1990 based on the relatively small number of missing data during this period.
To investigate the dynamical impact of the volcanic eruptions during boreal winter we make use of sea level pressure observations from the Hadley Centre Sea Level Pressure dataset (HadSLP2) (Allan and Ansell, 2006). This product does not contain missing data because of an applied interpolation procedure, which can cause uncertainties, especially in less covered regions like the Arctic, Antarctic or deserts. The spatial resolution of 5 • × 5 • is equal to the temperature product but the data include the period 1850 until 2004. However, an updated version is available using National Centers for Environmental Prediction (NCEP)/NCAR reanalysis fields from 2005 until present, named HadSLP2r (Kalnay et al., 1996). The mean values for both dataset are homogeneous but the variance is higher in HadSLP2r. Nevertheless we consider this adding of the dataset as justified since we use the updated period only for (a) the calculation of the climatology, (b) significance testing and (c) the empirical orthogonal function (EOF) analysis in order to calculate the NAO time series (Thompson and Wallace, 1998).

Reanalysis data
All reanalysis datasets span at least the period from 1979 until 2012 except two ECMWF products: ERA-20C (20th Century Reanalysis Product) (Poli et al., 2013), which ends in December 2010, and ERA-40 (Uppala et al., 2005), which ends in August 2002. To be able to compare all reanalysis over the same temporal region, we extend the ERA-20C and ERA-40 datasets until 2012 with data of the ERA-Interim reanalysis. ERA-Interim provides some advantages to the ERA-40 version including an improved data assimilation and a better representation of the stratospheric circulation (Dee et al., 2011). Since we only use these added data indirectly for calculations of climatology, anomaly fields and the EOF, the differences which would arise by using a full period consideration are negligible.
Both reanalysis products of the Japan Meteorological Agency (JMA) are used for the study. JRA-25 ends in 2004 but the data until January 2014 are available from the JMA Climate Data Assimilation System (JCDAS) with the same system as JRA-25 (Onogi et al., 2007). The subsequent JMA product is called JRA-55 and covers a longer period beginning from 1958. Several improvements have made in comparison to the previous product such as a significant reduction of the large temperature bias in the lower stratosphere by using a new radiation scheme (Kobayashi et al., 2015).
The MERRA reanalysis obtained from the National Aeronautics and Space Administration (NASA) is focused on the correct simulation of the hydrological cycle and is the only reanalysis used which does not represent the analysis field with spectral coefficients (Rienecker et al., 2011).
The first major reanalysis project was operated by the NCEP, called NCEP-R1 (Kalnay et al., 1996), and has been updated with the NCEP-R2 product (Kanamitsu et al., 2002). A more sophisticated and recent reanalysis product of NCEP is the Climate Forecast System Reanalysis (NCEP-CFSR) (Saha et al., 2010).
Most products include a period beginning from the late 20th century, which means that just the last two large tropical eruptions of Mt Pinatubo and El Chichón can be considered for the analysis. Due to the large number of individual ensemble members and the long time series beginning in 1851, the NOAA-CIRES 20th Century Reanalysis version 2c (NOAA-20CR) is the main product used in this study (Compo et al., 2011). The assimilation scheme of the NOAA-20CR reanalysis product uses an ensemble Kalman filter in streams of 5 years. Each stream has 56 members.
The ERA-20C product starts in 1900 and includes most of the period which we investigated. The assimilation of these long datasets only includes surface observation data, in contrast to the other products assimilating also satellite and radiosonde measurements. The reanalysis datasets are generally in good agreement with surface observation data, especially for sea level pressure and near-surface temperature data, used in this study (Simmons et al., 2004;Makshtas et al., 2007;Lindsay et al., 2014). NCEP-R1 and JRA-25 show differences in the sea level pressure field over Greenland and MERRA generally over mountain areas likely because of different surface pressure reduction methods (Lindsay et al., 2014). Since we use just anomaly fields for our calculations, this should not affect the results significantly.

CMIP5 models
The model data are obtained from the historic simulations of the CMIP5 models, which include simulations with just volcanic forcing (Table 2). All models except the coupledphysical model of the Geophysical Fluid Dynamics Laboratory (GFDL-CM3) start from the beginning of 1850 and have at least three ensemble members. The advantages of the GFDL-CM3 model is a sophisticated interaction scheme between aerosols and clouds and a focus on coupling between the troposphere and stratosphere (Donner et al., 2011). The Community Climate System Model 4 (CCSM4) and the Community Earth System Model version 1 with Community Atmospheric Model version 5 (CESM1-CAM5) models use the aerosol optical depth (AOD) description of Ammann et al. (2007). All other models use the updated version of the Sato et al. (1993) description. CESM1-CAM5 includes the direct and indirect effects of aerosols (Meehl et al., 2013), while CCSM4 just provides the direct effects. These models show a good reproduction of the El Niño-Southern Oscillation (ENSO) due to an improved deep convection scheme in the atmosphere component (Gent et al., 2011).  Sato et al. (1993). c Stenchikov et al. (1998).
The ModelE2 version of the NASA Goddard Institute for Space Studies (GISS-E2) provides four different simulations with just volcanic forcing (Schmidt et al., 2014). They differ by using distinguished ocean models and whether the models include interactive chemistry and parametrization of indirect aerosol effects. GISS-E2-R uses the Russell ocean model (Hansen et al., 2007) and GISS-E2-H uses the Hybrid Coordinate Ocean Model (HYCOM) (Sun and Bleck, 2006). Both realizations are available in a version with non-interactive chemistry (NINT), comparable to the prior CMIP3 simulation, but with a tuned aerosol indirect effect following Hansen et al. (2005) and a version with tracers of chemistry, aerosols and their direct and indirect effects (TCADI) including interactive chemistry and a parametrization of the first indirect aerosol effects (Menon et al., 2010).

Choice of volcanoes
Not all studies which concentrate on the large-scale impact of volcanic eruptions use the same criteria for the choice of which eruption should be considered for a composite analysis. The volcanic explosivity index (VEI) introduced by Newhall and Self (1982) is a very frequently used measurement for the strength of the eruption (Robock, 2000). The calculation of the index is restricted to volcanic measurements but does account for the height of the eruption column. In Fig. 1 all volcanic eruptions since the 19th century with a VEI of at least 5 are shown. Additionally, volcanoes are included which in other studies are considered to have an impact on the climate but just reached a VEI of 4. The size of the triangle of each eruption indicates the respective VEI. The colour shows the phase of ENSO in the first winter after the eruption. The last three tropical eruptions of Mt Agung in 1963, El Chichón in 1982 and Mt Pinatubo in 1991 were followed by an El Niño event, suggesting a possible connection between large tropical volcanic eruptions and ENSO (Adams et al., 2003). Since most studies assume that the events are coincidental, we remove the effects of ENSO with a linear regression (see Sect. 2.4).
Climate models represent volcanic eruption by an increase of the atmospheric aerosols due to the ejected material. Most models use the updated version of the so-called Sato Index (Sato et al., 1993). This index shows the AOD at wavelength 550 nm and is available as a zonal mean with global coverage and a meridional resolution of around 8 • . In Fig. 2 the tropical (30 • S-30 • N) AOD is plotted. As expected from the chosen region the values of low latitude eruptions are generally more pronounced than the extratropical eruptions. The tropical region is characterized by rising air in the stratosphere which lifts the aerosols into higher levels. The residual stratospheric meridional circulation transports the aerosols to high latitudes (Trepte and Hitchman, 1992). A volcanic eruption in higher latitudes is expected to have less influence on the climate system because the downward flow in the stratospheric extratropics avoids rising aerosols in higher levels. Nevertheless some studies show that also extratropical volcanic eruptions can have a significant large-scale impact on climate, but usually just on the hemisphere where the eruption took place (Graf and Timmreck, 2001;Oman et al., 2005). Since we focus on both the particular impact of the eruption on the NAO and the global temperature response, we consider the AOD in the tropical middle and upper stratosphere. In Fig. 2 we choose a threshold at 0.05 nm in order to distinguish which volcanoes likely have a strong impact on the global climate. Therefore we choose the eruption of Krakatau in August 1883, Santa María in October 1902, Mt Agung in March 1963, El Chichón in April 1982 and Mt Pinatubo in June 1991 for our analysis.

Methods
For the calculation of the NAO we use the leading EOF between 0-70 • W and 35-80 • N during the period 1979-2012 (Thompson and Wallace, 1998;Baldwin et al., 2008). We consider the SLP anomalies of the two winters following the eruption as volcanically influenced, except in the case of Santa María, which erupted so late during the year that the full influence on the winter circulation is unlikely. The  anomalies are calculated with respect to the mean for the years 1979-2012, excluding the following 2 years after the eruptions of El Chichón and Mt Pinatubo. The analysis of the temperature fields will be compared to the HadCRUT4 dataset, consisting of anomalies relative to the 1961-1990 reference period (Morice et al., 2012). Since the HadCRUT4 dataset contains missing data, we just use the grid points and time steps of the reanalysis products with non-missing values in the observational record. This ensures that all datasets can be compared directly. As some reanalysis products start in 1979 we only consider the five reanalysis datasets which include this reference period (ERA-20C, ERA-40, JRA-55, NCEP-R1 and NOAA-20CR).
After removing the mean seasonal cycle we subtract the data with a 10-year running mean to remove any further trend. To make sure that the running mean is not influenced by the considered year, we average over the 60 months before the year and the 60 months after the considered year.
Both the tropical temperature and NAO responses can be influenced significantly by the occurrence of ENSO (Lehner et al., 2016). To remove the effects of ENSO we calculate a monthly linear regression on the Niño 3.4 sea surface temperature index. This seasonally dependent ENSO pattern is then subtracted from the data. The significance of the resulting anomaly field is calculated with a Monte Carlo test assuming independence between the volcanic eruption events.

Pressure and NAO response
To analyse the NAO response to large volcanic eruptions we use surface pressure data of all available reanalysis products. Since MERRA, ERA-Interim, JRA25, NCEP-R2 and NCEP-CFSR start in 1979 with the beginning of continuous satellite F. Wunderlich and D. M. Mitchell: The observed response to volcanic eruptions observations, we use the years 1979-2012 as our reference period. Two tropical eruptions had a significant impact on the climate system during this period: El Chichón in 1982 and Mt Pinatubo in 1991. Figure 3 shows the mean SLP anomaly in the first two post-volcanic winters after El Chichón and Mt Pinatubo over the extratropical Northern Hemisphere in observations and for the multi-reanalysis mean.
The response pattern is captured well in all reanalysis products. Since the assimilation of surface pressure data is essential for reanalysis products, the difference between the individual products and the observations is expected to be small (Kalnay et al., 1996). Over the Arctic region a higher level of observational uncertainty is apparent due to the decreased number of assimilated measurements than in the midlatitudes. Low-pressure anomalies over the Arctic region are observed, bordered by a positive signal with the centre over Europe. Differences between the reanalysis products are weak but can be seen mainly above mountain areas, e.g. Rocky Mountains (Supplement, Figs. S1, S2), likely due to different pressure reduction techniques (Lindsay et al., 2014). Also the NOAA-20CR reanalysis product is in good agreement with other reanalysis datasets. NOAA-20CR is found to capture well the stratospheric temperature response to volcanic eruptions but with a slightly lower amplitude in comparison to other reanalysis products (Mitchell et al., 2014;Fujiwara et al., 2015). At the surface the NOAA-20CR does not show major differences in the pressure response after the eruption of Mt Pinatubo and El Chichón (Figs. S1, S2).
In comparison to an expected positive NAO, in the first winter the high SLP anomaly is shifted towards central Europe and therefore the negative centre is shifted northwards (Hurrell and Deser, 2009;Hurrell, 1995). Fujiwara et al. (2015) showed that the eruption of Mt Pinatubo influenced the stratospheric temperature and circulation stronger than the eruption of El Chichón. Therefore it is expected that the surface pressure signal is dominated by the response to the eruption of Mt Pinatubo. The strong Mt Pinatubo response was from the second winter (Fig. 3b, d), which shows an NAO-like pattern over the North Atlantic with significant positive anomalies in the region of the Azores, whereas the response from the first winter was diluted from interactions with the quasi-biennial oscillations (QBO) (Stenchikov et al., 2002). A sample of only two single events is not a robust data basis to make general conclusions. Hence, we use the NOAA-20CR product which agrees well with the other reanalysis products, considering the SLP response after the eruption of Mt Pinatubo and El Chichón (Figs. S1, S2), and contains a longer period between 1851 and 2012. During this period five major tropical eruptions took place: Krakatau, Santa María, Mt Agung, El Chichón and Mt Pinatubo. The mean SLP response to these five eruptions is shown in Fig. 4. The pattern in the North Atlantic region is similar to the mean over the last two eruptions but the pressure anomaly gradient is smaller, with a predominantly significant response pattern   (Table 1).
in the first winter. In general, the average pressure anomaly field indicates a shift towards a positive NAO in the first winter after the eruption. A significant mean NAO signal could be found by Christiansen (2008) and Driscoll et al. (2012), taking into account tropical and extratropical eruptions. We calculated the monthly winter NAO index of just the tropical eruptions and show the index individually for every eruption in Fig. 5. The observations are shown by crosses and all reanalysis data are shown by blue ranges. The orange ranges indicate the NAO index of the 56 NOAA-20CR ensemble members (95 % ensemble spread). The spread of the reanalysis products is small, showing high agreement between the individual reanalysis products. The observations agree well with the reanalysis data but in some cases the reanalysis NAO response can differ from the observational signal (e.g. in December and January after the eruption of El Chichón). Especially with higher uncertainty  of the observation data in the pre-satellite era, the differences of the NAO index between the individual members of NOAA-20CR are increased, indicated by a wider spread of the error bar. Therefore some ensemble members can differ with the corresponding observations by strength and even sign of the NAO phase (e.g. Santa María in December). Consistent with Jones et al. (2003) there is no evidence for a positive NAO shift due to the volcanic eruptions in particular months and in the second winter after the eruption (Fig. S3). The response in the mean winter NAO is not clear. In the first winter after the eruptions of El Chichón and Mt Pinatubo a strong positive NAO phase could be found. Relative to the distribution shown with grey bars, a moderate positive NAO was present following the eruption of Krakatau. For the eruption of Santa María we do not use the first winter after the eruption because the volcano erupted so late during the year that an influence of the injected aerosols on the stratospheric circulation is unlikely. In the winter 1903/04 an NAO index around zero was found. This agrees with the results of Christiansen (2008) and Driscoll et al. (2012). In the winter directly after the eruption (1902/03) a strong positive NAO was found (Christiansen, 2008). In the winter after the eruption of Mt Agung we found a negative NAO. Most of the aerosols after the eruption of Mt Agung were concentrated in the Southern Hemisphere (Fig. 2), which reduces the impact on the bo-real stratosphere in winter. Therefore we conclude that we do not find a significant positive NAO response to volcanic eruptions with taking just the strongest five tropical eruptions from the end of the 19th century until present. We confirm that the NAO generally shifts towards a positive state in Figs. 4 and 5d but exceptions like Mt Agung or Fernandina are found (Fig. S4). Extreme stratospheric variability linked with waves sources not associated with the volcanic eruptions might well be responsible for these exceptions as these can influence the NAO over a short timescale of weeks (Baldwin and Dunkerton, 1999). The positive NAO response to volcanic eruptions leads to a positive temperature anomaly over northern Europe in winter (Robock andMao, 1992, 1995;Fischer et al., 2007). This warming disagrees with the radiative-driven cooling of the troposphere following the eruptions. Therefore just dynamically driven effects like circulation changes could explain this response. It is important to understand these dynamically driven effects in order to understand the total volcanic signal. The studies of Stenchikov et al. (2006) and Driscoll et al. (2012) show that general circulation models are generally not able to reproduce these secondary effects and hence it is questionable if they reproduce the temperature response well. Only a subset of models show an associated warming over northern Eurasia but much weaker than the observations (Driscoll et al., 2012;Gillett and Fyfe, 2013). Figure 6 shows the TAS anomaly composites over the last three large volcanic eruptions of the reanalysis products in comparison to the observation data. The temperature anomalies are calculated over the whole first or second year after the eruption. We note here that, due to the observational anomalies respective to the period 1961-1990, we only consider the five reanalysis products which include this period. Influences due to ENSO and temperature trends are removed. Removing the ENSO effect is important due to the sampling of El Niño events after large volcanic eruptions. Without this step there is a warm bias in the tropics, especially over the eastern Pacific Ocean (Fig. S5). This warming counteracts the general tropical cooling signal after large volcanic eruptions and therefore would weaken the radiative-driven temperature signal. Only small temperature effects can be attributed to ENSO in the extratropics (Fig. S5). Another advantage of removing the ENSO signal is the reduction of the overall temperature variance. This allows for a smaller volcanic signal to be separated from other internal variability sources. We find a significant warming over northern Europe which is consistent with the expected winter warming over this region. The warming over Siberia and the cooling in the Middle East suggests a positive NAO response following volcanic eruptions (Thompson and Wallace, 1998). The cooling over Alaska and Northern Canada and the warming over the Southeastern United States and North Pacific suggest a negative Pa-F. Wunderlich and D. M. Mitchell: The observed response to volcanic eruptions (a) December Mt. Pinatubo, 1991El Chichon, 1982Mt. Agung, 1963Santa Maria, 1902 Krakatau, 1883 -4 -3 -2 -1 0 1 2 3 4 0 (b) January Mt. Pinatubo, 1991El Chichon, 1982 Mt    El Chichón (1982) andMt Agung (1963). Panels (a) and (b) show the mean anomaly for the HadCRUT4 observations and panels (c) and (d) show the mean anomaly of the five reanalysis datasets, containing this period (Table 1). Anomalies are calculated with respect to the average over the years 1961-1990, consistent with the Had-CRUT4 dataset. Single diagonal lines correspond to the 90 % and double diagonal lines to the 95 % confidence level obtained by a Monte Carlo test of three independent samples. cific/North American teleconnection pattern (PNA). This is consistent with a stronger polar vortex in boreal winter due to less wave propagation into the stratosphere (Garfinkel and Hartmann, 2008). The reanalysis mean shows a general cooling over sea areas, which is significant over the tropical Pacific and Southern Indian Ocean. The HadCRUT4 dataset shows a similar temperature anomaly pattern in the high northern latitudes. The cooling over the oceans is less pronounced in the observation dataset. The second year after the eruptions shows cooling over most continental and ocean areas with exceptions over Northern Siberia and the eastern Pacific.

Temperature response
To assess a larger sample of eruption events we expand the considered period and include the early eruptions of Krakatau and Santa María. During the period after this eruptions much less observation data are available, especially over continental areas. In Fig. 7a and b we average over at least 4 of 5 years after the volcanic eruptions. The results are similar to those of Fig. 6. To be able to evaluate the effect of the volcanic eruptions over less observed areas we calculated the mean temperature response for all five eruptions with the NOAA-20CR product without missing data consideration (Fig. 7c, d). Also the consideration of five independent eruptions confirms the significant warming over northern Europe due to the shift of the NAO towards a positive state. The strongest cooling signals are found in areas with sparse observations. Therefore, these signals highly depend on the reanalysis configuration. In the tropics a general cooling is found over land and sea areas with the exception of the eastern Pacific in the second year after the eruptions. This would confirm that the cooling after large volcanic eruptions is strongest in the tropics (Driscoll et al., 2012). For this reason we focus on the TAS signal in the tropical region (30 • N-30 • S). Figure 8 shows the mean anomalies of the tropical temperature in the first and second year after the eruptions, relative to the climatology of . The reanalysis and observational data are generally in good agreement, except in the first year after the eruption of Krakatau and Mt Pinatubo, with a separation of the reanalysis anomaly spread and the mean temperature response in the observational record. This shows that the spread of the reanalysis products does not account for the whole range of uncertainty present in the observations. A general cooling of the tropical TAS response is not always visible in the observational and reanalysis data. After the eruption of Mt Agung negative and positive anomalies in the first winter after the eruption are found. The second year after the eruption of El Chichón was characterized by a warming of the tropics. The study of Fujiwara et al. (2015) did not show a clear temperature signal in the tropical troposphere after the eruption of El Chichón. The eruption of Mt Pinatubo was the strongest in the 20th century and caused a cooling around 0.1 to 0.2 K, according to the observation and reanalysis data. In the first year after the eruption of Krakatau the cooling was just weak, regarding the TAS anomaly distribution, shown with grey bars. We generally find a cooling of the tropical temperature following large volcanic eruptions. The global mean response is less robust due to the impact of dynamical warming especially over Europe and the high uncertainty range of the NOAA-20C dataset after the eruption of Krakatau and Santa María (Fig. S7).
The historical simulations of the CMIP5 models involve a far higher range of TAS responses than the observations. Most simulations with just volcanic forcing show a negative temperature signal. However, due to the large spread of simulated TAS response just the eruption of Mt Pinatubo shows a significant tropical cooling. For consistency we calculated the anomalies relative to the same period as the observation and reanalysis data. This period  was highly influenced by volcanic activity. By calculating anomalies relative to a non-volcanic influenced period (e.g. 1931-1960) we found that all simulations with just volcanic forcing show a negative TAS response after both years of all five eruptions (not shown). The mean TAS response to the eruptions of the historical simulations including also anthropogenic and other natural forcing agrees well with the simulations with just volcanic forcing. Due to the high uncertainty range of the historical and volcanic simulations this result is not valid for all model runs. Most of this uncertainty is due to the difference in the response between the individual models rather than the internal climate variability.

Conclusion and discussion
In this study we investigated the uncertainty in surface climate response to strong volcanic eruptions. The most up-todate available reanalysis products, general circulation models and the newest observation datasets are used to best evaluate the radiative-driven tropical temperature response and the dynamically driven NAO response following eruptions. Given the availability of these new datasets, it is timely to revisit the surface response to test the robustness of past studies of volcanic influences on climate. A summary of volcanic eruption intensity and occurrence of El Niño events is presented in Fig. 1 from 1880 to present.
The shift of the NAO towards a positive state in boreal winter due to an intensification of the polar vortex was noted in some observational studies (e.g. Shindell et al., 2004) but only some models with a good representation of the stratosphere are able to reproduce the associated winter warming over Eurasia (Kirchner et al., 1999). The CMIP5 models generally fail or underestimate the impact of the volcanic eruptions on the northern hemispheric circulation (Driscoll et al., 2012). This shows that the dynamical mechanism is still not fully understood (Graf et al., 2007). Conditional on the injected material into the stratosphere, we selected the strongest five tropical volcanic eruptions from the end of the 19th century until present for our analysis. They are expected to have the biggest impact on the atmosphere.
We confirmed that a positive NAO phase is likely to be present during the first post-volcanic winter, similar to the results of Christiansen (2008) and Driscoll et al. (2012), and that this signal is very persistent throughout the individual months (while not necessarily being particularly extreme in each month). Nevertheless, uncertainties still re-  main because not all winters following large tropical volcanic eruptions show a positive NAO (e.g. Mt Agung and Fernandina). Also none of the particular winter months show a significant shift of the NAO after the eruptions. By taking into account all available reanalysis datasets and the Had-SLP2 observation data we have seen that there is a general agreement between the datasets. The sea level pressure in the North Atlantic varies between the reanalysis datasets, resulting in some uncertainty of the calculated NAO response, and these differences are larger for the more recent volcanic eruptions. The differences in temperature between the reanalysis datasets depend on the considered region. Poorly observed regions such as the poles show higher uncertainty (Fig. S6). The mean tropical response to volcanic eruptions is similar in all reanalysis datasets. While some reanalyses were in better agreement with observations, for the different diagnostics considered, the JRA-55 reanalysis performed consistently well. It is known that the atmospheric condition after the volcanic eruption is a big source of uncertainty for the impact of eruptions on the NAO. The QBO phase (Holton and Tan, 1982;Stenchikov et al., 2004), an El Niño or La Niña event (Moron and Plaut, 2003;Manzini et al., 2006;García-Herrera et al., 2006;Calvo et al., 2009) as well as the solar variability (Lean et al., 1995;Haigh, 2002;Gray et al., 2013) can influence the NAO phase directly or indirectly by modulating the stratospheric winter circulation in the Northern Hemisphere. A more robust indicator for the strengthening of the dynamically driven influence of the volcanic eruption is the characteristic winter warming over northern Europe (Shindell et al., 2004;Fischer et al., 2007). By averaging over the whole first year after the eruption we still could find a significant positive signal in northern Europe due to the winter warming. This means that the general decrease of the surface temperature due to the injected aerosols is overwhelmed by this dynamical effect at midlatitudes. This contributes to the fact that the strength of the radiative-driven cooling is still not completely determined. The less dynamically influenced tropical response to volcanic eruptions shows a general cooling mainly over ocean areas (Stenchikov et al., 2006). Nevertheless, considering the first 2 years after the eruption of Mt Agung and El Chichón no clear tropical cooling in both reanalysis and observation data is found. The model response to the eruptions shows a large uncertainty spread, suggesting a high dependence of the climate response on the atmospheric state during and after the eruption.
The Supplement related to this article is available online at doi:10.5194/acp-17-485-2017-supplement.