Supplement to An empirical model of global climate − Part 1 : A critical evaluation of volcanic cooling

Supplement to An empirical model of global climate − Part 1: A critical evaluation of volcanic cooling T. Canty , N. R. Mascioli , M. D. Smarte , and R. J. Salawitch 1,2,3 1 Department of Atmospheric and Oceanic Science, University of Maryland, College Park, MD 2 Department of Chemistry and Biochemistry, University of Maryland, College Park, MD 3 Earth System Science Interdisciplinary Center, University of Maryland, College Park, MD * Now at Department of Earth and Environmental Sciences, Columbia University, New York, NY ** Now at Division of Chemistry and Chemical Engineering, California Institute of Technology, Pasadena, CA


Introduction
It is well established that both natural and anthropogenic factors influence climate.Natural factors include variations in the intensity of sunlight driven by the ∼11 yr cycle of solar activity, variations in exchange of heat between the atmosphere and Pacific Ocean following the shift in ocean circulation recorded by the El Niño-Southern Oscillation T. Canty et al.: A critical evaluation of volcanic cooling (ENSO), and periodic volcanic eruptions with enough energy to dramatically enhance the stratospheric aerosol burden (e.g., Mass and Portman, 1989;Sects. 2.7 and 3.6.2 of IPCC, 2007;Lean and Rind, 2008).Anthropogenic factors include increases in the radiative forcing at the tropopause (RF) due to rising levels of greenhouse gases (GHGs) that cause the lower atmosphere to warm as well as the industrial release of precursors of tropospheric aerosols that can cause the lower atmosphere to either cool or warm, depending on a myriad of factors (Sect. 2.4 of IPCC, 2007).Throughout, we use RF to refer to the stratospheric-adjusted RF described in Sect.2.2 of IPCC (2007).
Atmosphere ocean general circulation models (GCMs) are commonly used to quantify the relative importance of natural (i.e., volcanoes and ocean circulation) and anthropogenic (i.e., GHGs and tropospheric aerosols) factors on global climate.Soden et al. (2002) used a GCM to show that the 0.5 • C cooling of the global lower troposphere, measured by the Microwave Sounding Unit (MSU) after the eruption of Mt.Pinatubo (hereafter, Pinatubo), was well simulated provided: (a) there is a significant positive feedback due to changes in atmospheric H 2 O, in response to the perturbation of the shortwave (SW) solar and longwave (LW) thermal radiation fields induced by Pinatubo; (b) the MSU data record is adjusted for a ∼0.2 • C warming due to ENSO.The need for a significant climate feedback to quantitatively account for the temperature perturbation following the Pinatubo eruption was also discussed by Hansen et al. (1993), Lacis and Mischenko (1995), Forster and Collins (2004), and Wigley et al. (2005).
Multiple linear regression (MLR) of the global surface temperature anomaly ( T ) has also been used to quantify the relative importance of natural and anthropogenic factors on climate (Lean andRind, 2008, 2009;Foster and Rahmstorf, 2011;Kopp and Lean, 2011;Tung and Zhou, 2013;Zhou and Tung, 2013).Typically, coefficients are found that relate a time series of T to the temporal variation of proxies that represent RF due to total solar irradiance (TSI), volcanoes, ENSO and anthropogenic radiative forcing of climate (AF) due to GHGs and aerosols.Time series of stratospheric optical depth (SOD), available for the past century from one of two independent analyses (Sato et al., 1993;Ammann et al., 2003;see Fig. 2.18 of IPCC, 2007) are used to represent the volcanic term.Lean and Rind (2008) used an MLR analysis to estimate that Pinatubo caused a ∼0.3 • C cooling of global surface temperature, considerably smaller than the canonical 0.5 • C cooling of global mean surface temperature commonly attributed to Pinatubo (e.g., Crutzen, 2006;IPCC, 2007).
Here we conduct an MLR analysis of the global temperature record from 1900 to present.Our model uses as input globally averaged mixing ratios for greenhouse gases from the Representative Concentration Pathways (RCP) database provided for the upcoming IPCC report (Meinshausen et al., 2011;van Vuuren et al., 2011).Mixing ratios of GHGs are es-sentially identical in the four RCP scenarios, for our period of interest, 1900 to present.Abundances from the RCP 8.5 scenario (Riahi et al., 2007(Riahi et al., , 2011) ) are used.Results would be unchanged had a different scenario been chosen; the abundance of GHGs differs very slightly between the scenarios starting in 2005 (CH 4 and N 2 O) and 2008 (CO 2 ).In our companion paper, Mascioli et al. (2012), we examine the sensitivity of future climate to the four RCP scenarios, which diverge strongly for GHGs midway through this century.
The start date of year 1900 for our analysis allows examination of perturbations to global climate following the major eruptions of Santa María (October 1902), Mount Agung (March 1963), El Chichón (April 1982), and Mt.Pinatubo (June 1991) and all minor eruptions strong enough to affect stratospheric optical depth, but precludes examination of perturbations due to eruptions of Tambora in 1815 and Krakatoa in 1883.Data needed for our analysis becomes more scarce and uncertain prior to 1900.Also, some of the key figures in IPCC ( 2007) important for our work, such as Figs.TS.23 and 9.14, begin around 1900.Our time period covers the same set of major volcanic eruptions examined by Wigley et al. (2005) and Thompson et al. (2009).
In addition to the commonly used regressor variables SOD, TSI, Anthropogenic RF and ENSO, we introduce to the regression proxies representing variations in the strength of the Atlantic Meridional Overturning Circulation (AMOC), the Pacific Decadal Oscillation (PDO) and the Indian Ocean Dipole (IOD).The detrended Atlantic Multidecadal Variability (AMV) index is used as a proxy for the strength of the AMOC (Andronova and Schlesinger, 2000;Knight et al., 2005;Stouffer et al., 2006;Medhaug and Furevik, 2011;Srokosz et al., 2012).Throughout, we provide extensive discussion of AMOC and AMV, and little discussion of the PDO and IOD, because we compute small contributions of the PDO and IOD to variations of global temperature.Indeed, this is the basis upon which Schlesinger and Ramankutty (1994) first identified the global, climatic significance of multi-decadal variations of sea surface temperature (SST) in the North Atlantic basin.We use the term AMV (e.g.Semenov et al., 2010) to represent North Atlantic SST, rather than the more commonly used Atlantic Multidecdal Oscillation (AMO), because the purely periodic nature of North Atlantic SST has been questioned (Vincze and Jánosi, 2011) and is of no consequence to our study.Below, we show that the AMV index accounts for a strong component of the variations in the century long temperature record provided by four data centres.
We find that global cooling attributed to volcanoes ( T VOLCANO ) declines sharply, by almost a factor of 2, when the AMV index is introduced into the regression.Much of this paper focuses on the robustness of this result.We show that the precise value of T VOLCANO depends on how SST in the North Atlantic, the basis for the AMV index, is detrended.Details of the various empirical parameters used in this analysis are provided in Sect. 2. A description of the model is given in Sect.3. Results of the regression analysis are provided in Sect. 4. Discussion of these results, including the physics of the AMOC as well as implications for geoengineering of climate, is provided in Sect. 5. A brief conclusion follows.Many abbreviations and symbols are used; although each is defined, a Glossary of terms is provided in Appendix A. Appendix B contains web addresses (URLs) for the many sources of data used in the analysis.Appendix C details how we have arrived at an estimate for the empirical range of net anthropogenic aerosol radiative forcing (NAA RF) from IPCC (2007).

External sources of data
This study uses many sources of data.We describe here data used in the manner provided: i.e., data obtained from external websites without further processing.Data records that require internal processing, such as the index for the AMV and the terms used to define anthropogenic RF of climate, are described in Sect. 3 (Model description).

Global temperature
Our regression model, described in Sect.3, uses as input monthly mean near surface air (hereafter, surface) temperature anomalies (either global or land) ( T OBS i ) and the 1sigma uncertainties of each monthly measurement (σ OBS i ).Some of the papers and/or data files provide 2-sigma uncertainties: if so, we have multiplied these values by 0.5 to obtain an estimate of the 1-sigma measurement uncertainty.Throughout, the use of italics, as for CRU4 below, denotes that a web-link for this data source is provided in Appendix B.

Surface
CRU4: the Climate Research Unit (CRU) of the University of East Anglia together with the Hadley Centre of the UK Met Office provide a global, monthly mean surface temperature (Morice et al., 2012) and land temperature (Jones et al., 2012) record.We use the most recent version of each: Had-CRUT4 for global and CRUTEM4 for land.Below, we refer to both HadCRUT4 and CRUTEM4 as CRU4, with clear notation for global or land.The land record is based on data from 5583 stations.The global record combines this information with a SST record based on ship and buoy observations from the International Comprehensive Ocean-Atmosphere Data Set (ICOADS).
Both HadCRUT4 and CRUTEM4 time series represent anomalies relative to the mean value of T OBS from 1961 to 1990.The uncertainties for HadCRUT4 are provided in data files accessible from the website noted in Appendix B. Uncertainties for CRUTEM4 were obtained from Jones et al. (2012).
The CRU4 record is the only dataset used for our "ladder plots" that compare modelled and measured surface T .We have chosen CRU4 for these plots due to the prominence placed on this record by IPCC (2007).Other datasets for global and land T are represented in summary figures that illustrate the conclusions of this study are similar for all long term climate records.
GISS: the Goddard Institute for Space Studies (GISS) provides a climate record based on SST from a combination of the Hadley Centre analysis (HadISST1) for 1880 to 1982 and satellite observations for 1982 to present (Hansen et al., 2010) and land temperature from over 700 surface meteorological stations that are part of Global Historical Climatology Network-Monthly (GHCN-M).Data from GISS are available from prior to 1900 until the end of 2011.Uncertainty estimates are from Hansen et al. (2006) for global temperature and from Hansen et al. (2010) for land.Values of T from GISS are presented as the anomaly with respect to 1951 to 1980.The use of a different time period for the anomaly, compared to CRU4, affects only the constant term (variable C 0 in Eq. ( 2) below) in the regression and is, therefore, of no significance for the present analysis.
NCDC: the National Climate Data Center also provides global and land T , from prior to 1900 until the end of 2011 (Smith et al., 2008).The land temperature anomalies are based on GHCN-M and SST is from ICOADS.Uncertainty estimates are from Smith et al. (2008) for global T and from Smith and Reynolds (2005) for land.This temperature record is presented as an anomaly relative to mean temperature from 1901 to 2000.We only show T VOLCANO versus SOD for NCDC Land; the T VOLCANO versus SOD relation for NCDC Global is virtually indistinguishable from similar relations using the global record from CRU4 and GISS.
BEG: the Berkeley Earth Group provides an estimate of T over land, from prior to 1800 to May 2010, based on measurements from 39 390 unique meteorological stations (Rohde et al., 2013).The estimate for T from BEG uses data from many more sources than considered by other data centres.GHCN-M has strict criteria for record length, completeness and establishment of a station baseline before data from a particular station becomes part of their record.The BEG has developed a methodology for the use of all data (Rohde et al., 2013).This record is available only for land at the present time.Uncertainties are provided in the BEG data file.The BEG temperature anomaly is relative to 1950 to 1980.(Christy et al., 2000) and the land TLT product provided by Remote Sensing Systems (Appendix B).The MSU anomalies are relative to the mean temperature over the January 1979 to April 2002 period.Uncertainty estimates for MSU are from Christy et al. (2003).The global MSU lower troposphere record shown below agrees well with data in Fig. 2 of Soden et al. (2002).

Regression variables
Here we describe the origin of variables used in the regression model that require no processing.SOD, GISS: we use a monthly mean, globally averaged time series of stratospheric optical depth (SOD), available from 1850 to present from GISS (Sato et al., 1993) as a proxy for the volcanic perturbation to the stratospheric sulfate layer.This dataset is based on ground, balloon-borne, and satellite observations.Satellite observations are available only from late 1978 to present.The scarcity of observing stations early in the record requires assumptions to be made regarding the geographic distribution of volcanic aerosols.We use the GISS record for SOD in the main body of our paper because it is the only SOD record regularly updated.
SOD, NOAA: a time series of monthly mean, globally averaged stratospheric optical depth (SOD) is also available from NOAA (Ammann et al., 2003).This time series is based on a 4-member ensemble simulation of volcanic eruptions, within a GCM that resolves the troposphere and stratosphere.This record is available from 1890 to 2008. Figure 2.18 of IPCC (2007) compares the Sato et al. (1993) and Ammann et al. (2003) records for SOD.Generally, the peak SOD from Ammann exceeds the peak SOD from Sato after major volcanic eruptions.As shown in the Supplement, use of SOD from NOAA rather than SOD from GISS in our regression has no bearing on our finding regarding the sensitivity of T VOLCANO to consideration of the AMV index.TSI: the total solar irradiance (TSI) time series used in our regression model is from the Naval Research Laboratory reconstruction of Lean (2000) and Wang et al. (2005).This dataset is based on measurements from a variety of spaceborne sensors starting in 1978, such as the Solar Stellar Irradiance Comparison Experiment on the Upper Atmosphere Research Satellite.For earlier periods of time, the reconstruction uses information such as number, location, and darkening of sunspots as well as time series of Mg-II and Ca-II Fraunhofer lines recorded by ground based instruments.There has been recent debate over the absolute value of TSI (e.g., Kopp and Lean, 2011) as well as the variation of solar output in the ultraviolet (UV) at different phases of the 11 yr cycle (e.g., DeLand and Cebula, 2012;Lean and De-Land, 2012).Neither affects our study.The 11 yr periodicity of TSI is well established and the timing of the peaks and valleys are known.The regression model results are insensitive to the absolute value of TSI as well as variations of solar irradiance in the UV that have little consequence for TSI.
ENSO: we use the NOAA Multivariate El Niño-Southern Oscillation Index (MEI) of Wolter and Timlin (2011).This index is derived from observations of cloud fraction, sealevel pressure, surface wind, sea surface temperature and surface air temperature.Unlike other ENSO indices, such as those based on surface pressure, this ENSO index is dimensionless.
PDO: the Pacific Decadal Oscillation (PDO) represents the temporal evolution of specific patterns of sea level pressure and temperature of the Pacific Ocean, poleward of 20 • N (Zhang et al., 1997), that have been shown to correlate with salmon, anchovy and sardine populations in the Pacific (Chavez et al., 2003).We use a dimensionless index based on analysis of Empirical Orthogonal Functions of the SST anomaly conducted by the University of Washington (Zhang et al., 1997).The PDO is caused by the response of the ocean to spatially coherent atmospheric forcing (Saravanan and McWilliams, 1998;Wu and Liu, 2003).This may explain why the PDO has little influence on the global climate record: the PDO is a response to local wind patterns, rather than an indicator of major release (or uptake) of oceanic heat at a magnitude important for global climate.
IOD: the Indian Ocean Dipole (IOD) index represents the temperature gradient between the Western and Southeastern Equatorial Indian Ocean.We use an index, with units • C, provided by the Japan Agency for Marine-Earth Science and Technology (Saji et al., 1999).We have decided to show results using the IOD so that all three major ocean basins are represented.There is little effect of the IOD on global climate, probably due to the size of the Indian Ocean as well as the fact oceanic deep water does not originate in the Indian Ocean basin.

Atmospheric radiation
ERBE: we show satellite observations of perturbations to Earth's radiation budget following the eruption of Mt.Pinatubo, at shortwave (solar) wavelengths and longwave (thermal) wavelengths, as measured by Earth Radiation Budget Experiment (ERBE) instruments on three satellites: ERBS, NOAA-9, NOAA-10.The ERBE instrument measures incoming solar radiation, reflected shortwave radiation and outgoing thermal radiation (Wielicki et al., 2002).We use Edition 3, Revision 1 ERBE data provided by the National Aeronautics and Space Administration, Langley Research Center, Atmospheric Science Data Center.

Model overview
Our MLR model builds on the work of Lean and Rind (2008) and Kopp and Lean (2011).However, our approach differs in four important manners.
First, we explicitly represent the observed increase in ocean heat content (OHC) (Domingues et al., 2008;Carton and Santorelli, 2008;Church et al., 2011).This term was neglected in prior MLR studies, including all of the MLR studies cited in Sect. 1.
Second, we quantify the sensitivity of the regression coefficients to uncertainty in NAA RF.Tropospheric aerosols that drive RF of climate can either cool (sulfate, dust, ammonium nitrate, organic carbon) or heat (black carbon and biomass burning) the lower atmosphere.Prior MLR studies, as well as many climate models, examine only a single scenario for NAA RF.Uncertainties in this term are large.We allow the regression to determine the best value of the climate feedback parameter (variable λ described in Sect.3.2) for prescribed values of NAA RF and OHC.Quantification of the sensitivity of MLR climate simulations to uncertainty in NAA RF, even in a highly parameterised fashion, constitutes an important step forward.
Third, we conduct a weighted MLR that accounts for uncertainties of the climate record.Knowledge of T has become better over time.As a result, the output of the MLR model tends to follow the climate record more closely during the latter half of the century, which is the proper interpretation of the climate record upon consideration of time dependent uncertainties.We also compute the 95 % confidence interval for each regression coefficient using the method outlined in Sect.8.4 of von Storch and Zwiers (2002).
Finally, in our regression model, we do not allow for multiple ENSO indices that are offset in time (Kopp and Lean, 2011) or a decade-long delay in the response of T to Anthropogenic RF (Lean and Rind, 2008).While most of our simulations are based on the use of a single multivariate ENSO index, we also show results using an independent estimate of the response of global-mean surface temperature induced by the spatial variability in ENSO (Thompson et al., 2009).The decade-long delay used by Lean and Rind (2008) is a surrogate for ocean heat uptake, which is explicitly represented in our approach.
Section 3.2 describes our regression model.The external forcing of global climate due to variations in the abundance of GHGs, black carbon (direct and indirect effects), TSI, stratospheric aerosols and sulfate (direct and indirect effects) from 1900 to present are in close agreement with the representation of these terms illustrated in Fig. 1a of Hansen et al. (2005) (see the Supplement).

Model details
Our regression model minimises: where T OBS i and T MDL i represent time series of observed and modelled global, monthly mean temperature anomalies, and σ OBS i is the 1-sigma uncertainty associated with each temperature observation.We write T MDL i as: where Bony et al., 2006) and i denotes month.The term λ p , the Planck feedback parameter (or Planck response), represents the response of RF to perturbations in surface temperature in the absence of any feedbacks, for Earth's present-day overall albedo.The numerical value corresponds to an Earth effective temperature of 241 K (Sect. 1.4.4 of McGuffie and Henderson-Sellers, 2005) and is similar to values commonly used in other empirical analyses (e.g., Bony et al., 2006;Forster and Gregory, 2006;Soden and Held, 2006;Murphy et al., 2009).
Model parameter γ represents the sensitivity of climate to feedbacks that occur in response to a perturbation of RF at the tropopause due to GHGs and aerosols: where The relation between RF and T MDL i in Eqs.
(2) and ( 3) is the same model framework used by Bony et al. (2006) andin Sect. 8.6 of IPCC (2007).Numerical values of λ given in our figures are comparable to values of λ "ALL" given in Fig. 8.14 of IPCC (2007) as well as values of x =p λ x described in Appendix A of Bony et al. (2006) (i.e., the sum of all feedbacks other than the Planck response).Even though upon first inspection it may appear that climate feedbacks are not allowed to operate on ocean heat export, based on the treatment of Q OCEAN i in Eq. ( 2), we note that climate feedback is inherent in the model formulation of Sect. 3.2.4.Our model framework is consistent with the representation of Earth's energy balance used in Raper et al. (2002) and Bony et al. (2006).
We work exclusively in a global, monthly mean framework.Values of the regression coefficients (C j, j =0 to 6 ) and the climate sensitivity parameter (λ) are found such that the cost function is minimised, for specified NAA RF i and OHC (via model variable Q OCEAN i , described in Sect.3.2.4)over the time period of consideration.Below, we show results for 1900 to the end of 2011 (ground-based observations of surface temperature) and 1978 to the end of 2011 (satellite observations of lower tropospheric temperature).
The index i − 6 for SOD represents the 6 month delay between volcanic forcing and surface temperature response.This lag is the same as used by Lean and Rind (2008) as well as Foster and Rahmstorf (2011) and agrees with the 6.8±1.5 month lag estimated by Douglass and Knox (2005).Since the volcanic perturbation occurs in the stratosphere and our model is based on stratospheric-adjusted RF (second panel, Fig. 2.2 of IPCC, 2007), this 6 month delay represents the time needed for the stratosphere to respond to a volcanically induced perturbation in sulfate aerosol loading.Various volcanic eruptions could exhibit different lags due to the latitude of the eruption, which could drive hemispheric and/or latitudinal asymmetries in SOD, especially soon after the eruption (Wigley et al., 2005).Our central result is unchanged if we vary the lag time by ±1 month but begins to vary with larger shifts.For simplicity, we assume all eruptions exhibit the same 6 month lag.
We impose a 1 month delay between variations in TSI and the related temperature response.A delay of 1 month yields the largest value of C 2 , the solar irradiance regression coefficient.This is the same delay for the response of surface temperature to variations in solar irradiance used by Lean and Rind (2008) as well as Foster and Rahmstorf (2011); the latter paper discusses the physical basis of this delay.
A 2 month delay is used for the response of T to ENSO.This delay is based on the lag that yields the largest value for the correlation coefficient of MEI versus the time series of T ENSO calculated by Thompson et al. (2009).Here, T ENSO represents the simulated response of global mean surface temperature to ENSO variability, taking into consideration the effective heat capacity of the atmospheric-oceanic mixed layer.Our 2 month delay is shorter than the 3 and 4 month delays used by Foster and Rahmstorf (2011) and Lean and Rind (2008), respectively.Lean and Rind (2008) and Foster and Rahmstorf (2011) obtained these delays based on the lag that yielded maximum values of their ENSO regression coefficients.Their analyses did not consider the AMV as a proxy for the strength of the AMOC, which could affect their estimate of the ENSO delay.For analysis of the 111 yr surface temperature record there is little variation in the value of C 3 , the ENSO regression coefficient, for either a 3 or 4 month delay compared to the value found using a 2 month delay.Interpretation of the shorter MSU temperature record is somewhat sensitive to specified ENSO delay.However, our overall scientific conclusions are robust for any reasonable choice of the delay in the response of T to TSI or ENSO.
The product of each regression coefficient (C j, j =0 to 6 ) and its associated regressor variable represents the contribution of this term to the global, monthly mean temperature anomaly.For instance, a time series of the effect of volcanoes on global temperature, T VOLCANO i is given by C 1 × SOD i−6 .The maximum cooling due to the eruption of Pinatubo, T PINATUBO , is given by the absolute value of C 1 × 0.15, where 0.15 is the maximum value of SOD observed after this eruption.
The origin of the SOD, TSI, ENSO, PDO and IOD terms, which are all based on external sources of data, has been described in Sect. 2. The anthropogenic RF terms are described in Sects.3.2.1 (GHGs) and 3.2.2(Aerosols), the treatment of the AMV index is detailed in Sect.3.2.3, and Ocean Heat Export is the focus of Sect.3.2.4.

Anthropogenic radiative forcing: GHGs
Figure 1 shows time series of direct RF due to GHGs, a model input represented as GHG RF i in Eq. ( 2).Monthly values of this forcing are specified from the RCP 8.5 scenario.We use global, annual mean mixing ratios of CO 2 , CH 4 and N 2 O provided on the RCP website (see Appendix B), which we interpolate to a monthly time grid.We then compute RF relative to year 1750 based on formula given in Table 6.2 of IPCC (2001).These relations were found using a model that allows for adjustment of stratospheric temperature in response to the radiative perturbation (Myhre et al., 1998).
The RF attributed to tropospheric O 3 is obtained directly from a radiative forcing file provided on the RCP Potsdam website.The RF due to tropospheric O 3 given by RCP 8.5 compares reasonably well to the Shindell et al. (2006) value: the RCP values exceed the Shindell values from 1900 to about 1950, and the Shindell values exceed the RCP estimate for the last few decades.We have run some of our simulations using the Shindell et al. (2006) estimate for RF due to O 3 , and the results are essentially identical to those shown here using O 3 RF from RCP 8.5.
The RF due to halocarbons shown in Fig. 1 is the sum of 30 compounds.For each compound, we have used global, yearly mixing ratios from either RCP 8.5 (Lamarque et al., 2011), Table 5A-3 of WMO (2011), or Velders et al. (2009) (updated by G. Velders, personal communication, 2011).Mixing ratios for HFC32, HFC125, HFC134a, HFC143a, HFC152a, HFC245fa, and HFC365mfc are from Velders.Mixing ratios for CFC11, CFC12, CFC113, CCl 4 , HCFC141b, HCFC142b, Halon 1301 and Halon 2402 are from WMO (2011).All other halocarbons are from RCP 8.5.For each halocarbon, the RF term has been found using the formula given in Table 6.2 of IPCC ( 2001) combined with the radiative efficiency tabulation given in Table 2.14 of IPCC (2007).The use of mixing ratios from Velders et al. (2009) and from WMO (2011) has no bearing on the present study because the differences, with respect to RCP values, are quite slight.However, use of HFC mixing ratios from Velders and the aforementioned halocarbons from WMO have a modest bearing on our companion paper, Mascioli et al. (2012), which projects T to 2060.Future mixing ratios of these species are projected by Velders and WMO to be higher than future mixing ratios in any of the RCP scenarios.Also, radiatively active HFC152a was overlooked in the RCP database.

Anthropogenic radiative forcing: tropospheric aerosols
Figures 2, 3 and 4 detail our treatment of radiative forcing due to anthropogenic tropospheric aerosols, the model input represented by NAA RF i in Eq. ( 2).The estimate of NAA RF i is tied to values of direct RF of mineral dust (Dust), ammonium nitrate (NH x ), fossil fuel organic carbon (OC), fossil fuel black carbon (BC), and biomass burning organic and black carbon (biomass) emissions given by RCP Potsdam (parenthetical terms refer to Fig. 3).First, we describe how we obtain direct and total RF for sulfate.Then parameters α COOL and α HEAT used to define NAA RF i are described.
The direct RF for sulfate aerosols (RF SULFATE−DIR ) given by RCP Potsdam exhibits a steady rise until about 1992, followed by a modest decline until about 2000, then a second peak occurring about around 2010.This time series does not follow the temporal evolution of sulfur emissions (S EMISS ) given by either Stern (2006b) or Smith et al. (2011) (Fig. 2a)., 2007).The curve labelled Biomass refers to emissions of OC and BC due to biomass burning, and the curves labelled OC and BC refer to fossil fuel burning emissions of these components.(c) total RF of aerosols that cool (blue, found for α COOL = 2.4), of aerosol that heat (red, found for α HEAT = 2.4), and their difference that defines NAA RF (black curve).The value of NAA RF in year 2005 is marked.It is coincidence that the two scaling parameters both require a value of 2. RCP provides values for S EMISS that are quite similar to those given by Smith et al. (2011), but only on a decadal time scale that does not reflect known, interannual variations (Fig. 2a).
We have, therefore, formed our own estimate for RF SULFATE−DIR , labelled Smith * in Fig. 2b, that is tied to the Smith et al. (2011) estimate of S EMISS and the Stern (2006a) estimate of total RF due to sulfate (RF SULFATE−TOT ).We have scaled RF SULFATE−TOT from Stern (2006a) by the ratio of S EMISS from Smith et al. (2011) divided by S EMISS from Stern (2006b).We scale by this ratio because the emissions from Smith et al. (2011) are an update to those from Stern (2006b).The resulting, scaled curve is multiplied by a constant factor, at all times, such that the value of RF SULFATE−DIR in 2005 equals −0.4 W m −2 , the best estimate of this quantity given in Table 2.12 of IPCC (2007).
The next step in the calculation of NAA RF is to scale RF SULFATE−DIR to RF SULFATE−TOT .Our best estimate for Figure 3 shows time series for total RF of anthropogenic aerosols that cool (Fig. 3a), aerosols that heat (Fig. 3b), and net anthropogenic aerosol RF (Fig. 3c). Figure 3a, b also show components that contribute to the cooling and heating terms, respectively.All time series shown in Fig. 3 are based on direct RF from RCP Potsdam for specific types of aerosols, except for RF SULFATE−DIR (described above).For aerosols that cool, the direct RF terms from RCP are all multiplied by α COOL and these components are summed to arrive at total RF for aerosols that cool (Fig. 3a).The same procedure is used for aerosols that heat, except a different scaling parameter, α HEAT , is used (Fig. 3b).
While use of two scaling parameters may seem overly simplistic, these parameters capture the essence of the temporal variation of NAA RF in a tractable manner.Chapter 2 of IPCC ( 2007) establishes that total RF due to aerosols is much larger than direct RF due to aerosols.Our scaling parameters represent the various aerosol cloud indirect effects that occur in the atmosphere, which have a myriad of names such as the first indirect effect, the second indirect effect, cloud albedo effect, the Twomey effect, the Albrecht effect and the cloud lifetime effect (Fig. 2.10 and Sect. 2.4.1 of IPCC, 2007).Aerosol cooling is dominated by sulfate particles and aerosol heating is dominated by black carbon.Aerosol cloud interactions that occur, for aerosols that cool, will be dominated by the interactions of sulfate and clouds.Interactions that occur for aerosols that heat will be dominated by effects related to black carbon.While it is possible aerosol cloud interactions may have changed over time, for example, due to a change in the height of power plant smokestacks, movement of meteorological fronts relative to point sources, or a change in the ratio of sulfate to nitrate emissions that alters the chemical composition of aerosols that cool, the simplest assumption is that the effect of aerosols on clouds is constant over time.The present state of knowledge regarding aerosol cloud interactions is so uncertain that, as noted above, quantification of uncertainty in NAA RF in this highly parameterised manner constitutes an important step forward.
Figure 3 shows time series of total RF from aerosols that cool, total RF from aerosols that heat, and the net affect (NAA RF i ) for specific values of α COOL and α HEAT .We choose year 2005 as a benchmark for NAA RF i due to the large number of tables and figures in IPCC ( 2007) that quantify RF of anthropogenic aerosols between 1750 (when NAA RF i was essentially zero) and 2005 (e.g., Figs.FAQ 2.1, 2.20 and 2.21 as well as Table 2.12 of IPCC, 2007).The value of NAA RF i at the end of 2005, denoted NAA RF 2005 , is marked on Fig. 3c.The value for α HEAT = 2.4 used in Fig. 3b, which coincidentally is the same value used for α COOL , was chosen such that a value of −1.0 W m −2 for NAA RF 2005 is obtained (Fig. 3c).This matches the IPCC (2007) central value for NAA RF 2005 (Appendix C).
There is a final important detail regarding the scaling parameters α COOL and α HEAT .There are infinitely many combinations of α COOL and α HEAT that yield the same value of NAA RF 2005 , as shown in  2.12 of IPCC ( 2007) (see Appendix C).The lines marked "High Road", "Middle Road", and "Low Road" show three ways that α COOL and α HEAT can be combined to arrive at the same values of NAA RF 2005 .The region of Fig. 4 bounded by the two limits of the empirical range as well as the High and Low Roads represents our estimate of realistic limits for α COOL and α HEAT.
The value for NAA RF 2005 of −1.0 W m −2 (red line, Fig. 4) can result from many combinations of α COOL and α HEAT .However, simulations of climate using Eq. ( 2 Simulations of T MDL are nearly identical.This model behaviour occurs because the RF terms, for aerosols that both cool and heat, are all tied to precursor emissions.Precursor emissions of all aerosol types have generally risen over time, driven by population growth, economic productivity, and technology (e.g., Myhre et al., 2001;Stern, 2006b;Smith et al., 2011).Values of S EMISS peaked in 1980 according to Smith et al. (2011) (Fig. 2b), but the deviation of present day emissions from the peak value is too small to discern whether the "High Road", "Middle Road" or "Low Road" analysis provides a better simulation of climate.
The key factor for simulating T from 1900 to present is the value of NAA RF for the modern epoch (represented as NAA RF 2005 on Fig. 4).When NAA RF 2005 is towards the upper limit of the empirical range (close to −0.4 W m −2 ), small values of the climate feedback parameter are found following minimisation of the Cost Function.When NAA RF 2005 is towards the lower limit of the empirical range (close to −2.2 W m −2 ), larger values of the climate feedback parameter result.Model parameters NAA RF 2005 and λ cantilever in a similar manner within GCMs (Kiehl, 2007).In Sect.4, we show this cantilevering has little effect on our primary result: volcanic cooling term is much more sensitive to whether the AMV index is included in the regression and how this index is detrended than it is to the value of NAA RF 2005 .

Atlantic multidecadal variability
The Atlantic Multidecadal Variability (AMV) index is a measure of the time evolution of SST in the North Atlantic Ocean, generally between the equator and 60 • N. Schlesinger and Ramankutty (1994) first described the relation of the modern climate record to the temporal variations of SST in the North Atlantic (Fig. 5a, top panel).Analysis of output from numerous oceanic GCMs shows that variations in the strength of the Atlantic Meridional Overturning Circulation (AMOC), also called the thermohaline circulation, are reflected by changes in North Atlantic SST (Knight et al., 2005;Stouffer et al., 2006;Zhang et al., 2007;Medhaug and Furevik, 2011;Srokosz et al., 2012).However, the temporal lag and magnitude of the correlation coefficient between North Atlantic SST and AMOC vary considerably between different ocean GCMs, with some showing considerably stronger relations than others (Medhaug and Furevik, 2011;Srokosz et al., 2012).There is considerable debate regarding the physical processes that drive variations in the strength of the AMOC, which is described in Sect. 4.3. Sutton and Hodson (2005) have shown that boreal summer climate in North America and Europe is modulated by the variations in the strength of the thermohaline circulation on multidecadal time scales.
Most pertinent to our analysis is whether variations in the strength of the AMOC exert an influence on global climate.Artificially induced variations in the strength of the thermohaline circulation, within an ensemble of GCMs, have primary influence on the North Atlantic, but there is global propagation of the signal (Stouffer et al., 2006).Dima and Lohmann (2007) show empirical evidence for the influence of variations in the strength of the AMOC in the North Pacific as well as the Atlantic basin.Further evidence for the influence of variations in the strength of the AMOC extending beyond the North Atlantic is provided by Semenov et al. (2010), DelSole et. al. (2011), Wu et al. (2011), Chiang and Friedman (2012), and Srokosz et al. (2012).Recently, Zhou and Tung (2013) and Tung and Zhou (2013) published MLR analyses of the global temperature record using an AMV index similar to the work reported here.The focus of both of these papers was improved quantification of the amount of global warming that can be attributed to anthropogenic activity upon the consideration of an index for variations in the strength of the AMOC.Neither paper discussed the effect of the AMOC on volcanic cooling.
There is also debate whether the strength of the AMOC has changed monotonically over time, due perhaps to rising GHGs (e.g., Box 5.1 of IPCC, 2007;Willis, 2010).Here, we exclusively use an AMV index that is detrended for the time period of the regression.Our sole focus is quantification of the impact of variability in the strength of the AMOC on T VOLCANO .
There are several groups that provide an index for North Atlantic SST.These groups commonly use the abbreviation AMO for this index.We use AMV because it is unimportant whether the SST signal repeats in a purely periodic manner.The NOAA AMV is based on SST measurements in the Atlantic from the equator to 70 • N, detrended using a linear regression (Enfield et al., 2001).The Royal Netherlands Meteorological Institute (KNMI) provides multiple AMV indices: one based on Atlantic SST from the equator to 60 • N detrended using near global SST (60 • S to 60 • N) (Trenberth and Shea, 2006) and another based on Atlantic SST from 25 • N to 60 • N detrended using a regression against global temperature (van Oldenborgh et al., 2009).Guan and Nigam (2009) compute an AMV based on principal component analysis of the SST poleward of 20 • N. Ting et al. (2009) also report an AMV based on principal component analysis, combined with a low pass filter, using output from six ocean GCMs.
We have examined the impact of these AMV indices in our model framework.The most important detail is how each index is detrended.Below, results are presented using three methods for detrending the AMV, performed internally based on SST data from two data centres.
The SST records are from the Hadley Centre (HadSST3) and NOAA (Kaplan Extended SST V2, hereafter Ka-planSST2).In the main paper, we show results using HadSST3.In the Supplement, we show that our overall conclusions are unaffected if KaplanSST2 is used.
Figure 5 shows various representations of the AMV index from HadSST3.The top panel shows area weighted SST in the Atlantic (equator to 60 • N).The bottom three panels show different representations of the AMV, based on how the index has been detrended.Figure 5b shows use of near global SST (60 • S to 60 • N), as suggested by Trenberth and Shea (2006).This results in an AMV index with less interannual variability (in an absolute sense) than the other detrending methods and, most importantly, a larger value in the 1900 to 1930 time period.Figure 5c shows use of a linear regression, as described by Enfield et al. (2001).We believe both of these methods for detrending can be critiqued.If the expression of changes in the strength of the AMOC truly extends beyond the North Atlantic, as suggested by Dima and Lohmann (2007), Zhang and Delworth (2007), Zhang et al. (2007) and Srokosz et al. (2012), then use of near-global SST to detrend the AMV index removes some of the physical signal.As noted by Trenberth and Shea (2006) as well as Ting et al. (2009), use of a linear regression to detrend the AMV index could result in the aliasing of a global warming signal into the resulting index.The intent of the AMV index is to arrive at a proxy for the internal variability of a component of the climate system independent of anthropogenic RF, which is known to have varied in a nonlinear manner over time.
Figure 5d shows a method for detrending the AMV index that uses the anthropogenic radiative forcing (AF) of climate.We import, to the model, the record of North Atlantic SST shown in Fig. 5a.For each iteration of the model, prior to computation of regression coefficients, North Atlantic SST is used to form a detrended AMV index by regression against (1+γ )(GHG RF i +NAA RF i ). Figure 5d shows converged results for GHG RF from RCP 8.5, NAA RF 2005 = −1.0W m −2 and λ = 1.22 W m −2 • C −1 (i.e., model results shown in Fig. 6d).An AMV index computed in this manner accounts for the fact that growth in background SST over time was nonlinear, which provides a more realistic means to infer variations in the strength AMOC from North Atlantic SST than detrending using a linear fit.

Ocean heat export
As atmospheric levels of GHGs rise, the associated RF perturbation leads to an increase in the temperature of the atmosphere as well as the world's oceans (e.g., Raper et al., 2002;Church et al., 2011;Schwartz, 2012).There has been a concerted effort within the oceanographic community to define the heat content of the upper 700 m of the global ocean (Ocean Heat Content, or OHC) from 1950 onwards, based on a variety of oceanic temperature measurements and data assimilation products (e.g., Carton and Santorelli, 2008).The ocean also responds to the GHG-induced RF perturbation at a depth below 700 m.We simulate rise in OHC below 700 m using a scaling parameter, described below.
We have based our analysis on the OHC record of Church et al. (2011).Their data span the 1950 to 2009 time period (bottom panel, Fig. 6a).We have conducted a linear least squares fit of their data to determine that OHC rose by 21.3 × 10 22 J from 1950 to 2009.To represent this rise in ocean heat within our model, we first convert OHC to heat flux (power), termed Ocean Heat Export (OHE), which represents the flow of energy from the atmosphere to the ocean.Using the surface area of the ocean considered in the measurement (3.3 × 10 14 m 2 ) and the time interval of the data record (59 yr or 1.86 × 10 9 s), the OHE has averaged 0.347 ± 0.0221 W m −2 over the observational period.Schwartz (2012) examined the time constant for the upper 700 m of the ocean to respond to an atmospheric RF perturbation within 5 GCMs and reported a median value of 6.3 yr, which we round to 6 yr.
We have designed our model to match observed OHE over the 59 yr period of observation, taking into consideration this 6 yr time lag.We make no attempt to model the ups and downs in the observational record, because the uncertainties in the reconstruction of OHC are the same magnitude as the fluctuations.Furthermore, these fluctuations are often not coherent in time for various OHC estimates (Carton and Santorelli, 2008).We assume a fixed fraction of the anthropogenic radiative forcing of climate (AF) at the tropopause can either flow into the ocean or can remain in the atmosphere, where it heats the land surface and ocean skin.If the transfer of heat from the atmosphere to ocean is proportional to thermal gradients between these components of the climate system and if T is linearly proportional to AF as commonly assumed in models of this nature (see for instance Eq. (3) of Raper et al., 2002), then the approach below provides an accurate representation of the physical process that drives OHE.
The export of heat from the atmosphere to the ocean in Eq. ( 2) is given by: Equations ( 4) and ( 5) represent the export, to the oceans, of a fixed fraction of the anthropogenic RF of climate.The index i − 72 in Eq. ( 4) represents the 6 yr (72 month) lag between an atmospheric RF perturbation and heat export to the upper ocean (Schwartz, 2012).The term OHE in Eq. ( 5) represents the rise in OHC over the period of observation in the upper 700 metres of the world's oceans, which accounts for 70% of the rise in OHC for the entire ocean (Sect.5.2.2.1 of IPCC, 2007; see also Gregory, 2000 andLevitus et al., 2005).The increase in OHC for ocean depths below 700 is accounted for by multiplying OHE by 1/0.7 in the numerator of Eq. ( 5).2011) estimate of OHC.Clearly the measurement of OHC is matched, on average.The representation of heat flow into the ocean as a fixed fraction of the anthropogenic RF perturbation allows us to simulate Q OCEAN from 1900 to present in a physically consistent manner.Figure 6a shows results for a simulation that excludes the AMV, PDO and IOD indices.When these terms are included in the regression, modelled OHC is indistinguishable from the red line on the bottom panel of Fig. 6a.For the converged model results shown in Fig. 6, 30 % of the anthropogenic RF perturbation has gone into the world's oceans (i.e., = 0.30).Gouretski and Reseghetti (2010)  is much larger than the estimate of Church et al. (2011).The value of OHC from Gouretski and Reseghetti (2010) implies 50 % of the anthropogenic RF perturbation has flowed into the oceans.Use of the Gouretski and Reseghetti (2010) value of OHC has no bearing on the results of this paper; λ adjusts to match T and all of our conclusions involving T VOLCANO are insensitive to OHC.On the other hand, future temperature is sensitive to λ, especially if RF due to aerosols declines as expected.The sensitivity of equilibrium climate sensitivity to OHC is quantified in our companion paper (Mascioli et al., 2012).

Global surface temperature
Figure 6 compares the monthly, global temperature anomaly ( T ) reported by CRU4, from 1900 to the end of 2011, with the modelled value of T (top panel of each "ladder plot").Figure 6a shows a regression that includes indices for volcanoes, TSI, humans (GHG RF and NAA RF) and ENSO.The various rungs of the ladder show the product of the regression coefficient and each index.The red line in the top rung of Fig. 6a is the sum of the four lines shown below, plus the constant term (not shown).The ∼1 • C rise in temperature over the 111 yr period is attributed to anthropogenic radiative forcing of climate ("Human" panel, Fig. 6a).All of the simulations shown in Fig. 6 use GHG RF from RCP 8.5 as well as the time series for NAA RF shown in Fig. 3c: i.e., NAA RF 2005 = −1.0W m −2 found using values of α COOL and α HEAT along the Middle Road of Fig. 4. High frequency variations superimposed on the long-term record are primarily attributed to ENSO.Volcanoes account for short-term decreases in T , with Pinatubo associated with a 0.31 • C drop in global mean surface temperature.Variations in total solar irradiance lead to an 11 yr cycle in T : note the y-axis for the solar rung of Fig. 6a has a different scale than the other rungs, to accentuate the solar cycle.
The cooling attributed to Pinatubo shown in Fig. 6a is quite similar to results of Lean and Rind (2008).They reported a maximum ∼0.3 • C cooling due to Pinatubo based on a MLR analysis of an earlier version of the CRU global surface temperature record.Thompson et al. (2009) (thermodynamic model of T from 1900 to 2009) as well as Foster and Rahmstorf (2011) (MLR analysis of T from 1979 to 2010) both conclude that Pinatubo led to a ∼0.3 • C reduction in global temperature.All three of these studies focused on the removal of volcanic cooling from the T record so that the human influence could be better quantified.None of these papers drew attention to the difference between their computed value of T PINATUBO and the 0.5 • C reduction of global surface temperature attributed to Pinatubo by Crutzen (2006), IPCC (2007), and several other studies.
Figure 6b, c and d show comparisons similar to Fig. 6a, except the AMV, PDO and IOD indices have been added to the regression.The comparison of measured and modelled OHC, as noted above, is only shown once because this plot is nearly identical for all simulations.Figure 6b considers the AMV from HadSST3 (abbreviated Had3) detrended using near-global SST (AMV Had3 SST ). Figure 6c shows results for Had3 AMV detrended using a linear regression (AMV Had3 LIN ) and Fig. 6d shows Had3 AMV detrended using anthropogenic radiative forcing of climate (AMV Had3 AF ).The effect of the Pacific Decadal Oscillation and Indian Ocean Dipole on global climate is minimal.These proxies are not discussed further even though for completeness we continue to use the PDO and IOD in other regressions.
The panels labelled Atlantic in Fig. 6 suggest, depending on how the AMV is detrended, variations in the strength of the AMOC may have an impact on global climate comparable in magnitude to ENSO.There is steady improvement in the ability to simulate global surface T reported by CRU4 as the AMV index is added to the model (Fig. 6b versus Fig. 6a), then as the method for detrending the AMV index is changed from SST to LIN (Fig. 6c versus Fig. 6b), and from LIN to AF (Fig. 6d versus Fig. 6c).The top panel of each ladder plot includes the numerical value of reduced chi-squared, defined as: where N FITTING PARAMETERS equals 5 for Fig. 6a (4 regression coefficients plus γ ) and equals 8 for Fig. 6b, c and d (3 additional regression coefficients).T OBS j and T MDL j represent observed and modeled annual mean global temperature anomalies (average of the twelve monthly values for each year), and σ OBS j is the uncertainty of the annual global mean temperature either provided for annual temperature (i.e., CRU4) or as the annual average of monthly uncertainties.The formulation of χ 2 is based on annual temperature anomalies, rather than monthly anomalies, because the monthly observed temperature anomaly displays an autocorrelation that differs significantly from that of modeled temperature whereas the CRU4 annual temperature anomaly displays an autocorrelation that closely resembles modeled temperature (see Supplement).Physically, a value of χ 2 ≤ 2 indicates that the model agrees with the observations, within the measured uncertainty.All four simulations meet this criterion.The drop in the value of χ 2 , from Fig. 6a to d, quantifies the steady improvement in the ability to simulate T as the AMV index is first considered, then as the method for detrending this index is varied.The increase in the amplitude of the green lines in the Atlantic rung, from Fig. 6b to d, demonstrates that increased importance of the AMV index is responsible for the steady decline in χ 2 .
Figure 1 of Lean and Rind (2008) includes notation for the pre-World War I (WWI) and World War II (WWII) time periods, during which global climate has been traditionally difficult to simulate.The AMV Had3 AF model (Fig. 6d) is able to simulate pre-WWI and WWII particularly well, due in part to changes in the SST record during these periods of time compared to the CRU3 record (Thompson et al., 2008), but mainly because of significant cooling during pre-WWI and large heating during WWII attributed to variations in the strength of the AMOC within this regression.The AMV Had3 AF regression simulates pre-WWI cooling and WWII warming of the CRU4 land temperature anomaly particularly well (Fig. 7d), supporting the possibility that variations in the strength of the AMOC have truly exerted influence on global climate and also that the AMV Had3 AF index is a valid proxy for variations in the strength of the AMOC.The AMV Had3 LIN and AMV Had3 AF simulations yield larger values of C 4 than the AMV Had3 SST regression.Larger values are found because for the 1900 to 1930 portion of the record, the AMV index is deemed to be strongly negative when detrended using either LIN or AF (Fig. 5).
As the influence of the variations in the strength of the AMOC on global climate becomes more prominent, the volcanic regression coefficient declines (volcanic rungs, Fig. 6).The maximum cooling attributed to Pinatubo is 0.31 • C without consideration of the AMV index.Pinatubo cooling drops to 0.14 • C for both the AMV Had3 LIN and AMV Had3 AF simulations (numerical values given in Table 1).Global T after the eruption of Pinatubo is simulated just as well in Fig. 6c and d as in Fig. 6a: in Fig. 6c and d, the model attributes a portion of the observed decline in temperature to variations in the strength of the AMOC rather than to SOD.
Further support for the validity of the use of the AMV index, as a proxy for the effect of variations in the strength of the AMOC on climate, are provided by an analysis that examines the CRU4 temperature record in each hemisphere (see the Supplement).When the CRU4 record for each hemisphere is examined, the numerical value of C 4 (AMV coefficient) is much larger for the Northern Hemisphere (NH) than the Southern Hemisphere (SH).Srokosz et al. (2012) state that "a weakened AMOC is typically accompanied by a slight warming of the Southern Hemisphere, though details differ between models".In our case, the value of C 4 is positive for each hemisphere.However, the contribution of variations in the strength of the AMOC to hemispheric T , as reflected by time series of C 4 ×(AMV index), is much stronger for the NH than the SH.Further discussion is given in the Supplement.
A factor of 2 reduction in the cooling attributed to Pinatubo, upon consideration of the AMV index, challenges conventional wisdom.Many prior studies have accounted for ENSO-related influence on temperature for quantification of volcanic cooling.None have considered an AMOC-related influence, which could potentially be comparable in magnitude to the ENSO influence (e.g., panels labelled El Niño and Atlantic in Fig. 6).Potential criticisms of the analysis presented in Fig. 6 are: (1) the model is based on indices that describe forcing of the climate system, rather than the thermodynamic response; (2) if major volcanic eruptions induce a drop in North Atlantic SST, then the AMV index may be responding to both volcanic forcing as well as variations in the strength of the AMOC; if so then a portion of the climate forcing we have attributed to the AMOC may indeed be due to volcanic eruptions; (3) the robustness of the coefficients resulting from a multiple regression, especially if the regressor variables are not orthogonal (i.e., exhibit a nonzero correlation coefficient).These criticisms are addressed in Sects.4.1.1,4.1.2,and 4.1.3.Finally, Fig. 6 demonstrates the steady rise in T over the course of the 111 yr record is due to the increase in RF of climate from the combination of GHGs and anthropogenic aerosols.The near monotonic rise of T over the past four decades in the panels labelled "Human" in Fig. 6 is quite similar in magnitude to observed T .We find a slight variation in the rise of T attributable to human activity over the past 32 yr when the AMV index is considered: The slopes of the orange lines in the Human rung are 0.118 ± 0.0011, 0.114 ± 0.0011, 0.110 ± 0.0011, and 0.122 ± 0.0012 (units of • C per decade) for panels a, b, c, and d, respectively, for 1979 to 2010 (the same time interval used by Zhou and Tung, 2013).The variation of the slopes is much smaller than the factor of 2 reduction in the component of global warming attributed to human activity during the past 32 yr reported by Zhou and Tung (2013) and Tung and Zhou (2013)  We have repeated the analysis shown in Fig. 6d using T ENSO and T VOLCANO from Thompson et al. (2009) rather than the MEI index for ENSO and the time series for SOD in Eq. ( 2).Thompson et al. (2009) used a simple thermodynamic model of the global atmospheric response to anomalous heating due to ENSO and volcanic eruptions.The T ENSO and T VOLCANO terms represent the response to these forcings.The response tends to be broader than the forcing and offset in time, due to the heat capacity of the atmosphere and ocean.When T ENSO and T VOLCANO are used in Eq. ( 2), no additional delay is assumed within the MLR model, since delay is inherent in the calculation of these response functions.Figure 7a is identical to Fig. 6a, except for the use of T ENSO and T VOLCANO .For Fig. 7, the C 1 (volcanic) and C 3 (ENSO) coefficients are dimensionless because the regressor variables, T ENSO and T VOLCANO , have units of • C. The maximum absolute value of T VOLCANO after the Pinatubo eruption from Thompson et al. (2009) is 0.32 • C. Using their response functions in MLR framework leads to T PINATUBO of 0.28 • C (0.879 × 0.32 • C).Our value for T PINATUBO differs slightly from that of Thompson et al. because we have handled radiative forcing due to aerosols in a different manner and we have included ocean heat export.Our value of 0.28 • C and Thompson's value of 0.32 • C for T PINATUBO both suggest Crutzen (2006) and IPCC ( 2007) have over-estimated the amount of global cooling induced by this eruption, even without consideration of possible effects of the AMOC.
Figure 7b is identical to Fig. 6d, except for the use of T ENSO and T VOLCANO .We find a value for T PINATUBO of 0.11 • C using the thermodynamic approach and the AMV index, which is close to the value of 0.14 • C found using indices for ENSO and SOD within the MLR model.Values for C 4 , the AMV regression coefficient, are similar.These comparisons demonstrate that our conclusions are robust and do not result from an artifact due to the use of indices (forcing terms) rather than response functions in our model framework.

Fourier analysis
In this analysis, there is a risk of double counting the volcanic signal; i.e., the improper attribution of volcanic cooling because a major eruption can potentially drive a decline in North Atlantic SST that could be reflected in the AMV index.We have repeated the calculations shown in Fig. 7, using Fourier analysis to remove the high frequency component of the AMV index, to address this concern.The power spectrum of the AMV index is shown in the Supplement.There is a broad peak at periods of 60 to 90 yr per cycle, which may or may not be indicative of purely periodic behaviour.There are numerous other, smaller peaks; the minor peak with highest frequency is at 1/9 yr −1 .
We have constructed the low frequency component of the AMV index by filtering the Fourier transform of the index, removing the high frequency component, then computing an inverse transform.Presumably any volcanic influence on North Atlantic SST would be realised soon after each major eruption, certainly within the first 8 yr.Therefore, any volcanic contribution to the AMV index should be removed by application of a filter that restricts frequencies higher than 1/9 yr −1 .
Figure 7c and d show results of a regression using two filtered reconstructions of the AMV index.For the filtered reconstruction based on removal of frequencies higher than 1/9 yr −1 , a value for T PINATUBO of 0.09 • C is found (Fig. 7c).A frequency cut off of 1/73 yr −1 results in T PINATUBO of 0.16 • C. The resulting AMV index in this case is nearly identical to the multidecadal variability (MDV) of climate attributed to changes in the strength of the AMOC by Wu et al. (2011).If indices for ENSO and SOD are used in the model, rather than T ENSO and T VOLCANO , values of 0.12 and 0.18 • C are found for T PINATUBO for the 1/9 yr −1 and 1/73 yr −1 filters, respectively (supplemental material).We have also conducted regressions for which the AMV index is held constant, at the value observed prior to each of the four major volcanic eruptions, for the duration of elevated T VOLCANO (thermodynamic approach) or enhanced SOD (index approach).In this case, T PINATUBO of 0.19 • C and 0.22 • C are found.All of these values for T PINATUBO are smaller than found when the AMV index is neglected in the regression.These calculations demonstrate that the diminution of T VOLCANO upon introduction of the AMV index is driven by the low frequency, high amplitude component of North Atlantic SST and is not the result of double counting.

Conditional regression
We use conditional regression analysis to quantify the possibility of a misleading result being caused by co-linearity of the regression variables (Denters and van Puijenbroeck, 1989).Using MLR, values of C 1 to C 6 in Eq. ( 2) are found in a single computational step.Significant correlation between regressor variables can potentially lead to misattribution of the cause of a variation in the signal, T MDL i in Eq. ( 2).
Conditional regression provides a means to assess whether the results of an MLR analysis have been compromised due to co-linearity.In this case, the value of one particular regression coefficient is first determined using linear regression.The product of the regression coefficient (i.e., C 1 ) and the regressor variable (i.e., SOD) is formed; then, this product is subtracted from the original signal, resulting in the first residual.The next regression coefficient is found, based on a linear regression of the first residual.The process is repeated until all regression coefficients have been determined.
Results from the conditional regression analysis illustrated in Fig. 8 are based on the use of indices for ENSO and SOD within the MLR model.Similar results are found for the thermodynamic approach (not shown).We have introduced one simplification: only four regression coefficients, C 1 to C 4 (i.e., indices for SOD, TSI, ENSO and AMV) are used.There are 24 possible combinations of the order of these four coefficients.All calculations in Fig. 8 are based on AMV detrended using anthropogenic RF (i.e., analogs of Fig. 6d).
Figure 8 shows results of the conditional regression analysis.The figure displays values of C 4 (AMV coefficient) and T PINATUBO for the 24 conditional regressions (gray points).The figure also shows C 4 and T PINATUBO from the MLR simulation (black points).Error bars represent 95 % confidence intervals, found based on the approach described in von Storch and Zwiers (2002).
The mean values of C 4 and T PINATUBO from the 24 conditional regressions (red horizontal bars; error bars show 2σ standard deviation) are close to the values of these quantities found using MLR.When C 1 is the first regression coefficient found, T PINATUBO equals 0.22 • C (first six gray points).The largest values of T PINATUBO are found when C 3 (ENSO) is the first regressor index and C 1 (SOD) is computed prior to C 4 (AMV).When C 4 is the first computed coefficient, values of T PINATUBO always fall below 0.14 • C. Proper interpretation of Fig. 8 is that, upon consideration of the possible co-linearity of regressor variables, the true value of T PINATUBO lies between 0.05 • C and 0.3 • C, with a most likely value of either 0.14 • C (MLR) or 0.18 • C (mean of conditional regressions).These values of T PINATUBO are all much less than the 0.5 • C global cooling given by Crutzen (2006) and IPCC (2007).

Land surface temperature
Since T over land is completely independent of SST, the determination of a value for C 4 that is significantly larger than zero provides further observational evidence that the variations in the strength of the AMOC do indeed exert an effect on global climate.Figure 9 is identical to Fig. 6, except Fig. 9 is based on analysis of the global land surface temperature anomaly ( T LAND ) from CRU4.Greater volcanic cooling over land is readily apparent.Large, high frequency variations in T LAND are also apparent, resulting in higher values of χ 2 compared to the simulation of global temperature.There is a steady improvement in the ability to simulate T LAND as the AMV is first added to the model (Fig. 9b versus Fig. 9a), then as the method used to detrend the AMV is changed from SST to either LIN (Fig. 9c) or AF (Fig. 9d).A maximum cooling of 0.51 • C over land is associated with Pinatubo when the AMV index is neglected in the regression, consistent with the commonly quoted value for global cooling (ocean and land) following Pinatubo.
The computed cooling due to Pinatubo falls to 0.28 • C and 0.32 • C, respectively, when AMV Had3 LIN or AMV Had3 AF are used in the regression.As apparent in the top panels of Fig. 9c and d 2)) (a) and maximum cooling due to the eruption of Mt.Pinatubo (absolute value of C 1 × 0.15) (b) calculated using the MLR model (large black circles) and from a series of conditional regressions (small gray circles; see text).Black and gray error bars denote the 95 % confidence intervals (von Storch and Zwiers, 2002).
The red horizontal bar denotes the mean of the 24 conditional regressions, and the red error bar denotes the 2σ standard deviation about the mean.For each panel, the first six conditional regression points represent calculations for which C 1 was first computed.The second sextet show results for which C 4 was first computed.For the third sextet, C 3 was first computed and for the last, C 2 was first computed.
time of peak SOD following the eruption of Pinatubo.However, when either AMV Had3 LIN or AMV Had3 AF is used in the regression, a significant portion of this observed cooling is attributed to the AMOC and not Pinatubo.The CRU4 climate record indicates a perturbation to T LAND similar to that reported by Hansen et al. (1993) and Lacis and Mischenko (1995): our interpretation of the cause of a portion of this cooling is all that differs.A literal interpretation of the χ 2 values given in the top panels of Fig. 9 could lead one to conclude that only the simulations shown in panels Fig. 9c and d represent a model consistent with observations, to within the uncertainty of measurement.We are not suggesting such a literal interpretation because values of χ 2 are affected by the high frequency noise of the data record, as well perhaps the invariably subjective nature of specification of measurement uncertainty (i.e., our analysis of the BEG land temperature record never achieves anywhere close to χ 2 = 2 because of the vanishingly small uncertainties associated with this data record; see Table 2).For the CRU4 land record, when the AMV index is neglected (Fig. 9a), T OBS exhibits a sharper increase since the mid-1990s than represented in the model.Simulations that include the AMV index are able to represent T OBS since the mid-1990s better than the simple model, suggesting that a change in the strength of the AMOC may indeed be responsible for part of the recent temperature rise.And, similar to the simulation of global T , variations of T LAND at the time of pre-WWI and WWII are better represented when the AMV index is included in the regression.compared to the simulation of global T , indicative of more rapid warming over land than ocean.

Cause and effect
The reduction of T PINATUBO upon consideration of the AMV index begs the question regarding cause and effect of the variations in the strength of the AMOC.Kuhlbrodt et al. (2007) and Srokosz et al. (2012) have written detailed overviews of the AMOC.There are two distinctly different theories regarding AMOC variability: the ocean salinity/sea ice theory and the atmospheric aerosol/volcano theory.Dima and Lohmann (2007) suggest variations in the strength of the AMOC are caused by the export of sea ice through the Fram Strait, driven by atmosphere-ocean patterns of sea level pressure throughout the Arctic.Temporal variations in sea ice export affect salinity and hence the rate of deep water formation.The freshening of the North Atlantic due to larger flux of sea ice inhibits deep water formation, causing a cooling (negative AMV index).
The reconstruction of Fram Strait sea ice export published by Schmith and Hansen (2003) exhibits temporal variations quite similar to North Atlantic SST (our AMV index).This conceptual picture is supported by the 230 yr long sortable silt grain size record of Boessenkool et al. (2007), invoked as a proxy for the strength of the Iceland-Scotland overflow, that also mirrors our AMV index.Zhang et al. (2007) have suggested the AMOC drives multidecadal variability of temperature throughout the Northern Hemisphere.Delworth and Zeng (2012) present results of a 4000 yr GCM simulation that exhibits internal variability of the AMOC, due to propagation of salinity anomalies, that lead to 0.1 to 0.3 • C hemispheric-scale temperature anomalies (i.e., the same magnitude as the green curves in Figs.6c, d and 9c, d).Meehl et al. (2011) recently presented a five-member ensemble GCM simulation of future climate that shows ∼0.5 • C variations in global temperature they identify as being due internally generated, decadal timescale variability in OHC.
The decadal time scale variability for future global temperature shown in Fig. 1a of Meehl et al. (2011) looks remarkably similar to the decadal time scale variability shown in the Atlantic rungs of our Figs.6c, d and 9c, d.
On the other hand, a number of recent studies have focused on AMOC variability driven by tropospheric aerosols and volcanoes.Church et al. (2005) suggested the eruption of Pinatubo resulted in rapid reductions in both OHC and global mean sea level.Stenchikov et al. (2009) present a GCM simulation that showed volcanic cooling could strengthen the AMOC.Zanchettin et al. (2012) suggested the strengthening of the AMOC would maximise about 10 yr after a major volcanic eruption.Murphy et al. (2009) point to a relation between enhanced volcanic aerosol and a brief cessation of flow of energy into the ocean, based on an 8 yr linear fit smooth-ing of the derivative, with respect to time, of OHC reported by Domingues et al. (2008).Booth et al. (2012) implicate an optical depth anomaly driven by temporal variation of tropospheric aerosols (pollution) and stratospheric sulfate (volcanoes) as a primary driver of SST variability in the North Atlantic.
Linkages between volcanic eruptions, AMOC and OHC are not well established.Iwi et al. (2012) do not find strong evidence for an effect of Pinatubo on the strength of the AMOC and Carton and Santorelli (2008) state that their analysis of OHC does not reflect an impact from the eruption of Pinatubo.As shown in the Supplement, the derivative of OHC from Church et al. (2011), an update to the record of Domingues et al. (2008), bears no relation to SOD (hence, the conclusion of Murphy et al. (2009) seems highly dependent on which OHC record is used, and possibly how the data are smoothed).There is extensive literature on this subject, nearly entirely focused on the debate of whether or not major volcanic eruptions affect OHC and the strength of the AMOC.Our study seems to be the first to suggest that variations in the strength of the AMOC may have compromised prior estimates of volcanic cooling.
Figure 10 examines time series of SOD and the AMV Had3 AF in an attempt to probe cause and effect.Data obtained after perturbations to SOD from the four major eruptions since 1900, Santa María, Agung, El Chichón and Pinatubo, are indicated using specific colours.Figure 10 shows that the AMV index was in a negative phase (i.e., North Atlantic SST tended to be lower than average) at the time SOD was perturbed by the eruptions of Santa María, El Chichón, and Pinatubo.The AMV Had3 AF was neutral at the time of the Agung eruption.As shown in the Supplement, this plot looks similar for AMV Had3 SST and AMV Had3 LIN , except the AMV index is much closer to neutral after the eruption of Santa María when detrended using global SST.
The middle panel of Fig. 10c shows a scatter plot of SOD versus AMV Had3 AF .This is purely a combination of the SOD record as reported by Sato et al. (1993) and the AMV index deduced from HadSST3.Following the eruptions of Santa María, Agung, El Chichón, and Pinatubo, there were 176 months when SOD exceeded 0.01, which we consider to be the volcanic threshold based on visual inspection of the SOD record.For 137 of these months, AMV Had3 AF was negative.The left panel of Fig. 10c shows a similar scatter plot, except here the SOD signal has been moved earlier in time by 6 months.Individual data points move because SOD is now aligned with a different value of AMV Had3 AF .Remarkably, the numerical breakdown is nearly unaltered: 6 months prior to volcanic perturbation of SOD, the AMV index tended to be in a negative phase 80 % of the time.The SOD time series must be shifted backwards in time by about 2 yr for there to be a 50:50 split between negative and positive phases of the AMV index in the scatter plot.If the SOD signal is moved in the opposite direction, later in time, the AMV index tends to be slightly less negative (right panel, Fig. 10c).These scatter plots suggest that variations in North Atlantic SST, as reflected in the AMV index, occurred prior to the volcanic perturbation to SOD.Also, there is no suggestion that when SOD achieves a local maximum (i.e., peak volcanic perturbation), the AMV index is at a local minimum either coincident in time with peak SOD (as suggested by the GCM of Booth et al., 2012) or ∼2 to 3 yr after peak SOD (as suggested by the GCM studies of Stenchikov et al., 2009 andZanchettin et al., 2012) (see the Supplement).
Our empirical examination of AMV from HadSST3 (Fig. 10) and the OHC record (in the Supplement) does not present evidence for a volcanic signature.If a future consensus emerges that major volcanic eruptions truly do impose a significant imprint on the AMV index, such that the indices used in our analysis are flawed (i.e., do not represent variations in the strength of the AMOC) either at the time SOD was elevated or soon (∼6 months) after each eruption, then clearly our primary finding may not be valid.We conclude by noting three crucial points in support of this finding: (1) the regression coefficient for the AMV index is driven by the low frequency, high amplitude component of the signal, which can not possibly be driven by volcanic eruptions; (2) if there is a 2 yr delay in the imprint of SOD on the AMV index, as suggested by Stenchikov et al. (2009) and Zanchet-tin et al. (2012), then our analysis should be valid because the AMV index at the time of SOD perturbation would not be strongly affected by the volcano; (3) if the low frequency temporal variation of the AMV index is truly caused by a tropospheric aerosol optical depth anomaly, as suggested by Booth et al. (2012), our finding is still valid, because this finding is dependent only on the notion that North Atlantic SST variability reflects an external forcing of the climate system with global influence that has been overlooked in prior empirical estimates of volcanic cooling, regardless of physical origin (provided it is not volcanic).

Consideration of other climate records and uncertainty in aerosol RF
Table 1 and Figs.11 and 12 are designed to further assess the robustness of our suggestion that variations in the strength of the AMOC may have compromised prior estimates of volcanic cooling.Table 1 shows values of the maximum cooling attributed to Pinatubo, denoted T PINATUBO , found using all available long-term climate records (Sect.2).Table 2 shows minimum values of χ 2 found from each regression.Numerical values are based on the absolute value of the product of 0.15 (maximum SOD after Pinatubo) and C 1 , the SOD   Kiehl (2007).Figure 11a and b show the cantilevering of λ and NAA RF 2005 described by Kiehl (2007).Most importantly, the relation of these two model parameters is insensitive to how the AMV index is treated.Our model indicates a more compact relation between these terms than found within the GCMs examined by Kiehl (2007), probably because all of our simulations are constrained to match the same value of OHC.The sensitivity of λ to OHC is explored in Mascioli et al. (2012).
Figure 11c and d show the sensitivity of T PINATUBO (same quantity reported in Table 1) to NAA RF 2005 .The computed value of T PINATUBO is moderately sensitive to NAA RF 2005 .As the cooling attributed to tropospheric aerosols drops (as NAA RF 2005 approaches −0.4 W m −2 ), cooling attributed to Pinatubo rises.Figure 11c and d also show the dependence of T PINATUBO on the AMV index.Cooling attributed to Pinatubo depends strongly on whether the AMV index is considered as well as the method used to detrend this index, for all values of NAA RF.Bond et al. (2013) have provided a comprehensive assessment of RF due to BC aerosols, published in preliminary form on 15 January 2013.Their best estimate for the change in global RF of climate due to the direct BC effect is 0.71 W m −2 warming (lower limit of 0.08, upper limit of 1.27), between 1750 and 2005.This estimate is considerably larger than the value of 0.19 W m −2 warming for year 2005 given in the direct RF files provided by RCP Potsdam.Proper interpretation of the Bond et al. (2013) assessment will require quantification of the impact of this new understanding of BC on the total cloud albedo effect, considering all aerosols.Such analysis has not yet been conducted.Once this analysis is complete, it is likely the best estimate for NAA RF 2005 will be a larger numerical value (less cooling) than the baseline value of −1.0 W m −2 used throughout.The calculations for NAA RF 2005 equal to −0.4 W m −2 in Fig. 11 may serve as a useful guide for anticipating the impact, in our model framework, of the Bond et al. (2013) assessment.More cooling is attributed to volcanoes as net anthropogenic aerosol cooling drops, but the importance of the AMOC remains and the cooling attributed to Pinatubo is still well below the 0.5 • C reduction of global surface temperature stated by Crutzen (2006), IPCC (2007), and several other studies.
Figure 12 is complementary to Table 1 and Fig 2)).The thin blue lines are a root sum of squares combination of the statistical uncertainty and the variation in C 1 due to uncertainty in NAA RF 2005 (see text).All regressions begin at January 1900.The end date is driven by availability of data.Regressions using data from GISS and NCDC extend to December 2011, regressions using data from CRU4 extend to December 2010, and the regression using data from BEG runs to May 2010.
SOD regression coefficient.We have also computed sensitivity of the SOD regression coefficient to uncertainty in NAA RF 2005 (Fig. 11) and combined this value using root sum of squares with the confidence interval to find the total uncertainty, shown by the thin, light blue error bars.The six long-term climate records considered in Fig. 12 show a similar pattern: after the AMV index detrended using global SST is introduced into the regression, T VOLCANO falls.If detrending of the AMV index using either LIN (shown) or AF (not shown but nearly identical to LIN) is the correct way to represent the impact of variations in the strength of the AMOC on global temperature, then inferred T VOLCANO is about a factor of 2 smaller than has been previously reported.Table 1 as well as Figs.11 and 12 support our contention that neglect of variations in the strength of the AMOC may have compromised prior estimates of volcanic cooling.

Atmospheric temperature
We now turn our attention to the record of lower tropospheric temperature provided by MSU.Hansen et al. (1993), Lacis and Mischenko (1995) and Soden et al. (2002), often cited in reference to the 0.5 • C cooling due to Pinatubo, also examined lower tropospheric temperature from MSU.
Figure 13 shows an analysis of the global lower tropospheric T from MSU, which covers the time period December 1978 to December 2011.Figure 14 is identical to Fig. 13, except for land.We use the same approach as applied to the 111 yr record of surface T except, since the MSU data record is relatively short, the values of C 4 (the AMV index regression coefficient) have been fixed at values determined from the regression of the CRU4 global temperature record (Fig. 6) and the CRU4 land temperature record (Fig. 9), respectively, for the Figs. 13 and 14.The results in Figs. 13 and 14 must be interpreted with caution because the 32 yr period covers only 2 volcanoes, less than 3 full solar cycles, and we have imposed an external constraint on the effect of variations in the strength of the AMOC.Nonetheless, cooling attributed to volcanoes falls as the AMV index is first introduced into the regression, then as the method used to detrend the AMV index is altered.It has long been known that cooling associated with the eruption of El Chichón in spring 1982 was moderated by ENSO warming (e.g., Thompson, 1995;Wigley et al., 2005).Warming due to ENSO at the time of peak SOD due to El Chichón is picked up well by our model.Maximum cooling due to Pinatubo is found to be 0.50 • C globally and 0.61 • C over land when the AMV index is not included in the regression, which are close to values of volcanic cooling reported by the MSU analyses of Hansen et al. (1993), Lacis and Mischenko (1995) and Soden et al. (2002).
Our estimates of the cooling attributed to Pinatubo drop to 0.35 • C (global) and 0.45 • C (land) when AMV Had3 AF is used in the regression.Model parameters T VOLCANO and λ are both larger for land than global, consistent with the no-tion that lower tropospheric temperature is more sensitive to radiative perturbations over land than over ocean.These estimates for the cooling attributed to Pinatubo are smaller than the values reported by Hansen et al. (1993), Lacis and Mischenko (1995) and Soden et al. (2002).As apparent in our figures, the lower atmosphere did indeed cool by ∼0.5 • C at the time of peak SOD following the eruption of Pinatubo.However, when AMV Had3 AF is used in the regression, a significant portion of this cooling is attributed to variations in the strength of the AMOC rather than Pinatubo.The quality of fit indicated by numerical values of χ 2 improves considerably as the AMV index is added to the regression, and moderately as the method used to detrend AMV index is varied.However, none of the regressions approach χ 2 = 2.This is likely due to a combination of small values for σ OBS of MSU (error bars for MSU T shown in Figs. 13 and 14 are much smaller than error bars for CRU4 T shown in Figs. 6 and 9) as well as climate variability.For instance, the cloud height anomaly reported by Davies and Molloy (2012) shows a sharp drop in late 2007/early 2008 (their Fig. 1) that is well aligned with a drop in MSU T : a decline in cloud height, all else being equal, should cause lower tropospheric cooling.This feature is not picked up by our regression because the climate feedback parameter is assumed to be constant throughout the simulated time period.Lower tropospheric temperature will be more sensitive to short term forcings, such as cloud height anomaly, than global surface temperature.The numerical values of T PINATUBO given above support our suggestion that volcanic cooling may have been over estimated, by about a factor of 2, due to prior neglect of variations in the strength of the AMOC.

Discussion
We have presented an analysis of the climate record, using multiple linear regression, that suggests variations in the strength of the AMOC may be responsible for a portion of the ∼0.3 to 0.5 • C cooling that followed the eruption of Mt.Pinatubo.If correct, then a recalibration of the use of Mt.Pinatubo as a proxy for geoengineering of climate is needed.Here, we focus on the physical understanding of volcanic cooling and the implications of our study for geoengineering of climate.

Physical understanding of volcanic cooling
Many studies of volcanic cooling, including Hansen et al. (1993), Lacis and Mischenko (1995), Soden et al. (2002), Forster and Collins (2004), Wigley et al. (2005) and Stenchikov et al. (2009), use GCMs far more sophisticated than our simple regression analysis.This section is motivated by the fact, illustrated in Fig. TS.23 of IPCC (2007), that most GCMs overestimate the observed perturbation to global mean surface temperature following the eruption of Mt.Pinatubo.The thick red line of Fig. TS.23 of IPCC (2007) shows modelled T from many GCMs exhibited a decline that is much larger than the observed T following the eruption of Pinatubo.
We begin by examining the perturbation to Earth's solar (shortwave, SW) and thermal (longwave, LW) radiation budget after Pinatubo.Figure 15 shows ERBE observations for the tropics (20 • S to 20 • N) and 60 • S to 60 • N. The tropical panel is nearly identical to Fig. 2 of Trenberth and Dai (2007) and the 60 • S to 60 • N panel is similar to Fig. 1 of Soden et al. (2002).The eruption of Mount Pinatubo led to a dramatic increase in stratospheric sulfate aerosol loading, causing a large rise in the reflection of solar radiation due to the optical properties of sulfuric acid droplets.The effect of volcanic aerosols was clearly seen in surface solar radiation measurements: after the eruption of Pinatubo, the amount of direct SW radiation hitting the surface fell dramatically, the amount of diffuse SW radiation rose in tandem, resulting in a net maximum ∼6 W m −2 drop in total SW radiation reaching the surface at Mauna Loa (Dutton and Bodhaine, 2001).However, volcanically induced aerosols also trap thermal radiation.This is less well documented, particularly at the surface.In the tropics, the ERBE record shows that increased reflection of solar radiation following Pinatubo dominates trapping of thermal radiation, resulting in a peak net perturbation to the radiative budget of ∼6.0 W m −2 .In contrast, for 60 • S to 60 • N, the solar and thermal terms are close in magnitude resulting in a peak net perturbation to the radiative budget of less than 3.0 W m −2 .Nearly all of this perturbation occurs in the tropics: the net radiative effect of Pinatubo poleward of 20 • latitude was small in the Northern Hemisphere and essentially zero in the Southern Hemisphere (see the Supplement).Our analysis of ERBE data is in good agreement with Harries and Futyan (2006).Soden et al. (2002) showed that the GFDL GCM was able to simulate, remarkably well, the perturbation to the solar and thermal radiation fields induced by Pinatubo over the 60 • S to 60 • N region.The ERBE measurements represent the radiative perturbation at the top of the atmosphere, not the tropopause.Matching ERBE is a necessary, but not sufficient, condition for accurate simulation of the affect of volcanic aerosols on climate.We discuss below the need for accurate simulation of the stratospheric response.Also, as discussed in Sect. 1, Soden et al. (2002) (and other studies) require a significant positive feedback due to changes in tropospheric H 2 O to match the MSU lower troposphere temperature anomaly.They state "without the strong positive feedback from water vapour, the model is unable to reproduce the observed cooling".Soden et al. (2002) show quite favourable comparisons of modelled and measured upper tropospheric H 2 O that support the notion that a strong, positive water vapour feedback did occur.Forster and Collins (2004) reached similar conclusions as Soden et al. (2002), also based on GCM simulations.Both studies accounted for the effect of ENSO on T and neither study considered the possibility that variations in the strength of the AMOC could have affected T .Wigley et al. (2005) present an analysis of the major volcanic eruptions since 1900, using the National Center for

T. Canty et al.: A critical evaluation of volcanic cooling
Atmospheric Research (NCAR)/US Department of Energy (DOE) GCM.They provide extensive discussion of the effect of ENSO on the inference of volcanic cooling, but do not discuss possible effects due to variations in the strength of the AMOC.Wigley et al. (2005) concluded maximum global surface cooling associated with Pinatubo was 0.61 ± 0.1 • C, consistent with a climate sensitivity of 3.03 • C (range 1.79 to 5.21 • C).Their climate sensitivity refers to equilibrium warming for a doubling of CO 2 ( T 2×CO 2 ).In our notation, T 2×CO 2 is expressed as: where the RF due to CO 2 , 5.35 ln(CO FINAL 2 /CO INITIAL 2 ) W m −2 , is the IPCC (2007) expression originally published by Myhre et al. (1998).Using values for γ from Figs. 6d and 9d, our model would imply T 2×CO 2 of 1.89 • C for global temperature and 2.58 • C for surface land temperature.Our global value is above the lower limit of IPCC ( 2007) and just below the lower limit of Wigley et al. (2005).The land value is within the Wigley et al. (2005) range, although Wigley focused exclusively on global surface temperature.The model value of γ , which drives T 2×CO 2 , is highly dependent on both NAA RF (Fig. 11) and the treatment of OHE (Mascioli et al., 2012).Our model run for NAA RF 2005 = −1.0W m −2 along the "Middle Road" of Fig. 4 and constrained to match the OHC measurement of Gouretski and Reseghetti (2010), for AMV Had3 AF , yields T PINATUBO = 0.14 • C and γ = 1.08, implying T 2×CO 2 = 2.41 • C. The relation between γ , NAA RF, and OHE is further quantified by Mascioli et al. (2012).We include these numerical details to emphasise our MLR results are within the range of previously published GCM results.Douglass and Knox (2005) conducted a regression analysis of MSU lower tropospheric temperature measurements and concluded the atmosphere exhibited a negative feedback following the eruption of Pinatubo.This paper has been discussed in a series of published comments and replies following initial publication, concluding with another paper, Douglass et al. (2006), that reinforces their notion of a negative feedback within the climate system in the short time period following the Pinatubo eruption.Our estimates of climate feedback, given above, represent a best fit to the entire 111 yr temperature record, without direct focus on the time period immediately following major eruptions.We resist the temptation to assess climate feedback immediately after the four major eruptions since 1900 because our framework is based upon stratospherically adjusted RF, and relating the ERBE measurements (top of atmosphere) to the tropopause is beyond the scope of this study.If the true influence of Pinatubo on global cooling is as small as suggested by our lower limits, and there was indeed a strong, global, negative RF anomaly at the tropopause, then perhaps there was a negative feedback following the eruption of Pinatubo, as suggested by Douglass and Knox (2005) and Douglass et al. (2006).Stenchikov et al. (2009) state; "In this study, Pinatubo aerosols globally decrease the incoming net radiative flux at the top of the atmosphere by about 3 W m −2 at maximum that is consistent with most of the IPCC-AR4 models . . .this radiative perturbation dominated all other forcings for at least two years".The radiative perturbation of Pinatubo was enormous for the tropics (Fig. 15a).The radiative perturbation is much smaller when examined over the 60 • S to 60 • N region.Here, the SW perturbation peaks at 3.7 W m −2 , which occurs 2 months after the eruption and averages 1.7 W m −2 for the first 2 yr.ERBE shows a drop in outgoing LW radiation that counteracts the increased reflection of SW radiation.The net perturbation (LW-SW) at the top of the atmosphere, from 60 • S to 60 • N, peaks at −2.5 W m −2 about 2 months after the eruption and has a mean value of −0.67 W m −2 for the first 2 yr.The GCM simulations of Stenchikov et al. (2009) represent this LW trapping of heat by volcanic aerosol.Within their model, LW trapping tends to heat the lower stratosphere, reinforcing the notion that correct modeling of the stratospheric response and the downward influence on the troposphere is vital.
We now turn to the stratosphere.The eruption of Pinatubo induced major changes in stratospheric dynamics and temperature.Tropical stratospheric temperature rose and tropical upwelling increased (e.g., Kinne et al., 1992).Stratospheric ozone fell, first in the tropics (Schoeberl et al., 1993) and then at mid-latitudes (Kinnison et al., 1994).Quantitative understanding of the amount of stratospheric ozone depletion following the eruption of Pinatubo requires realistic representation of the dynamical change (Kinnison et al., 1994) as well as accounting for the influence of both natural and anthropogenic halogen sources (Salawitch et al., 2005).
Figure 9.14 of IPCC (2007) suggests the stratospheric perturbation following the eruption of Pinatubo may not be handled particularly well by GCMs that contributed to IPCC (2007).This figure shows a dramatic drop in the height of the tropopause associated with enhanced SOD within GCMs that appears to be distinctly different than the observed, monthly mean tropopause height anomaly based on the European Center for Medium-Range Weather Forecasts (ECMWF) 40 yr re-analysis (ERA-40).The modelled tropopause height anomaly shows GCMs tend to "collapse the tropopause" when SOD is enhanced.There is a sense of similar behaviour in the observed tropopause after the eruption of Agung (although the observed perturbation is broader in time than the modelled response), opposite behaviour after El Chichón (i.e., observed tropopause rose at the time GCMs suggest it should have fallen), and a confused situation after Pinatubo.In early 1992, the reanalysis shows a well timed but smaller drop in the height of the tropopause than modelled; however, there are tropopause height transients before and after the eruption of Pinatubo in the reanalysis that are not apparent in the models.Perhaps these comparisons are affected by the low-pass filter applied to the data and model; the purpose of Fig. 9.14 is to document the importance of rising GHGs for the proper simulation of the longterm rise in tropopause height.Nonetheless, given the importance of proper quantification of the stratospheric response, future comparisons of modelled and measured tropopause height would be useful for evaluating and perhaps improving GCM simulations of the response of surface temperature to volcanic perturbations.

Implications for geoengineering
The possibility of geoengineering of climate via injection of sulfur to the stratosphere has a long history, starting with Budyko (1974), followed by a remarkably detailed and prescient chapter entitled "Geoengineering" in a National Academy of Sciences report (NAS, 1992), as well as papers by Dickinson (1996) and Schneider (1996).The suggestion by Crutzen (2006) that "if sizable reductions in greenhouse gas emissions will not happen and temperatures rise rapidly, then climatic engineering" by artificial enhancement of stratospheric sulfate could be "the only option available to rapidly reduce temperature rises and counteract other climatic effects" led to a widespread renewal of interest in geoengineering, undoubtedly due to the prominence of Paul Crutzen, but also because cooling after the eruption of Pinatubo had been so well studied by the time his paper was written.Indeed, Crutzen (2006) wrote "enhanced reflection of solar radiation to space by the particles cooled the Earth's surface on average by 0.5 • C in the year following the eruption" (of Pinatubo).Since 2006, many studies (e.g., Rasch et al., 2008a, b;Ammann et al., 2010) have estimated cooling due to stratospheric sulfate injection from first principles (i.e., optical properties of aerosols), rather than using T VOLCANO inferred from Pinatubo.Nonetheless, Pinatubo as an analogy for geoengineering of climate is pervasive in the modern literature: for instance, Robock et al. (2008) suggest the equivalent of one Pinatubo every 4-8 yr would be required to stop global warming (see also Robock et al., 2009).
If further studies support our suggestion that the cooling after the eruption of Pinatubo and other major volcanoes has been over estimated by about a factor of 2, there are three important implications for geoengineering of climate via injection of stratospheric sulfate.First, numerical values from past major volcanoes used as a proxy for cooling by geoengineering would need to be revised.This revision would be straightforward to implement if consensus on the true value of T VOLCANO is achieved.
Second, the behaviour of GCMs used to assess the response of climate to major volcanoes and the response of climate to geoengineering will have to be critically appraised.As described in Sect.5.1, GCMs tend to strongly and consistently collapse the tropopause after major volcanic eruptions, a response not borne out by observation.The atmospheric sciences community has placed enormous ef-fort towards evaluating the chemical, dynamical and radiative behaviour of chemistry-climate general circulation models (CCMs) (Eyring et al., 2010).However, none of the primary GCMs that have quantified volcanic cooling participated in the Eyring et al. (2010) effort.Chapter 4 of Eyring et al. (2010), entitled Stratospheric Dynamics, includes evaluative metrics for 13 CCMs, but does not include the GFDL GCM used by Soden et al. (2002) and Stenchikov et al. (2009), the NCAR PCM model of Wigley et al. (2005), the HadCM3 model used by Forster and Collin (2004), or the HadGEM2-ES model of Booth et al. (2012).Given the importance of stratospheric dynamics for the quantification of both volcanic cooling and geoengineering, we suggest the upcoming Geoengineering Model Intercomparison Project (GeoMIP) (Kravitz et al., 2011) adopt some of the evaluative metrics developed by Eyring et al. (2010) and apply these metrics to all GeoMIP GCMs.Representation of the response of stratospheric ozone to increased sulfate loading (Tilmes et al., 2008(Tilmes et al., , 2009;;Heckendorn et al., 2009) is also needed to properly forecast the response of surface temperature to geoengineering, because ozone is the primary absorber of SW radiation in the stratosphere.Model representation of natural, very short lived (VSL) organic halogen sources may be necessary to properly treat the response of ozone to geoengineering of climate, because ozone loss due to decomposition products of VSL halogen sources is acutely sensitive to aerosol loading in the lowermost stratosphere (Tilmes et al., 2012).
The third implication is subtle.Trenberth and Dai (2007) examined a 55 yr record of global land precipitation, river discharge and the Palmer Drought Severity Index and concluded the eruption of Pinatubo led to precipitation and discharge anomalies much larger than those observed for other years, resulting in widespread drought.However, the AMV index has also been associated with 20th century drought (e.g., Nigam et al., 2011).Figure 10 indicates the AMV index was in a negative phase at the time Pinatubo erupted and that the transition from positive to negative started 2 yr prior to the eruption.If a consensus emerges that the cold SSTs in the North Atlantic during 1992 were driven in part by a process such as sea ice export through the Fram Strait (Dima and Lohmann, 2007) or internal variability of salinity (Delworth and Zeng, 2012), then future analyses of the effects of Pinatubo on the hydrological cycle may have to isolate volcanic and oceanic induced perturbations to serve as a realistic proxy for geoengineering of climate via stratospheric sulfate injection.

Conclusions
The climate-record from 1900 to present exhibits a monotonic rise driven by increasing levels of anthropogenic GHGs (IPCC, 2007), declines after major volcanic eruptions (Hansen et al., 1993;Lacis and Mischenko, 1995; T. Canty et al.: A critical evaluation of volcanic cooling Soden et al., 2002), as well as low frequency, high amplitude variations that have been attributed to changes in the strength of the Atlantic Meridional Overturning Circulation (AMOC) (Schlesinger and Ramankutty, 1994;Andronova and Schlesinger, 2000;Knight et al., 2005;Stouffer et al., 2006;Medhaug and Furevik, 2011;Srokosz et al., 2012;Zhou and Tung, 2013;Tung and Zhou, 2013).Prior studies that quantified volcanic cooling have not considered variations in the strength of the AMOC.We have shown, using multiple linear regression (MLR) analysis of many climate records, that cooling attributed to volcanic eruptions is reduced by about a factor of 2 when sea surface temperatures in the North Atlantic are used as a proxy for the effect of AMOC on global climate.
The surface temperature records from CRU4, GISS, NCDC and BEG as well as lower tropospheric temperature from MSU, analysed with the MLR model and excluding the Atlantic Multidecadal Variability (AMV) index, yield values for the cooling associated with the eruption of Mt.Pinatubo that range from 0.31 • C to 0.61 • C.This range is broadly consistent with values for Pinatubo cooling reported by Hansen et al. (1993), Lacis and Mischenko (1995), Soden et al. (2002), Lean and Rind (2008), Thompson et al. (2009) and Foster and Rahmstorf (2011).The cooling associated with Mt.Pinatubo drops nearly a factor of 2, ranging from 0.14 • C to 0.45 • C, when these same datasets are analysed allowing for variations in the strength of the AMOC in the regression, based on an AMV index detrended by anthropogenic RF.
All of our modelled temperature anomalies drop by ∼0.5 • C at the time of peak SOD following the eruption of Mt.Pinatubo.However, a significant portion of this cooling may be due to ocean circulation and not volcanoes, based on results of our regression analysis.The timing of stratospheric optical depth (SOD) and North Atlantic SST anomalies suggest changes in ocean circulation preceded the four major volcanic eruptions that have occurred since 1900.Hence, the data record suggests ocean circulation has affected prior estimates of volcanic cooling, rather than volcanic eruptions drive ocean circulation.We use a detrended AMV index as a proxy for the AMOC: i.e., this study considers variability in the strength of the AMOC rather than possible long-term monotonic change.Precise determination of volcanic cooling is sensitive to the manner in which the AMV index is detrended.If global SST is used to detrend AMV, as suggested by Trenberth and Shea (2006), the influence of the AMOC on volcanic cooling is moderate.If AMV is detrended using either a linear regression (Enfield et al., 2001) or, as we suggest, anthropogenic RF of climate (which varies in a nonlinear manner over time), then the influence of the AMOC on volcanic cooling could be strong.The regression based on AMV detrended using anthropogenic RF simulates pre-WWI cooling and WWII warming of the CRU4 land temperature anomaly particularly well (Fig. 9d), giving credence to the possibility that variations in the strength of the AMOC truly have exerted influence on global climate.If a consensus emerges that variations in the strength of the AMOC do indeed exert a major influence on global climate, as suggested throughout this paper, then this may provide an important new opportunity to improve estimates of temperature and climate forecasts on decadal time scales, as suggested by Knight et al. (2005).
If the factor of 2 reduction of volcanic cooling suggested here is borne out by future studies, the implications for geoengineering of climate by artificial enhancement of stratospheric sulfate are immense.It would be straightforward to "recalibrate" Pinatubo as a proxy for geoengineering.Of greater concern is the fidelity of atmospheric general circulation models (GCMs) used to assess the response of the atmosphere to major volcanoes as well as geoengineering.Accurate GCM representation of the response to either perturbation requires realistic treatment of a myriad of physical processes, including the trapping of longwave thermal radiation by stratospheric sulfate aerosols, the dynamical response of the stratosphere, and changes in stratospheric ozone.The tendency of the GCMs used by IPCC (2007) to collapse the tropopause following major perturbations to SOD, which is not borne out by the ERA-40 reanalysis, suggests the physical response to stratospheric sulfate aerosol injection is not properly represented in modern climate models.

Appendix A
Glossary of symbols and terms α COOL : parameter that scales the direct RF of climate due to aerosols that cool the atmosphere (negative RF) to total RF; units: dimensionless α HEAT : parameter that scales the direct RF of climate due to aerosols that heat the atmosphere (positive RF) to total RF; units: dimensionless -T : Global, monthly mean temperature anomaly for either the surface or lower troposphere.Here, we consider four records of T : the global surface (land and oceans), the land surface, the global lower troposphere (land and oceans), and the lower troposphere over land; units: • C -T : Global, annual mean temperature anomaly used here for computation of χ 2 (Eq.7); units: • C -T 2×CO 2 : Equilibrium climate sensitivity for a doubling of CO 2 relative to pre-industrial levels, which in our modelling framework is equal to  et al. (2009), using a first order differential equation that takes into consideration the spatial variability of variations in SST as well as the effective heat capacity of the atmospheric-oceanic mixed layer; units:

Fig. 1 .
Fig. 1.Direct RF due to greenhouse gases (GHG RF) used as input for all model calculations.The coloured regions show individual contributions from CO 2 (red), CH 4 (blue), tropospheric O 3 (orange), halocarbons (green) and N 2 O (purple), all based on global, monthly mean mixing ratios from the RCP 8.5 scenario.
Fig. 2. (a) Global, annual anthropogenic emission of sulfur (S EMISS ) from RCP Potsdam compared to values published by Stern (2006b) and Smith et al. (2011).(b) Direct RF due to sulfate aerosols (RF SULFATE−DIR ) from RCP Potsdam (post 2005) and the estimate of RF SULFATE−DIR (1900 to 2005) used in our model, labelled Smith * , which is tied to S EMISS from Smith et al. (2011) and the IPCC (2007) estimate of −0.4 W m −2 for year 2005.(c) Total RF due to sulfate aerosols (RF SULFATE−TOT ) from David Stern, Australian National University, provided at the URL given in Appendix B, compared to the value of RF SULFATE−TOT used in our model, labelled Smith * , found by multiplying Smith * RF SULFATE−DIR by α COOL , for a value of α COOL = 2.4, chosen to yield match our estimate of RF SULFATE−TOT of −0.96 W m −2 for year 2005.Panels are all extended to 2060 for illustrative purposes, to support the statement in the paper that tropospheric aerosol forcing of climate will be diminishing over time.
Fig. 3. (a) Total RF of tropospheric aerosols that cool, as labelled, based on our Smith* estimate of RF SULFATE−DIR and RCP 8.5 estimates of direct RF for other components, all multiplied by α COOL = 2.4.The curve labelled Sum denotes total RF due to aerosols that cool.(b) Same as (a), except for aerosols that heat.Direct RF components, from RCP 8.5, have been multiplied by α HEAT = 2.4, chosen so that net anthropogenic aerosol RF (NAA RF) in year 2005 equals −1.0 W m −2 (IPCC, 2007).The curve labelled Biomass refers to emissions of OC and BC due to biomass burning, and the curves labelled OC and BC refer to fossil fuel burning emissions of these components.(c) total RF of aerosols that cool (blue, found for α COOL = 2.4), of aerosol that heat (red, found for α HEAT = 2.4), and their difference that defines NAA RF (black curve).The value of NAA RF in year 2005 is marked.It is coincidence that the two scaling parameters both require a value of 2.4 to match IPCC (2007) estimates of RF SULFATE−DIR and NAA RF in year 2005.Model results are shown for a variety of values of α COOL and α HEAT .
Fig. 3. (a) Total RF of tropospheric aerosols that cool, as labelled, based on our Smith* estimate of RF SULFATE−DIR and RCP 8.5 estimates of direct RF for other components, all multiplied by α COOL = 2.4.The curve labelled Sum denotes total RF due to aerosols that cool.(b) Same as (a), except for aerosols that heat.Direct RF components, from RCP 8.5, have been multiplied by α HEAT = 2.4, chosen so that net anthropogenic aerosol RF (NAA RF) in year 2005 equals −1.0 W m −2 (IPCC, 2007).The curve labelled Biomass refers to emissions of OC and BC due to biomass burning, and the curves labelled OC and BC refer to fossil fuel burning emissions of these components.(c) total RF of aerosols that cool (blue, found for α COOL = 2.4), of aerosol that heat (red, found for α HEAT = 2.4), and their difference that defines NAA RF (black curve).The value of NAA RF in year 2005 is marked.It is coincidence that the two scaling parameters both require a value of 2.4 to match IPCC (2007) estimates of RF SULFATE−DIR and NAA RF in year 2005.Model results are shown for a variety of values of α COOL and α HEAT .

Fig. 4 .
Fig. 4. Contours of net anthropogenic aerosol RF in year 2005 (NAA RF 2005 ) (black solid lines) as a function of α COOL and α HEAT , scaling parameters used to relate direct RF of aerosols to total RF.The contour for NAA RF 2005 = −1.0W m −2 , the best estimate from IPCC (2007), is shown in red.The dashed green lines denote the range of NAA RF 2005 inferred from data analyses given in Table2.12 of IPCC (2007) (Appendix C).The black dashed/yellow highlight lines denote various manners upon which values of NAA RF 2005 , ranging from the lower limit of −2.2 W m −2 to the upper limit of −0.4 W m −2 , can be sampled.
Fig. 4. The red line and black solid lines on this figure denote isopleths of NAA RF 2005 .The green dashed lines represent the empirical range for NAA RF 2005 , which we have computed as −0.4 W m −2 to −2.2 W m −2 based on Table Fig. 5. (a) Area weighted, monthly mean SST in the North Atlantic (equator to 60 • N) based on HadSST3 data.(b) AMV index (blue and red) found by detrending North Atlantic SST using area weighted SST between 60 • S and 60 • N (gray curve).(c) AMV index (blue and red) found by calculated by detrending North Atlantic SST using a linear regression (gray line).(d) AMV index (blue and red) found by detrending North Atlantic SST using anthropogenic RF of climate (gray line).The particular anthropogenic RF curve is an MLR simulation constrained to match the CRU4 global temperature record, NAA RF 2005 = −1.0W m −2 along the Middle Road, and OHC from Church et al. (2011).
For the Church et al. (2011) measurement of OHC, Eq. (5) becomes: γ ){GHGRF + NAA RF} 1944 to 2003 (6) The notation 1944 to 2003 denotes the mean value of the term enclosed within brackets over the 1944 to 2003 time period.Years 1944 to 2003 are used for the average of anthropogenic RF in the denominator of the expression because of the 6 yr delay between atmospheric perturbation and upper ocean response: OHC was measured by Church et al. (2011) from 1950 to 2009.For each iteration of the model involved with minimisation of the Cost Function, as λ adjusts to prescribed NAA RF, is updated.The red curve on the bottom panel of Fig. 6a compares the integral over time of the modelled value of Q OCEAN with the Church et al. (

Fig. 6 .
Fig. 6.Ladder plots showing the MLR simulation of the global temperature anomaly from CRU4.The top rung of each of the four ladder plots compares measured (black) and modelled (red) T , from the start of 1900 to the end of 2010.Values of χ 2 are given and 1-sigma uncertainty of T , available for every month, is shown periodically (blue error bars).The other rungs show contributions to T from volcanoes (based on SOD), solar (based on TSI), humans (sum of GHG RF and NAA RF terms in Eq. (2)), and the various ocean terms (see text).Values of the regression coefficients are given; for the volcanic coefficient, the statistical uncertainty from the regression is also presented.Specified values of NAA RF 2005 and output values of the sensitivity parameter are given in the rung labelled Human.(a) (Top) Results of MLR for a simulation where the regression coefficients for AMV, PDO and IOD have been set to zero.(Bottom) Modelled and measured ocean heat content of the upper 700 m of the world's oceans.Data are from Church et al. (2011) except 9.5×10 22 J has been added to the measurements at all times so that the mean value of OHC of our model matches the mean measured value of OHC for the time period of observation (measurement is an anomaly so this type of adjustment is valid).The red line is the modelled value of 0.7 × Q OCEAN × 3.3 × 10 14 m 2 × t, where t is the time evolved since the start of 1900.(b) Same as top part of (a) but also including rungs showing contribution to the regression from variations in SST within the Atlantic (based on AMV Had3 SST ), Pacific (PDO) and Indian (IOD) Oceans.(c) Same as (b) but for AMV Had3 LIN .(d) Same as (b) but for AMV Had3 AF .
Fig. 7. (a) Same as Fig.6a, except indices for ENSO and SOD have been replaced by values for the response of global temperature to ENSO and volcanic eruptions found by the thermodynamic model ofThompson et al. (2009), with no temporal delay.(b) Same as Fig.6d, except the response of global temperature to ENSO and volcanic eruptions found byThompson et al.(2009) are used.(c) Same as (b), except the AMV index has been filtered to remove all components with frequencies higher than 1/9 yr −1 .(d) Same as (b), except the AMV index has been filtered to remove all components with frequencies higher than 1/73 yr −1 .

Fig. 8 .
Fig. 8. Values of the AMV index regression coefficient (C 4 in Eq. (2)) (a) and maximum cooling due to the eruption of Mt.Pinatubo (absolute value of C 1 × 0.15) (b) calculated using the MLR model (large black circles) and from a series of conditional regressions (small gray circles; see text).Black and gray error bars denote the 95 % confidence intervals(von Storch and Zwiers, 2002).The red horizontal bar denotes the mean of the 24 conditional regressions, and the red error bar denotes the 2σ standard deviation about the mean.For each panel, the first six conditional regression points represent calculations for which C 1 was first computed.The second sextet show results for which C 4 was first computed.For the third sextet, C 3 was first computed and for the last, C 2 was first computed.

Fig. 9 .
Fig. 9. Same as Fig. 6, but for MLR simulation of the land temperature anomaly from CRU4.
Fig. 10.(a) Monthly mean stratospheric optical depth (SOD) fromSato et al. (1993).The four major volcanoes, Santa María (purple), Mount Agung (green), El Chichón (red) and Mt.Pinatubo (blue) are indicated.(b) Monthly mean Atlantic Multidecadal Variability index calculated by detrending HadSST3 using anthropogenic radiative forcing (AMV Had3 AF ).Line is the same as shown in Fig.5d, except AMV index during the times SOD exceeded 0.01, following the four major volcanic eruptions since 1900, are coloured as in (a).(c) Scatter plots of SOD versus AMV index.In the centre, SOD is plotted versus AMV index with no time shift.Data collected during the four major volcanic eruptions since 1900 are coloured as in (a).The numerical values on the top of the plot tabulate the number of months the AMV index was either positive or negative, when SOD exceeded 0.01 after these four major eruptions.The left hand panel shows a scatter plot for a shift of SOD six months earlier in time and the right-hand panel shows a scatter plot for a shift of SOD six months later in time.
Fig. 11.(a) Climate feedback parameter (λ in Eq. (2)) found for a regression of monthly mean global surface temperature anomaly from CRU4, as a function of NAA RF 2005 (using values of α COOL and α HEAT along the "Middle Road" of Fig. 4), for different treatments of the AMV index: none, AMV Had3 SST , AMV Had3 LIN and AMV Had3 AF .(b) Same as (a) except for a regression of the CRU4 monthly mean land surface temperature anomaly.(c) Maximum cooling due to the eruption of Mt.Pinatubo found by the regression of monthly mean global surface temperature anomaly from CRU4 as a function of NAA RF 2005 for different treatments of the AMV index.Here we show T PINATUBO as a negative numerical value to be consistent with the depiction T VOLCANO in other figures.We provide positive values of T PINATUBO in the context of "cooling", in the text, to avoid the use of a myriad of minus signs.(d) Same as (c) except for a regression of CRU4 monthly mean land surface temperature anomaly.
Figure 11a and b show how the climate feedback parameter λ (Eq. 3) varies as a function of NAA RF 2005 , for values of α COOL and α HEAT along the middle road of Fig. 4. End points of the line segments denote a range of NAA RF 2005 between −2.0 W m −2 and −0.4 W m −2 , similar to the empirical range presented in IPCC (2007) (Appendix C).Model results for NAA RF 2005 = −1.0W m −2 are included, as well as the range of NAA RF 2005 represented in GCMs described by

Fig. 12 .
Fig. 12. Cooling attributed to volcanoes ( T VOLCANO ) as a function of stratospheric optical depth (SOD), for regression of monthly mean global temperature anomalies ( T ) from various data centers, as indicated (after "T Data:").The black lines each represent C 1 × SOD versus SOD, where C 1 is the volcanic regression coefficient.Results for no AMV (solid lines), AMV Had3 LIN (dotted), and AMV Had3 SST (dashed), assuming NAA RF 2005 = −1.0W m −2 along the "Middle Road" of Fig. 4, are shown for each panel.Maximum SOD associated with the eruptions of Santa María, Mount Agung, El Chichón and Mt.Pinatubo from Sato et al. (1993) are denoted on each panel by the coloured arrows marked SM, A, C and P. The thick dark blue vertical error bars denote the 95 % confidence interval of the volcanic regression coefficient (C 1 , Eq. (2)).The thin blue lines are a root sum of squares combination of the statistical uncertainty and the variation in C 1 due to uncertainty in NAA RF 2005 (see text).All regressions begin at January 1900.The end date is driven by availability of data.Regressions using data from GISS and NCDC extend to December 2011, regressions using data from CRU4 extend to December 2010, and the regression using data from BEG runs to May 2010.

Fig. 13 .
Fig. 13.Same as Fig. 6, but for analysis of the global, monthly mean LT5.4 (lower troposphere) temperature anomaly measured by Microwave Sounding Unit (MSU) instruments.The analysis is from December 1978 to December 2011.

Fig. 14 .
Fig. 14.Same as Fig. 11, but for the MSU land, monthly mean lower troposphere temperature anomaly.
Fig. 15.(a) Deseasonalised times series of shortwave, longwave and net radiation anomalies for 20 • N to 20 • S, from the start of 1985 to end of 1999, with respect to a 1985 to 1989 (pre-Pinatubo) baseline.Data are from ERBE Edition 3 Revision 1, non-scanner, wide fieldof-view observations(Wielicki et al., 2002).The gray shaded region denotes the period of time from the eruption of Pinatubo (15 June 1991) until the end of 1992.For this latitude region, raw data are provided as 36-day means.(b) Same as (a), but for 60 • N to 60 • S. For this latitude region, raw data from the tropics (36-day means) have been combined with raw data from the extra-tropics (72-day means; extra-tropics refer to 20 • to 60 • in both hemisphere) to produce a 72-day mean, near global average (using latitudinal weighting).

--
RF CO 2 = 5.35 W m −2 ln(2) or 3.71 W m −2 ; units: • C -T PINATUBO : maximum cooling following the eruption of Mt.Pinatubo in our model framework, which equals the absolute value of C 1 × 0.15, where 0.15 is the largest value of SOD observed after this eruption; units: • C -T VOLCANO : the effect of the volcanic perturbation to stratospheric optical depth on surface temperature, which equals C 1 × SOD i−6 in our model framework; units: • C γ : sensitivity of climate to all feedbacks that occur in response to a GHG perturbation to RF at the tropopause; related to λ by Eq. (3); positive values represent amplification of RF due to internal feedbacks and negative values represent dampening; units: dimensionless λ p : Planck feedback parameter; response of surface temperature to a change in the emission to space of long-wave thermal radiation in the absence of any feedbacks for Earth's present-day overall albedo; value of 3.2 W m −2 • C −1 used throughout λ: climate feedback parameter; response of surface temperature to an external RF at the tropopause to internal feedbacks such as water vapour, clouds, albedo, or lapse rate; the total climate feedback parameter, shown in numerous figures, represents the sum of individual components; positive values of individual components represent amplification of RF and negative values represent dampening; units: W m −2 • C −1 -: Fraction of RF of climate due to anthropogenic GHGs and aerosols, modified by climate feedback, that is expended by heating the world's oceans; units: dimensionless -AF: Anthropogenic radiative forcing of climate; units: W m −2 -AMO: Atlantic Multidecadal Oscillation -AMOC: Atlantic Meridional Overturning Circulation -AMV: Atlantic Multidecadal Variability -AMV Had3 AF : AMV index found by detrending monthly mean SST in the Atlantic, from 0 to 60 • N, based on the HadSST3 record using anthropogenic forcing of climate (GHG RF + NAA RF); units: • C -AMV Had3 LIN : AMV index found by detrending monthly mean SST in the Atlantic, from 0 to 60 • N, based on the HadSST3 record, using a linear regression; units: • C -AMV Had3 SST : AMV index found by detrending monthly mean SST in the Atlantic, from 0 to 60 • N, based on the HadSST3 record, using global SST (also from HadSST3); units: • C -CRU: Climate Research Unit, University of East Anglia -CRU4: Latest surface temperature record from CRU; HadCRUT4 is the latest global temperature record and CRUTEM4 is latest land temperature record -GCM: Atmosphere ocean general circulation model -ENSO: El Niño-Southern Oscillation or El Niño-Southern Oscillation Index, units: dimensionless -GCM: General Circulation Model -GHG RF: Direct RF due to Greenhouse Gases, here taken to be CO 2 , CH 4 , Tropospheric O 3 , Halocarbons, and N 2 O -HadSST3: Latest sea surface temperature record from the Hadley Centre -KaplanSST2: Latest sea surface temperature record from NOAA, termed Kaplan Extended SST V2 -LW: Longwave or thermal radiation -MEI: Multivariate ENSO Index; units: dimensionless -MLR: Multiple Linear Regression -MSU: Microwave Sounding Unit -NAA RF: Net Anthropogenic Aerosol RF; units: W m −2 -OHC: Ocean Heat Content of the upper 700 m of the global ocean; units: J -OHE: Ocean Heat Export (the time derivative of OHC divided by the ocean surface area considered); units: W m −2 -Q OCEAN : Model representation of Ocean Heat Export; units: W m −2 -RCP: Representative Concentration Pathways -RF: stratospheric-adjusted radiative forcing of climate, as described in Sect.2.2 of IPCC, 2007; units: W m −2 RF SULFATE−DIR : Direct RF due to sulfate aerosols; units: W m −2 RF SULFATE−TOT : Total RF due to sulfate aerosols; units: W m −2 -S EMISS : Emission of sulfur, a precursor of sulfate aerosols; units: Tg yr −1 -SOD: Stratospheric Optical Depth; units: dimensionless -SST: Sea Surface Temperature; units: • C -SW: Short wave or solar radiation www.atmos-chem-phys.net/13etal.: A critical evaluation of volcanic cooling -T ENSO : Simulated response of global mean surface temperature to ENSO variability found by Thompson

Canty et al.: A critical evaluation of volcanic cooling measurements
of atmospheric temperature, from December 1978 to present, from a series of National Oceanic and Atmospheric Administration (NOAA) satellites.We use the global LT5.4 (lower troposphere) product provided by the University of Alabama, Huntsville MSU: the Microwave Sounding Unit (MSU) and the Advanced Microwave Sounding Unit (AMSU) provide www.atmos-chem-phys.net/13/3997/2013/Atmos.Chem.Phys., 13, 3997-4031, 2013 T.

Table 1 .
Maximum cooling due to Mt. Pinatubo for all model simulations.

Table 2 .
Values of reduced chi-squared (Eq.7) calculated for all model simulations.