%% Version 4.3.2, 25 August 2014 % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Template.tex -- LaTeX-based template for submissions to the % American Meteorological Society % % Template developed by Amy Hendrickson, 2013, TeXnology Inc., % amyh@texnology.com, http://www.texnology.com % following earlier work by Brian Papa, American Meteorological Society % % Email questions to latex@ametsoc.org. % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % PREAMBLE %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% Start with one of the following: % DOUBLE-SPACED VERSION FOR SUBMISSION TO THE AMS \documentclass{ametsoc} % TWO-COLUMN JOURNAL PAGE LAYOUT---FOR AUTHOR USE ONLY % \documentclass[twocol]{ametsoc} \usepackage{lineno} %\usepackage{ametsoc} \usepackage{graphicx}%%%\usepackage{graphicx} % CIG \usepackage{tabularx} \usepackage{rotating} \usepackage{amsmath}%%%\usepackage{amsmath} \usepackage{amssymb}%%%\usepackage{amssymb} \usepackage{epstopdf}%%%\usepackage{epstopdf} \usepackage{paralist} \usepackage{multirow} \usepackage{lineno} \usepackage{pdflscape} \usepackage{amssymb}%%%\usepackage{amssymb} \usepackage{epstopdf}%%%\usepackage{epstopdf} \usepackage[normalem]{ulem} \usepackage[font={small}]{caption} \linenumbers*[1] \definecolor{orange}{RGB}{255,127,0} \definecolor{black}{RGB}{0,0,0} \newcommand{\TO}[1]{\textcolor{orange}{#1}} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%% To be entered only if twocol option is used \journal{jcli} % Please choose a journal abbreviation to use above from the following list: % % jamc (Journal of Applied Meteorology and Climatology) % jtech (Journal of Atmospheric and Oceanic Technology) % jhm (Journal of Hydrometeorology) % jpo (Journal of Physical Oceanography) % jas (Journal of Atmospheric Sciences) % jcli (Journal of Climate) % mwr (Monthly Weather Review) % wcas (Weather, Climate, and Society) % waf (Weather and Forecasting) % bams (Bulletin of the American Meteorological Society) % ei (Earth Interactions) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %Citations should be of the form ``author year'' not ``author, year'' \bibpunct{(}{)}{;}{a}{}{,} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%% To be entered by author: %% May use \\ to break lines in title: \title{The effect of atmospheric transmissivity on model and observational estimates of the sea ice albedo feedback?} % Corresponding author mailing address and e-mail address: \authors{Aaron Donohoe\correspondingauthor{ Polar Science Center, University of Washington, Seattle, Washington/USA}} \affiliation{Applied Physics Laboratory, University of Washington, Seattle, Washington} \email{adonohoe@u.washington.edu} %---------------------------------------------------------------------------------------- \extraauthor{Ed Blanchard-Wrigglesworth} \extraaffil{Department of Atmospheric Sciences, University of Washington, Seattle, WA, USA} \extraauthor{Axel Schweiger} \extraaffil{Polar Science Center, Applied Physics Laboratory, University of Washington, Seattle, Washington, USA} \extraauthor{Philip J. Rasch} \extraaffil{Pacific Northwest National Laboratory, Richland, WA, USA} \abstract{The sea-ice-albedo feedback (SIAF) is the product of the ice sensitivity (IS) -- how much the surface albedo in sea-ice regions changes as the planet warms-- and the radiative sensitivity (RS) -- how much the top of atmosphere radiation changes as the surface albedo changes. We demonstrate that the RS calculated from radiative kernels in climate models is reproduced from calculations using the ``approximate partial radiative perturbation'' method that uses the climatological radiative fluxes at the top of atmosphere and the assumption that the atmosphere is isotropic to shortwave radiation. This method facilitates the comparison of RS from satellite-based estimates of climatological radiative fluxes with RS estimates across a full suite of coupled climate models and, thus, allows model evaluation of a quantity important in characterizing the climate impact of sea ice concentration changes. The satellite based RS is within the model range of RS that differs by a factor of two across climate models in both the Arctic and Southern Ocean. Observed trends in Arctic sea ice are used to estimate IS which, in conjunction with the satellite-based RS yields an SIAF of 0.16 $\pm$ 0.04 W m$^{-2}$ K$^{-1}$. This Arctic SIAF estimate suggests a modest amplification of future global surface temperature change by approximately 14$\%$ \TO{relative to a climate system with no SIAF}. We calculate the global albedo feedback (GAF), \TO{including changes in snow cover over land,} in climate models using model specific RS and IS and find a model mean feedback parameter of 0.37 W m$^{-2}$ K$^{-1}$ which is 40$\%$ larger than the IPCC AR5 estimate based on using RS calculated from radiative kernel calculations in a single climate model.} \begin{document} %% Necessary! \maketitle %\begin{itemize} %\item The surface albedo feedback calculated from radiative kernels in climate models is reproduced from the climatological radiative fluxes at the top of atmosphere and surface by assuming the atmosphere is isotropic %\item Climatological satellite based radiation provide an observational estimate of the surface albedo magnitude %\item The surface albedo feedback differs in climate models by a factor of two in both the Arctic and Southern ocean with ensemble mean values near the observational estimate %\end{itemize} \section{Introduction} Sea ice area is expected to decrease as the climate system warms, and this in turn will lead to a darker surface, and increase in solar radiation absorbed by the climate system. This additional radiative input reinforces the initial warming providing a positive climate feedback often termed the sea-ice-albedo feedback (SIAF). Early literature on climate stability in simplified models suggested that SIAF could cause abrupt and dramatic climate state transitions under smoothly varying external forcing \citep{north_84, budyko_69} or produce\ multiple equilbria in more comprehensive coupled climate models \citep{ferreira_multiple}. More modest estimates of the global albedo feedback \TO{(including changes associated in surface albedo over land)} were found in coupled climate models \citep{ipcc2013, bony_et_al, soden_feedback}, producing an IPCC AR5 ensemble mean global SIAF of 0.26 W m$^{-2}$K$^{-1}$ \citep{ipcc_models} leading to a \TO{~22$\%$} increase in the global climate response to external forcing \citep{roe_feedbacks} \TO{relative to system with no surface albedo feedback}. \citet{pistone, pistone2} used the co-variance of year-to-year sea ice anomalies and satellite radiation to produce an observationally based estimate of SIAF with a similar magnitude for the Arctic sea ice (0.31 W m$^{-2}$K$^{-1}$) and \TO{pointed out this additional radiative input to the climate system due to Arctic ice melt to date ~25$\%$ the anthropogenic forcing}. There is still a substantial ($\pm$ 0.1 W m$^{-2}$K$^{-1}$) inter-model spread in strength of the SIAF \citep{winton_albedo, hall_06} that is understood to be the leading cause of inter-model differences in the high latitude climate response (polar amplification) \citep{hall_04, kay_amp} \citep{holland_bitz_03}. SIAF measures how much additional radiative energy the Earth system gains due to sea ice loss as the planet warms, which amplifies the warming relative to a system with no SIAF. SIAF is quantified as the global (area weighted) average of RI$_{TOA, \alpha}$, the Radiative Impact of sea ice change (the local TOA radiative flux change due to surface albedo changes ($\alpha$) from sea ice loss per degree of global averaged surface temperature change): %The SIAF is defined as the (area weighted) \emph{global average} TOA radiative impact of sea ice changes per unit of global mean surface temperature change (): \begin{equation} SIAF = \left[RI_{TOA, \alpha}(x,y)\right] \end{equation} \noindent where $\left[ \right]$ brackets indicate a global average. Following \citet[][Eq. 1]{winton_albedo}, the spatial map of RI$_{TOA, \alpha}$(x,y) is the product of two quantities \citep{soden_kernel, shell_kernel}: 1.) the surface albedo change due to sea ice loss per unit of global mean surface temperature change, $\left[dT_{S}\right]$, ($\frac{d\alpha_{SI}}{\left[dT_{S}\right]}$), and 2.) the sensitivity of top of atmosphere (TOA) radiation to surface albedo ($\frac{\partial RAD_{TOA}}{\partial \alpha}$) that we hereafter refer to as radiative sensitivity (RS) : \begin{equation} RI_{TOA, \alpha}(x,y)= \underbrace{\frac{d\alpha_{SI}}{\left[dT_{S}\right]}}_{IS(x,y)} \underbrace{\frac{\partial RAD_{TOA}(x,y)}{\partial \alpha(x,y)}}_{RS(x,y)} . \end{equation} \noindent The normalization of $RI_{TOA, \alpha}(x,y)$ by global mean temperature (T$_{S}$)change is integrated into the IS term and RS is defined as the local radiative change at the TOA per unit of surface albedo change. This study considers only the radiative impact of $\alpha$ changes in high latitudes (poleward of 60$^{\circ}$N and 55$^{\circ}$S, in the northern (NH) and Southern (SH) hemispheres respectively) over oceans, and calculations of SIAF exclude the impact of changes in terrestrial snow cover. \TO{RS and $\alpha$ changes are calculated for each month and then their product is time averaged.} Changes in $\alpha_{SI}$ are calculated over the ocean and capture both the impact of sea ice loss and changes in surface albedo over sea ice (i.e. snow and melt ponds). \citet{hall_06} \TO{claim} that RS varies very little between climate models. As a result, much of the literature on SIAF uncertainty has focused on processes controlling sea ice albedo changes \citep{horvat_albedo} and the sensitivity of \TO{sea ice concentration (SIC)} to warming \citep{winton_albedo, qu_hall, curry_albedo} which both vary substantially between models. The IPCC estimate of SIAF \citep{ipcc_models, soden_kernel} used a RS calculated from a single model, neglecting inter-model differences and biases (relative to observations) and assuming RS does not contribute to SIAF uncertainty. We assess the validity of this assumption below. RS depends primarily on cloud reflectivity; clouds impede the amount of downwelling solar radiation reaching the surface and also reduce the amount of solar radiation reflected by the surface from reaching the TOA \citep{taylor_et_al_07, donohoe_battisti_11a} leading to a quadratic dependence of RS on cloud reflectivity. High latitude cloud properties vary substantially between models and exhibit many biases relative to observations \citep{gorod_cloud, varvus_clouds,tren_fasullo_10}. Cloud differences can contribute to model differences in RS that in turn influence: (i) the sensitivity of sea ice loss to future warming \citep{hwang_amp} via local positive radiative feedbacks and (ii) the impact of sea ice loss on the global energy budget and, thus, the global climate sensitivity to external forcing. This study assesses inter-model differences in RS and consistency compared to estimates from satellite observations. We also identify relative contributions of IS and RS to model spread and biases (relative to observations) in the amplification of global warming by SIAF, and evaluate the impact of using RS from a single climate model to calculate SIAF across models as was done in \citet{soden_kernel} and the IPCC AR5 estimate of surface albedo feedback. The manuscript is organized as follows: section 2, outlines how a simplified isotropic model often discussed in textbooks on radiative transfer, and further developed by \citet{taylor_et_al_07} can be used to calculate RS from standard climate model output and demonstrates that the method reproduces results from more computationally demanding radiative kernel techniques. This facilitates further evaluation of inter-model spread in RS in the coupled models participating in the Coupled Model Inter-comparison Project \citep[CMIP3 and CMIP5][]{meehl_et_al, cmip5}. Most importantly, this method also provides an observational estimate of RS from satellite data (Section 3). These estimates of RS along with the sea ice response over the historical period are used to calculate an observational SIAF (Section 4). The observational SIAF is compared to that in model simulations under historical forcing and 4XCO$_{2}$ and the model spread and biases are decomposed into contributions from RS and IS (Section 5). A summary and discussion follows. \section{The impact of surface albedo changes on TOA radiation in radiative kernels and a simplified model} \subsection{Radiative kernels} The impact of surface albedo changes on TOA radiation (RS) has been rigorously calculated using radiative kernel techniques in a small number of climate models \citep{smith_kernel,angie_kernel, shell_kernel, block,soden_kernel, previdi_10}. RS can be calculated directly from offline radiative model calculations by prescribing changes to the surface albedo ($\alpha$) at each grid point and then running the radiative code with all other fields unchanged -- a technique referred to as a radiative kernel calculation \citep{soden_kernel, shell_kernel}. Radiative kernels are generally calculated \TO{at each gridpoint over a global domain} by perturbing the surface albedo at each grid point by a specified amount (independent of whether that surface albedo change is feasible) using atmospheric models with prescribed historical climatological (seasonally varying) sea surface temperatures. We use kernel calculations (for specific models) provided by: 1.) Karen Shell, NCAR CAM3 \citep{shell_kernel}, 2.) Karoline Block, MPI ECHAM6 \citep{block}, 3.) Angie Pendergrass, NCAR CAM5 \citep{angie_kernel}, 4.) Chris Smith, UKMO HadGEM2 \citep{smith_kernel}, 5.) Brian Soden, GFDL AM2p12b \citep{soden_kernel}. RS is reported in W m$^{-2}$ $\%^{-1}$ where the $\%$ refers to a 0.01 unit change in surface albedo (independent of the climatological surface albedo). Summertime (MJJA) daily-averaged TOA insolation in the Arctic (defined as the region poleward of 60$^{\circ}N$) is of order 420 W m$^{-2}$, and a 4.2 W m$^{-2}$ $\%^{-1}$ RS would be expected in a completely transparent atmosphere. Radiative kernel calculations produce an RS Arctic average of 1.63 W m$^{-2}$ $\%^{-1}$ across the 4 different models (numbers in the upper right of each panel in Fig. 1) indicating that the atmosphere attenuates the surface contribution to reflected radiation at the TOA by a factor of $\sim$ 2.6 (4.2/1.63). Kernel estimates of RS in Arctic summer (May-June-July-August) are largest over Greenland (2-3.5 W m$^{-2}$ $\%^{-1}$), smallest in the Greenland Iceland Norwegian (GIN) Seas (0.5-1 W m$^{-2}$ $\%^{-1}$) with intermediate values in the central Arctic (1-2.5 W m$^{-2}$ $\%^{-1}$ -- Upper panels of Fig. 1). This spatial structure primarily reflects the climatological pattern of solar radiation reaching the surface in the Arctic \citep{lindsay_arctic}. Highest RS values are found where cloud cover and water vapor are low over the high topography of Greenland. Moderate RS values are seen in the central Arctic due to the thin but persistent cloud cover over the perennial sea ice. RS is smallest in the GIN seas due to abundant thick clouds. %consistent with the finding of \citet{donohoe_battisti_11a}. There is remarkable inter-model spread in Arctic RS across the different radiative kernel calculations, especially over the central Arctic where the models differ by a factor of \TO{two}. As shown below, the diversity of RS across the different kernel calculations is a consequence of inter-model differences in the mean state cloudiness and \emph{not} due to differences in radiative transfer code or the methodology used to calculate the kernels between the different groups. %Similar estimates of RS for the MPI ECHAM5 model \citep{previdi_10} are not shown because fields needed for the alternate RS calculation discussed in the next subsection were not available, but that model is also within the range shown in Fig. 1. %The spatial averages of RS are displayed in the corners of Fig. 1 and also differ by a factor of two across models over the whole Arctic (upper left), Arctic ocean (lower right) and sea ice domain (lower left) In the SO, RS during the Austral summer (NDJF) calculated from radiative kernels shows a zonally annular structure in all models with smaller values over the cloudy storm track region equatorward of the ice edge and larger values over the sea ice (upper panels of Fig. 2). However, the models differ to first order on the magnitude of RS over the open ocean and on the location and aerial extent of the region of larger RS adjacent to the Antarctic continent. In HADGEM2, the value of RS over the open ocean is 2 W m$^{-2}$ $\%^{-1}$ whereas in NCAR CAM3 RS is 1 W m$^{-2}$ $\%^{-1}$ over the same region. In NCAR CAM5, the region of high RS adjacent to the Antarctic coast extends substantially into the SO whereas in NCAR CAM3 and ECHAM6 the high RS region is confined to the coast itself with the exception of the Weddel and Ross Seas. The inter-model differences in the aerial extent of the high RS region roughly correspond to inter-model biases in summertime ice extent; the gradient in atmospheric tranmissivity is linked to the sea ice edge via cloud coverage and atmospheric water content although in some models the gradient in cloudiness is significantly poleward of the ice edge (i.e. NCAR CAM3) while in other models the cloud gradient is co-located with the ice-edge(i.e. NCAR CAM5). Overall, the SO domain average RS (excluding the Antarctic continent -- to focus on the sea ice) ranges from 1.29 to 1.75 W m$^{-2}$ $\%^{-1}$ (as shown by the values in the upper right corner of Fig. 2). \subsection{Isotropic single layer model} \citet{taylor_et_al_07} -- hereafter T07-- developed a model (hereafter the isotropic model) that can be used for approximating RS from the climatological radiative fluxes at the TOA and surface and some basic assumptions about shortwave radiative transfer in the atmosphere. Part of the T07 derivation is repeated here for clarity with a few modifications to variable names. Of the incident shortwave radiation at the TOA (S), assume a fraction (A) is absorbed in the atmosphere above cloud top and a fraction R of the radiation incident on cloud top is reflected back to space (Fig. 3). This resultant downwelling radiation at the surface is S(1-A)(1-R). A fraction ($\alpha$, equal to the surface albedo) of this downwelling radiation is reflected upwards. Of this surface upwelling radiation, R is reflected back (downward) to surface with the remainder (S[1-A][1-R]$^{2}$) transmitted to space. Reflections and transmissions are continued indefinitely subject to the three primary assumptions: (i) cloud optical properties can be represented by a single layer ,(ii) cloud reflection is isotropic -- the same fraction (R) of broadband shortwave radiation incident on the cloud layer is reflected independent of the direction (upwelling/downwelling) and how many previous interactions with the surface and cloud occur and (ii) all of the atmospheric absorption occurs above cloud top on the first downward pass which is apt for describing SW absorption by ozone in the stratosphere \citep{chou_and_lee}. \TO{We further analyze the limitations of these assumptions at the end of this subsection.} In the isotropic model, loss of shortwave radiative energy from the climate system due to surface albedo is a three step process: (i) insolation must be transmitted to the surface then (ii) reflected by the surface and finally (iii) transmitted from the surface to the TOA. Mathematically, upwelling SW radiation at the TOA that results from reflection off the surface is equal to the insolation (S) times the downwelling transmissivity ($\left[1-A\right]\left[1-R\right]$) times the upwelling transmissivity (1-R). The isotropic model also includes higher order reflections where the SW radiation reflected at the surface is reflected back to the surface off clouds and thereafter will contribute additional upwelling SW fluxes at the TOA with each subsequent reflection equal to the value of the previous order contribution times $\alpha$R. These terms form an infinite geometric series % \begin{equation} %\label{eq_iso} %SW \uparrow_{TOA} = SR\left(1-A\right) + \sum_{n=1}^{\infty} S \alpha^{n} R^{n-1} \left(1-A\right)\left(1-R\right)^2. %\end{equation} that converges to the expression: \begin{equation} \label{eq_iso} SW \uparrow_{TOA} = \underbrace{SR\left(1-A\right)}_{SW\uparrow_{TOA,atmos}} + \underbrace{S\alpha \frac{\left(1-A\right) \left(1-R\right)^{2}}{1-\alpha R}}_{SW\uparrow_{TOA,surf}}, \end{equation} \noindent where SW$\uparrow_{TOA, atmos}$ and SW$\uparrow_{TOA, surf}$ indicates the upwelling radiation at the TOA that was derived from atmospheric and surface reflection respectively. Thus, if the values of R and A along with $\alpha$ and S are known, the contribution of the surface to the SW flux at the TOA can be calculated. In our case, the isotropic model provides equations relating 3 satellite derived quantities ($SW \uparrow_{TOA}$, $SW \uparrow_{SURF}$ and, $SW \downarrow_{SURF}$) in terms of 3 unknown variables (A,R,$\alpha$) and the satellite measured S. The result is a determined set of 3 equations in terms of 3 variables. Thus, the climatological radiative fluxes allow the calculation of the single pass A and R for each climate model. We can then calculate the expected change of $SW \uparrow_{TOA}$ as $\alpha$ changes with all else being equal by taking the partial derivative of Eq. 3 with respect to $\alpha$. \begin{equation} \label{eq_RSiso} RS = \frac{\partial SW \uparrow_{TOA}}{\partial \alpha} = S\frac{\left(1-A\right) \left(1-R\right)^{2}}{1-\alpha R} \left(1+ \frac{R\alpha}{1-R\alpha}\right). \end{equation} \noindent This provides an alternate method for calculating RS that relies only on readily available model output at monthly resolution that can also be compared with the RS calculated from radiative kernel techniques. %In this framework, RS consists of two terms. The first term on the RHS of Eq. \ref{eq_RSiso} is the climatological ratio of the surface contribution to planetary albedo to the surface albedo. The second term includes the impact of the modified $\alpha$ on secondary and successive reflections and tends to be substantially smaller than the first term except over regions that have both very bright clouds and very bright surfaces. Specifically, the second term becomes comparable in magnitude to the first term when $\frac{R\alpha}{1-R\alpha} \approx 1$ that occurs for R and $\alpha$ both $>$ 0.7 which is not realized in the observed Arctic (Fig. 5) . The lower panels of Fig. 1 show the RS in the Arctic summer calculated from Eq. \ref{eq_RSiso} applied to the monthly climatological output from the same control simulations that were used to calculate the radiative kernels. % We note that climatological fields from the GFDL simulation used to calculate the radiative kernels were not saved and, thus, we cannot use the isotropic model to calculate RS in the GFDL climate model. The RS calculated from the isotropic model is in good agreement with that calculated from radiative kernels in terms of the spatial pattern of RS and inter-model differences. Spatial correlation between RS in the isotropic model and radiative kernel calculation for each model is high with an R$^{2}$ that exceeds 95$\%$ in all but NCAR CAM3. The inter-model differences in domain average of RS is within 10$\%$ in the absolute sense and captures the rank of RS in models (c.f. the adjacent upper and lower panels of Fig. 1 with R$^{2}$ listed in the middle). The isotropic model explains 94$\%$ of the variance in MJJA RS calculated from radiative kernels considered across models and over all Arctic grid points collectively with a root mean squared (RMS) error of 0.15 W m$^{-2}$ \TO{$\%^{-1}$} (top panel of supplemental Fig. A2). As a basis for comparison, if one used the spatial pattern of MJJA RS calculated using radiative kernels from one model to predict the kernel based RS in a different model -- as was done in the IPCC estimate of SIAF-- the RS variance explained is 21$\%$ with a RMS error of 0.67 W m$^{-2}$ \TO{$\%^{-1}$} (bottom panel of Fig. A2). Thus, the isotropic model offers a factor of 4 improvement on the practice of applying RS calculations from a single climate model. The isotropic model also captures the spatial pattern and inter-model spread of the kernel calculated RS in the SO (Fig. 2) although the absolute values of RS differs by as much as 20$\%$ (in the HADGEM2 model). The isotropic model explains 96$\%$ of the variance in NDJF kernel RS across models over the SO (top panel of supplemental Fig. A3) with a root mean squared (RMS) error of 0.23 W m$^{-2}$ \TO{$\%^{-1}$}. When radiative kernels from one model are used to predict the kernel based NDJF RS in a different model the variance explained is 71$\%$ with a RMS error of 0.47 W m$^{-2}$ \TO{$\%^{-1}$} (bottom panel of Fig. A3). Thus, the isotropic model offers a factor of 2 improvement on the practice of applying RS calculations from a single climate model in the SO. These results indicate that the isotropic model captures the essential SW radiative processes that determine the RS of surface albedo changes, and that the inter-model spread in RS is determined by the climatological cloud reflectivity which is adequately calculated from the modeled TOA and surface fluxes according to Eq. \ref{eq_iso}. The isotropic model tends to bias the RS high relative to the radiative kernel (c.f. the domain average values listed in the upper right of the map in the upper and lower panels of Figs. 1 and 2) and we speculate this results from the simplifying assumption that atmospheric absorption only occurs during the first pass as this allows more of the radiation reflected off the surface to be transmitted to space than would occur if the atmosphere absorbed upwelling solar radiation. Alternative formulations of similar isotropic models \citep{donohoe_battisti_11a} assume the atmospheric absorption occurs in the same layer as the cloud reflection and occurs on all passes through the atmosphere to account for shortwave absorption by water vapor that occurs throughout the troposphere \citep{donohoe_battisti_seas}. This model better matches the RS calculated by radiative kernels in the tropics and mid-latitudes but substantially underestimates RS relative to the radiative kernel derived value at high latitudes (Appendix Fig. A1). We speculate that in the dry Arctic, the atmospheric absorption is primarily by stratospheric ozone whereas in the lower latitudes water vapor also contributes. For this reason, we choose to assume the absorption occurs only on the downward pass and return to possible impacts and improvements of this method in the discussion section. \subsection{Causes of inter-model spread in RS} What processes are responsible for the factor of 2 spread in modeled RS in Figs. 1 and 2? The ability of the isotropic model to reproduce the kernel based RS calculated for each model demonstrates that the mean state atmospheric opacity is the primary determinant. Generally speaking, RS is determined by how much insolation is transmitted to the surface and thus how much impact surface albedo changes have on reflected solar radiation. More specifically, RS is proportional to the atmospheric transmissivity squared with higher order modifications due to the impact of multiple reflections (Eq. \ref{eq_RSiso}). What then causes the inter-model spread in atmospheric opacity? Clear sky surface albedo kernels (Fig. 4) have much larger magnitudes than their all-sky counterparts. The very similar spatial structures and absolute values in the four models with available kernel calculations have domain averages that differ by $~2\%$ from the multi-model mean, indicating that 1) clear-sky processes are not responsible for the inter-model spread in all-sky RS and 2) the different radiative transfer codes used in the climate models find a similar RS for a similar (clear-sky) mean-state. %The minimal role of differences in the radiative transfer code in the inter-model spread of RS is also supported by the large differences in full-sky RS between NCAR CAM3 and NCAR CAM5 which use the same radiative transfer model but have very different RS. The atmospheric opacity parameters -- reflectivity and absorptivity -- \TO{calculated by the} isotropic model applied to the mean states of the different climate models are shown in Fig. 5. The all-sky reflectivity is subdivided into a clear-sky and cloud component by applying the isotropic model to the clear-sky mean state radiative fields (as in T07) to define a clear-sky reflectivity, and the cloud reflectivity is then defined as the all-sky minus clear-sky reflectivity. All climate models have very similar and nearly spatially uniform clear-sky reflectivity and all-sky absorptivity with Arctic domain average absolute differences from the model mean of order 0.02 fractional units. The slight spatial structure in clear-sky reflectivity and absorptivity is consistent between climate models. Clear-sky reflectivity is larger near the North Pole consistent with enhanced Rayleigh scattering due to the shallower angle of incidence with latitude. Absorptivity is smaller over the thinner atmosphere above topography and drier continents consistent with reduced absorption by water vapor. In contrast to the consistency of absorption and clear-sky reflection between models, the cloud reflectivity differs substantially between models in both spatial structure and domain average values (which differ between models by over 0.20 fractional units). In general, regions of stronger cloud reflectivity have smaller RS values consistent with less downwelling solar radiation at the surface in cloudy regions. However, the anti-correlation between the spatial variability in RS and cloud reflectivity is significant but far from perfect (R $\approx$ -.60) within a given climate model due to the (comparable in magnitude) impact of the spatial structure of mean state albedo (Eq. \ref{eq_RSiso}) on the multiple reflection contribution to RS. On a broader scale, the Arctic domain average cloud reflectivity is very strongly anti-correlated (R=-0.99) with the domain average RS indicating that Arctic averaged RS is primarily determined by the mean state cloud reflectivity. \section{Observational estimate of radiative sensitivity to surface albedo changes and comparison to coupled models} Given the strong correspondence between RS calculated from radiative kernels and the isotropic model (Figs. 1 and 2), we can use the isotropic model to calculate RS from observational estimates of radiative fluxes at the TOA and surface and use these same fields (routinely available from model simulations) to assess model biases in RS and diagnose their role in the SIAF. Observational estimates of climatological radiative fluxes are taken from the CERES EBAF surface product version 4.0 \citep{ceres_4_toa, ceres_4_surf} between 2000 and 2018. Climate model RS is estimated using the isotropic model for the last decade (1995-2005) of historical CMIP5 \citep{cmip5} climate simulations forced. \footnote{Most of the radiative kernel calculations discussed in Section 2 used ``modern'', slightly differing time periods.} Maps of summer (MJJA) RS estimated from satellite products and models are shown in Fig. 6. Three spatial averages of RS are also provided: (i) the whole domain poleward of 60$^{\circ}$N (upper left corner in black) -- observational value of 1.79 W m$^{-2}$ $\%^{-1}$; (ii) the Arctic ocean excluding land masses (lower right in blue) -- observational value of 1.68 W m$^{-2}$ $\%^{-1}$ and; (iii) the spatial average over the sea ice (the spatial footprint and region varies between models, lower left in purple) -- observational value of 1.92 W m$^{-2}$ $\%^{-1}$. The observational RS is very similar to multi-model mean values ( 1.72, 1.65 and, 1.79 W m$^{-2}$ $\%^{-1}$ over the entire Arctic domain, Arctic ocean and sea ice regions respectively). The models and observations generally agree on the spatial pattern of RS over the Arctic with high values over the Greenland ice sheet where the reduced mass of the atmosphere above the high topography is associated with enhanced atmospheric SW tranmissivity, lower RS values over the GIN Sea and more spatially uniform RS values over the Central Arctic. The magnitude of RS differs substantially across models with domain average RS varying by almost a factor of two between the models, consistent with results from the radiative kernel based RS calculation (Fig. 1). The inter-model (2$\sigma$) spread in Arctic average RS is 0.57, 0.53 and, 0.64 W m$^{-2}$ \TO{$\%^{-1}$} over the full Arctic domain, Arctic ocean and climatological sea ice respectively. The SO observational estimate of summertime (NDJF) RS is similar but slightly lower (domain average excluding the Antarctic continent of 1.56 W m$^{-2}$ $\%^{-1}$) than the multi-model mean (1.71 W m$^{-2}$ $\%^{-1}$). All models and observations show an annular structure in RS with smaller values in the storm track region and larger values adjacent to the Antarctic continent over the sea ice (Fig. 7). RS differs substantially between models (order factor 2) in the storm track region and on the location and lateral extent of the high RS region adjacent to the continent. Some models (i.e. CSIRO MK5) also have zonal asymmetries in RS that are best characterized as a zonal wavenumber 1 pattern. The domain average RS values differ by less than the factor of 2 differences seen in the Arctic, but the local RS difference between models --especially in the storm track region -- are of order a factor of 2. The inter-model (2$\sigma$) spread in SO average RS is 0.54 W m$^{-2}$ \TO{$\%^{-1}$}, comparable in magnitude to that over the Arctic domain and Arctic ocean. These results collectively suggest that while CMIP5 ensemble average RS of high latitude ice loss is quite similar to that implied from observational constraints, models diverge substantially on the radiative impact of ice loss because of differences in atmospheric optical properties (i.e. clouds). %These inter-model differences in RS impact SIAF and the response to anthropogenic forcing. \section{Observational estimate of ice albedo feedback} The \TO{Arctic sea-}ice-albedo feedback (SIAF) is (the spatial average of) the product of the RS -- the TOA radiative impact of surface albedo changes -- and the ice sensitivity (IS) -- the surface albedo change due to \TO{Arctic sea} ice loss per unit of global warming (Eqs. 1,2). Thus, the RS calculated from the climatological radiative fluxes and the isotropic model in the previous sections along with estimates of IS from the observational record provide an observational estimate of the SIAF that can be compared to the SIAF calculated using the same methodology applied to CMIP5 simulations with historical and long term forcing. Furthermore, we can explicitly ask if the model spread (and potential bias relative to observations) in SIAF is explained by RS or IS spread. The observational estimate of IS is calculated from the changes in decadal surface albedo of the Arctic ocean from 1982 to 2016 (2007-2016 average minus 1982-1991 average -- Fig. 8) \TO{during each summer month} divided by the global mean surface temperature (T$_{S}$) change over the same time period. \TO{We use two different observational based data sets to calculate the change in surface albedo over this time period} : (i) sea ice concentration calculated by the National Snow and Ice Data Center \citep{nsidc} from passive microwave brightness measured by the Nimbus 7 satellite available from 1979-2016 and (ii) broadband \TO{(all-sky)} surface albedo measured by the Advanced Very High Resolution Radiometer (AVHRR) Polar Pathfinder (APP-x) extended data set \citep{appx} that covers the 1982-2017 time period. \TO{The central estimate of our observational based IS is the average of calculations from these two data sets (elaborated on below) and our uncertainty estimates account for differences across the two data-sets.} The NSIDC sea-ice concentration changes are converted to a surface albedo change record by multiplying the SIC changes by the albedo contrast between sea ice and open ocean ($\Delta \alpha$), which is assumed to be spatially and temporally invariant: \begin{equation} \label{eq_IS} IS = \frac{dSIC}{\left[dT_{S}\right]} \Delta{\alpha}. \end{equation} \noindent Eq. \ref{eq_IS} assumes that changes in $\alpha_{SI}$ are isolated to regions of sea-ice melt. NSIDC monthly maps of the decadal average change in sea ice concentration are multiplied by an assumed surface albedo contrast between the open ocean and sea ice ($\Delta \alpha$) of 0.54 -assuming a typical ice $\alpha$ of 0.6 \citep{hummel_reck} and an ocean albedo of 0.06 \citep{hansen_83}. This choice of typical ice albedo is an average of snow covered sea ice found during the late spring and sea ice with melt ponds in the late summer \citep[see Fig. 9 of][]{perovich_02}. This map of monthly NSIDC ice-concentration derived surface albedo change and those derived from the APP-x \TO{(also monthly)} data are averaged to produce the observational best estimate of change in surface albedo (Fig. 8C) -- hereafter referred to as the Observational Best Estimate (OBE). Both products produce similar estimates of surface albedo changes (Appendix Fig. A4). We use differences between the two surface albedo data sets as well as the intra-decadal variability within each data set to calculate the uncertainty in observational IS (Fig. 8D) as outlined in the Appendix. Observational IS is calculated by normalizing OBE surface albedo changes by a global surface temperature change of 0.7 $\pm$ 0.1K over the 1982-2016 time period. The central estimate and uncertainty in global mean surface temperature change come from the average and standard deviation of the mean across three different global surface temperature data sets: (i) National Centers for Environment Prediction (NCEP) reanalysis surface air temperature \citep{kalnay_et_al}, (ii) the Goddard Institute for Space Studies Surface Temperature Analysis (GISTEMP) \citep{gistemp}, and (iii) the modification by Cowtan and Way \citep{cowtan_way} of the Met Office Hadley Centre surface temperature dataset \citep{hadcrut} version 4 (HadCRUT4). The monthly IS is then multiplied by the monthly RS derived from CERES data, and then time averaged (over the summer months) to produce a map of radiative impact of sea ice changes (Fig. 8E). While the previous figures showed MJJA average in the NH and NDJF in the SH, figure 8E extends the summertime season to include the six months centered on the summer solstice (AMJJAS in NH and ONDJFM in the SH) since previous work \citep{flanner_albedo} found an appreciable contribution to the SIAF during the shoulder seasons especially in April. The uncertainty in the $RI_{TOA,\alpha}$ is (Fig. 8F) is assessed from a Monte Carlo simulation that takes into account 3 different uncertainties in the input data sets propagated (in quadrature) onto the calculation of $RI_{TOA,\alpha}$ (see Appendix for details): (i) the uncertainty in RS -- due to uncertainty in the climatological radiative fluxes, (ii) uncertainty in the surface albedo change -- due to both intra-decadal variability and differences between the APP-x and NSIDC ice concentrations data sets and (iii) uncertainty in the global mean temperature change that goes into the calculation of IS. We note that $RI_{TOA,\alpha}$ in Fig. 8E is, by definition, the radiative impact of sea ice changes normalized by global mean surface temperature change (=0.7K) and has a summertime (AMJJAS) Arctic domain average of 4.9 $\pm$ 1.4 W m$^{-2}$ K$^{-1}$ which translates to an absolute change in summertime radiation of 3.4 $\pm$ 1.0 W m$^{-2}$ over the Arctic. To convert this number to a global and annual mean radiative impact, one must weight this number by the ratio of summer months to the year ($\frac{6}{12}$) and the spatial area of the Arctic (poleward of 60$^{\circ}$N) divided by that of the globe (.065) resulting in a global TOA radiative change of 0.11 W m$^{-2}$ over the 1982-2016 period. This translates to a global radiative feedback \TO{(divide by 0.7 K global T$_S$ change)} of 0.16 $\pm$ 0.04 W m$^{-2}$ K$^{-1}$ given the observed global surface temperature change over the same period. The uncertainties cited above reflect 2 standard deviations. We do not estimate the observational based surface albedo feedback in the SO because the change in SO sea ice concentration over the observational period is not statistically significant above the year-to-year variability \citep{SO_trends}. We also note that this estimate is isolated to the Arctic ocean (we have masked the APP-x albedo changes over land) and, thus, does not include the impact of changes in snow cover over land. \section{Comparison of observational and model SIAF and decomposition of inter-model spread of SIAF into RS and IS} We now compare the observational Arctic SIAF derived above with that derived by the same methodology in historical CMIP5 simulations. The RS for each climate model that was calculated using the isotropic model in the previous section (from the climatology at end of the historical simulation -- 1995 to 2005) is multiplied by the decadal average surface albedo change --calculated as the ratio of upwelling to downwelling broadband shortwave radiation at the surface -- over the historical simulation (1995 to 2005 minus 1975 to 1985). We note that this time period was chosen to correspond to the end of the historical simulations and differs from the 1982 to 2016 period used for the observational calculations. The RS and surface albedo changes are calculated for each month and the product is spatially averaged over the Arctic ocean to calculate the SIAF; we exclude the impact of changes in snow cover over land from our calculations. For simplicity, we will only discuss the annual and global mean of the calculations normalized by the global mean surface temperature change over the same time period, as we did for the observations. The CMIP5 ensemble mean Arctic SIAF in the historical simulations is 0.12 W m$^{-2}$ K$^{-1}$ with a spread (2 standard deviations, $\sigma$) of 0.13 W m$^{-2}$ K$^{-1}$ (gray histogram in Fig. 9A with wide bars). The ensemble mean is slightly smaller than the observational estimate (c.f the solid and dashed vertical black lines in Fig. 9) but the large inter-model spread indicates that the models differ in either RS and/or IS. We now ask how much RS and IS contribute to the inter-model differences in Arctic SIAF. To estimate the IS contribution to the SIAF spread, the calculation of SIAF is repeated but the model specific RS is replaced with the observational based RS value. The resulting distribution of SIAF (blue histogram in Fig. 9A) shows the spread produced by biases and inter-model differences in IS. The mean value of SIAF in the fixed RS distribution (0.12 W m$^{-2}$ K$^{-1}$ -- Table 1) is nearly equal to that of the full SIAF calculation (c.f. the blue and black vertical lines). The CMIP5 ensemble average SIAF is quite insensitive to RS model biases, and it is lower than the observed estimate because the modeled IS is smaller than the observational estimate. Furthermore, the spread in the fixed RS distribution is only slightly smaller than that of the full SIAF calculation (2$\sigma$ = 0.12 W m$^{-2}$ K$^{-1}$) indicating that the majority of the inter-model spread in SIAF calculated from the historical simulation is a result of the IS differences between models. A similar analysis can be made to estimate the impact of biases (relative to observations) and inter-model RS differences on the calculated SIAF by replacing the model specific IS with that derived from observations (red histogram in Fig. 9A). The CMIP5 ensemble average SIAF of the fixed IS distribution (0.16 W m$^{-2}$ K$^{-1}$ -- Table 1) is larger than that of the full SIAF calculation (c.f. the red and black vertical lines in Fig. 9) indicating that the CMIP5 ensemble average IS is smaller than that observed (the OBE value) a result also found by \citet{rosenblum}. The inter-model spread in SIAF in the fixed IS experiment (2$\sigma$ = 0.04 W m$^{-2}$ K$^{-1}$) is smaller than that of the full calculation and fixed RS experiment indicating that inter-model differences in RS play a smaller but not insignificant role in the SIAF spread calculated over the historical simulations. A summary of the role of biases and inter-model differences in RS and IS in determining the model distribution of SIAF is provided in Table 1. This partitioning of SIAF differences in contributions from RS and IS takes spatial and temporal co-variances of ice loss and RS into account by weighting the ice loss to the RS at that location and time. Similar results for the impact of IS and RS on the total spread in SIAF are obtained by simply noting the fractional spread (relative to the ensemble mean) of summertime Arctic domain average RS and IS between models. The ratio of domain and summertime average inter-model spread (2$\sigma$) to the ensemble mean domain and summertime average of RS is 40$\%$ whereas that of IS is 107$\%$ roughly scaling with the fractional contribution to SIAF spread calculated above. This result suggests that inter-model differences in IS and RS are fairly spatially and temporally homogenous and the resultant inter-model spread in SIAF is independent of the spatiotemporal co-variability of RS and IS. Previous work has found similar large magnitude inter-model spread in IS in CMIP3 \citep{cmip3_ice_sensitivity} and CMIP5 \citep{cmip5_ice_sensitivity} linked to the spread in the magnitude of Arctic amplification. The sea ice retreat over the historical record represents the superposition of the response to climate forcing and natural variability and, thus, the inter-model spread in IS calculated over the 30 years of historical simulations is expected to exceed that in response to long-term sustained forcing. \TO{\citet{schneider_albedo} found that decadal trends in sea-ice during periods when global mean temperatures increased by more than 0.5K provided good estimates of the long-term SIAF in an ensemble of climate models. Other studies suggest that as much as 50$\%$ of the observed Arctic sea loss since 1979 could be a result of the natural variability of atmospheric circulation \citep{ding_ncc, ding_cesm,kay_trend}}. To reduce the amount of internal variability relative to the forced component, we also look at the contribution of RS and IS to the inter-model spread in SIAF in response to an abrupt and sustained quadrupling of atmospheric CO$_{2}$ where the forced climate change signal is expected to be larger than the natural variability. The IS in the CO$_{2}$ quadrupling simulations is calculated from the change in surface albedo and global mean surface temperature between the PI and the average over years 50-100 since CO$_{2}$ quadrupling. The RS used to calculate the \TO{SIAF} is calculated from the PI climatological fields in the same model. The ensemble average Arctic SIAF calculated from the 4XCO$_{2}$ simulations is 0.13 $\pm$ 0.09 W m$^{-2}$ K$^{-1}$ (uncertainty is 2$\sigma$) and is in close agreement with the ensemble average of the historical simulation (0.12 $\pm$ 0.13 W m$^{-2}$ K$^{-1}$) with reduced inter-model spread. The central estimate and range of SIAF from all model simulations -- calculated from 2$\sigma$ of the mean -- is 0.13$\pm$ 0.02 W m$^{-2}$ K$^{-1}$ and is slightly smaller than but not statistically different from the observational estimate (0.16 $\pm$ 0.04 W m$^{-2}$ K$^{-1}$). Because the 4XCO$_{2}$ ice response primarily reflects the forced response, the similarity of the ensemble average SIAF diagnosed from historical and 4XCO$_{2}$ simulations suggests that the same physics responsible for the long-term SIAF are evident in historical simulations despite the additional statistical noise from internal variability. When the model specific RS is replaced by the observational estimate of RS the resultant Arctic SIAF for the CO$_{2}$ quadrupling simulations is 0.13 $\pm$ 0.08 W m$^{-2}$ K$^{-1}$ and when the model specific IS is replaced by the observational estimate of IS the resultant SIAF is 0.16 $\pm$ 0.04 W m$^{-2}$ K$^{-1}$ (lower left panel of Fig. 9 and Table 1). These results suggest that in the long-term response to sustained anthropogenic forcing: 1.) the CMIP5 ensemble average RS (spatially and temporally weighted by the relevant regions of ice loss) is very near the observational estimate, 2.) the CMIP5 ensemble average IS (spatially and temporally weighted by structure of RS) is slightly smaller than the observational estimate and is responsible for the model SIAF being smaller than the observational estimate and 3.) inter-model differences in IS contribute twice as much to the inter-model spread in SIAF (63$\%$ of the ensemble average value) than do inter-model differences is RS (30$\%$ of the ensemble average value). We note that, the inter-model spread in IS and RS are significantly (R=0.54) correlated (at 95$\%$ confidence interval) and return to the implication of this result in the discussion section. A similar analysis can be performed for the 4XCO$_{2}$ simulations in the SO (poleward of 55$^{\circ}$S) to indicate an ensemble average SIAF of 0.08 $\pm$ 0.13 W m$^{-2}$ K$^{-1}$ (Fig. 9C, Table 1). The SO SIAF is negative in a single model (GFDL ESM2G) that simulates sea ice growth in the Weddell Sea under 4XCO$_{2}$. When the model specific RS is replaced by the observational estimate of RS, the calculated SO SIAF is 0.07 $\pm$ 0.11 W m$^{-2}$ K$^{-1}$ suggesting the the ensemble average RS is slightly larger than that estimated from the observations, consistent with Fig. 7. Because no observational estimate of SO IS is available, we probe the sensitivity of SO SIAF to RS by replacing the model specific IS with the ensemble average IS resulting in a calculated SIAF of 0.08 $\pm$ 0.04 W m$^{-2}$ K$^{-1}$; inter-model differences in RS result in inter-model differences in SO SIAF of magnitude 50$\%$ the ensemble mean estimate. However, the contribution of inter-model spread in RS to SIAF spread is dwarfed by the impact of inter-model differences in IS which produces inter-model differences in SO SIAF exceeding the central estimate by almost a factor of 1.5 (160$\%$). This result is consistent with the large inter-model differences in SO ice response to global warming reported by \citet{cmip_ice, polvani_SO}. \subsection{Global surface albedo feedback -- comparison to IPCC AR5 value} The IPCC AR5 estimated a global surface albedo feedback of 0.26 W m$^{-2}$ K$^{-1}$ based on the calculations of \citet{soden_kernel} which use a single RS -- derived from kernel calculations in the GFDL model (Fig. 1) -- applied the surface albedo change in each CMIP3 model. These calculations are global and include the impact of $\alpha$ changes over land (due to changes in snow cover) in addition to the sea-ice related changes considered up to this point and we term this combined contribution of land and sea ice changes the global albedo feedback (GAF).\TO{More recently, \citet{schneider_albedo} presented a CMIP5 ensemble mean \TO{GAF} 0.40 W m$^{-2}$ K$^{-1}$ using NCAR CAM5 based kernels. It is unclear if this discrepancy \TO{results} from the different RS used in these studies or the IS in different GCM ensembles.} Here, we compare the GAF produced using the (kernel based) RS from a single model to that calculated using a model specific RS derived from the isotropic model. Our \TO{GAF} calculations are based upon surface albedo change calculated from the 4XCO$_{2}$ simulations minus that in the pre-industrial simulation normalized by the global mean surface temperature (TS) change in that model -- a quantity akin to IS in Eq. 2 but including the albedo changes over land. This albedo change is multiplied by RS estimated two ways: (i) using the method introduced in this study, where RS is calculated from the isotropic model (Eq. \ref{eq_RSiso}) using radiative fluxes from appropriate model specific PI simulation and (ii) using the method introduced by \citet{soden_kernel} where the GFDL surface albedo kernel (Fig. 1) is used to estimate RS for all models. We separate the \TO{GAF} calculation into hemispheres. In the NH, the GAF calculated in this study is larger than that calculated using the GFDL kernel in all models (all the red dots fall below the 1:1 line in the upper left panel of Fig. 10) as would be expected from the GFDL RS being at the very low end of the model range especially over the ocean domain. In the CMIP5 ensemble average, the NH \TO{GAF} is 0.20 W m$^{-2}$ K$^{-1}$ using the GFDL kernel as compared to 0.27 W m$^{-2}$ K$^{-1}$ using the isotropic model methodology (35$\%$ greater -- Table 2). The NH \TO{GAF} has 64$\%$ more spread using the model specific RS because: (i) the ensemble mean RS is larger than the GFDL kernel RS and (ii) the inter-model spread in RS contributes to the \TO{GAF} spread as discussed in the previous subsection. If we restrict the calculation to the Arctic ocean poleward of 60$^{\circ}$N (as was done in Sections 4 and 5) we find a CMIP5 ensemble average SIAF of 0.09 W m$^{-2}$ K$^{-1}$ using the GFDL kernel compared to the 0.13 W m$^{-2}$ K$^{-1}$ (Table 1) using the isotropic model methodology (45$\%$ greater). This result suggests that approximately half of the \TO{GAF} is due to $\alpha$ changes over land as found by \citet{flanner_albedo}. In the Southern Hemisphere, the \TO{GAF} estimates from the two methods are in closer agreement; the dots cluster along near the 1:1 line in the upper right panel of Fig. 10 with the exception of the models producing the highest \TO{GAF}. This result is expected since the GFDL RS is near the ensemble mean over the SO (Fig. 2 and Fig. 7). The ensemble average \TO{GAF} in the SH is, therefore, very similar when using the methodology in this study (0.09 W m$^{-2}$ K$^{-1}$) as compared to that calculated using GFDL RS only (0.08 W m$^{-2}$ K$^{-1}$) -- Table 2. Globally, we calculate an \TO{GAF} of 0.37 W m$^{-2}$ K$^{-1}$ which is 30$\%$ greater than the same result found applying the GFDL RS to CMIP5 4XCO$_{2}$ simulations of 0.29 W m$^{-2}$ K$^{-1}$. We note the the IPCC AR5 cites a global \TO{GAF} of 0.26 W m$^{-2}$ K$^{-1}$ derived from the GFDL kernel and CMIP4 simulations and, thus, our estimate is 40$\%$ larger than the AR5 value. We attribute 30$\%$ of this increase to improved methodology of using model specific RS and 10$\%$ to the difference between CMIP4 and CMIP5 model characteristics. \TO{Importantly, the IPCC diagnosis of the overall climate sensitivity of climate models is unaffected by our revised more positive \TO{GAF}. Rather, our results suggest that the shortwave cloud feedback should be revised downward by the same amount because cloud feedbacks are diagnosed from all-sky minus clear-sky TOA radiation adjusted by all-sky minus clear-sky radiative kernel calculations.} \section{Summary and discussion} We have shown that the radiative impact of surface albedo changes (RS) calculated using offline radiative transfer models (radiative kernels) can be closely replicated using a single layer isotropic SW radiation model applied to the climatological radiative fluxes at the TOA and surface. This procedure allows estimates of SIAF to be conveniently calculated from observational data sets and standard model output without use of a kernel calculation, facilitating a comparison of observational and model estimates of SIAF. It also allows the differences between models and observations based calculations to be decomposed into contributions from RS and IS. The multi-model mean of RS is close to the observational estimate in the Arctic and only slightly larger than the observational estimate in the SO. However, the inter-model spread in RS (Figs. 6 and 7) is substantial, producing inter-model differences in SIAF estimates that are 30$\%$ and 50$\%$ the magnitude of the ensemble mean SIAF in the Arctic and SO respectively. In agreement with \citet{sledd_albedo}, high latitude clouds tend to mask the impact of surface albedo variations on the TOA albedo by a factor of 2-3 in observational estimates. Differences in climate model clouds influence the degree of cloud masking. Our results indicate that inter-model differences in IS are more important than RS in explaining the inter-model spread in SIAF. However, IS is not statistically independent of RS (R = 0.54). \TO{It is possible that inter-model differences in RS contribute to inter-model difference in IS because models that have a larger radiative response to sea ice loss will tend to have greater sea ice loss due to a stronger positive feedback between initial ice loss and radiative heating.} In this sense, the contribution of RS to inter-model differences in SIAF of 0.04 W m$^{-2}$ K$^{-1}$ \TO{both in the Arctic and SO} can be thought of as a lower bound on the contribution of mean state radiative biases to the SIAF. We hope to explore the impact of mean state radiative biases (RS) on IS and the persistence of sea ice loss events in future work. We estimate an observationally based global, and annually averaged increase in TOA radiation of 0.11 W m$^{-2}$ from Arctic sea ice changes over the 1982-2016 time period using observationally based estimates of sea ice changes and the CERES derived radiative sensitivity (RS) implying a SIAF of 0.16 $\pm$ 0.04 W m$^{-2}$ K$^{-1}$. \citet{flanner_albedo} found a Northern Hemisphere average ''crypospheric radiative forcing" of 0.45 W m$^{-2}$ over the 1979-2008 time period about half of which (0.22 W m$^{-2}$) was attributed to sea ice changes -- the other half was attributed to snow changes over land. Thus, the \citet{flanner_albedo} result converted to a global average (0.22/2 =0.11 W m$^{-2}$) agrees very well with our findings. Similarly, \citep{cao_albedo} found a Northern Hemisphere SIAF of 0.25 W m$^{-2}$ K$^{-1}$ using observed surface albedo change and RS estimated using model based kernels derived from GFDL \citep{soden_kernel} and CAM3 \citep{shell_kernel}. This result translates to a global feedback of Arctic changes of 0.12 W m$^{-2}$ K$^{-1}$ which is smaller than our central estimate and we speculate this result follows from the lower than observed RS in the CAM3 kernel (Fig. 1). \citet{pistone,pistone2} calculated a substantially larger SIAF (0.31 $\pm$ 0.04 W m$^{-2}$ K$^{-1}$) from the inter-annual covariance of sea ice concentration and TOA radiation measured by CERES. We speculate that some of the TOA radiative variability that coincides with ice loss events in \citet{pistone} is not directly a consequence of (i.e. geographically co-located with and/or a radiative consequence) surface albedo changes but, rather, is a consequence of atmospheric optical properties (i.e. clouds, water vapor, etc) that co-vary with Arctic sea ice concentration. A central question moving forward is whether the atmospheric changes (and the associated radiative anomalies) accompanying Arctic sea ice loss over the limited historical period result from natural variability of atmospheric circulation initiated by tropical and mid-latitude processes or are a direct result of sea ice loss and, thus, should be expected to also apply to future climatological changes. Additionally, how accurately does the observational IS calculated over the historic record represent the expected relationship between future changes in Arctic ice concentration and global mean temperature? \citet{pistone} suggest that the SIAF \TO{(Arctic ocean only)} alone results in a 25$\%$ enhancement of global warming via radiative feedbacks, a value they derive from the ratio of their calculated radiative impact of historic ice loss divided by the anthropogenic climate forcing to date. We offer two modifications as updates to their calculation: (i) a significantly lower estimate of the radiative impact of Arctic sea ice loss outlined above and (ii) consideration of how the implied feedback relates to equilibrium climate sensitivity, noting that the climate system is not currently in equilibrium with the anthropogenic forcing to date. For the latter reason, the feedback gain of the Arctic SIAF should be calculated by comparing the SIAF to the equilibrium radiative feedback of all other radiative processes as opposed to the ratio of the transient radiative impact of ice loss to date to the applied forcing. Given observational central estimates of the total equilibrium feedback parameter of -1.19 W m$^{-2}$ K$^{-1}$ \citep{kyle_lambda} and our observational estimate of the Arctic SIAF ($\lambda_{SIAF} =$ +0.16 $\pm$ 0.04 W m$^{-2}$ K$^{-1}$) the implied feedback parameter of all processes excluding the SIAF ($\lambda_{0}$) satisfies the equation -1.19 W m$^{-2}$ K$^{-1}$ = $\lambda_{0}$ + 0.16 W m$^{-2}$ K$^{-1}$. This implies that $\lambda_{0}$ (the reference climate feedback parameter of a system with no SIAF) is -1.35 W m$^{-2}$ K$^{-1}$. We note that the reference climate feedback parameter is more negative than that of a system with a SIAF implying a smaller climate sensitivity of the reference system relative \TO{to} the full system with a SIAF as is expected for the positive SIAF. The fractional amplification of global mean temperature changes -- the feedback gain, G$_{SIAF}$-- due to the SIAF is then \citep{roe_feedbacks}: \begin{equation} G_{SIAF} = \frac{1}{1+\frac{\lambda_{SIAF}}{\lambda_{0}}} = 1.14 \pm .04 . \end{equation} Thus, our analysis suggests that the Arctic SIAF amplifies global warming by 14$\%$ (2$\sigma$ range between 10 and 19$\%$) at the equilibrium timescale and is a more modest amplifier of global warming than the 25$\%$ suggested by \citet{pistone}. The IPCC AR5 report \citep{ipcc_models} points out a discrepancy between the observational based SIAF of \citet{flanner_albedo} and the model based estimate of \citet{soden_kernel} and speculates that models are biased toward low IS, but the role of inter-model spread and biases in RS were neglected. While we find no ensemble mean model bias in Arctic RS (Fig. 6), the model estimate of RS used in \citet{soden_kernel} is taken from radiative kernel calculations in a single (GFDL) model and then applied to the IS across models. The RS from that model (Fig. 1) is biased low relative to both the observational based RS (by 46$\%$ of the kernel RS in the Arctic average) and the CMIP5 ensemble mean. As a result, the AR5 estimate of the global surface albedo feedback of 0.26 W m$^{-2}$ K$^{-1}$ based on the calculations of \citet{soden_kernel} is substantially lower than our calculated value of 0.37 W m$^{-2}$ K$^{-1}$ which uses model specific RS estimates. This result suggests that at least some part of the low model bias identified in the IPCC AR5 is a consequence of using a RS that is inconsistent with some climate models. We recommend using model specific RS derived from the isotropic model as a better practice to applying radiative kernels across models. Additionally, our results identified no discernible model bias in the SIAF at least when considering like quantities over the Arctic ocean domain. \acknowledgments We thank Karen Shell and 3 anonymous reviewers for thoughtful critique of an earlier version of this manuscript. We thank Angeline Pendergrass, Chris Smith, Ryan Kramer, Karen Shell, Brian Soden, Block and Marotzke and Michael Previdi for providing radiative kernel calculation and their assistance providing further clarification on the appropriate climatological radiation fields to use for isotropic model calculations. We also thank Jeff Key and Xuanji Wang for providing the NOAA AVHRR Polar Pathfinder data. This work was funded by a Department of energy mini grant to the HLES team at PNNL, the NSF Antarctic Program Grant Number PLR 1643436 and the NOAA MAPP grant eGC1$\#$A127135. \appendix %\section{Supporting Information} We describe the methodology used to calculate the uncertainty in our observational estimates of the RS, IS and RI$_{TOA, \alpha}$ the spatial average of which gives the resultant SIAF (Eq. 1). We do so by first bootstrapping (random re-sampling with replacement) the original observational data into subsets half the temporal length of the original data to produce an ensemble of records. For example, in the CERES data used to calculate the RS, we produce an ensemble of radiative climatologies derived from random selections of 9 years of the 18 years of data. This procedure queries how sensitive the radiative climatologies are to the limited length of the CERES record. Similarly, the surface albedo changes are calculated from the difference of random selections of 5 year averages within the period 1982-1991 and 2007-2016. We then use the re-sampled data to calculate the RS -- using the isotropic model-- and IS in a Monte-Carlo simulation. We calculate 100 different estimates of RS and 100 different estimates of IS with 50 derived from re-sampled NSIDC ice concentration data and 50 derived from re-sampled APP-x data. Thus, our estimates of IS (Fig. 8D) account for two sources of uncertainty: (i) the impact of intra-decadal variability on calculating longer term changes in surface albedo and (ii) instrumental uncertainty. The within data set intra-decadal variability of surface albedo contributes more to the IS uncertainty than the differences between APP-x and NSIDC sea ice concentration data sets; the standard deviation in IS calculated from ensembles of just the 50 NSIDC or 50 APP-x data is similar to that derived from the 100 member ensemble considered collectively. Given that the NSIDC estimate of surface albedo change is derived from sea ice concentration changes only and does not account for changes in the albedo over ice, the similarity of the NSIDC and APP-x derived IS suggest that albedo changes are primarily associated with changes in ice area, in opposition to the findings of \citep{horvat_albedo}. The uncertainty in RS (taken as 2 standard deviation across the re-sampled ensemble) is approximately 10$\%$ of the mean RS with larger values in the vicinity of sea ice edge (Fig. 8B) suggesting that the cloud properties that determine the RS are fairly constant from year-to-year. In contrast, the uncertainty in the IS (Fig. 8D) is approximately 60$\%$ of the mean value with particularly large uncertainties in the Beaufort Sea suggesting that the intra-decadal variability and measurement uncertainty of sea ice changes substantially hinders the calculation of long term IS over the relatively short observational record. We now describe how we use the uncertainty in IS and RS to calculate the uncertainty in RI$_{TOA, \alpha}$, the spatial average of which gives the SIAF uncertainty. We diagnose uncertainty RI$_{TOA, \alpha}$ by convoluting the 100 estimates of RS and the 100 estimates in IS to produce 10,000 estimates of RI$_{TOA, \alpha}$. This procedure accounts for the spatial co-variance of IS and RS uncertainty and central estimates. For example, the uncertainty in IS will have a larger impact in the regions and seasons where RS is largest. The uncertainty in the RI$_{TOA, \alpha}$ looks like and is comparable in fractional magnitude to that in surface albedo change with a slight modification by the spatial pattern of the mean RS. The spread in the spatial average of these 10,000 RI$_{TOA, \alpha}$ is combined with the uncertainty in global mean temperature changes -- propagated in quadrature since both quantities are scalars-- to produce a probability distribution function of SIAF (dark black distribution in left panels of Fig. 9). These calculation give an Arctic SIAF of 0.14 $\pm$ 0.4 W m$^{-2}$ K$^{-1}$ where the uncertainty is taken as 2$\sigma$ . The uncertainty in the observational global SIAF can be decomposed into contributions from the RS and IS uncertainty as follows: (i) the contribution of RS is calculated as 2$\sigma$ of the distribution derived from the 100 estimates of RS and multiplied by the OBE IS and (ii) the contribution of IS is calculated as 2$\sigma$ of the distribution derived from the 100 estimates of IS and multiplied by the mean RS. The uncertainty in the observational SIAF is almost entirely ($\pm$0.04 W m$^{-2}$ K$^{-1}$) due to uncertainty in the IS (dark blue narrow distribution in Fig. 9) whereas the uncertainty in the RS contributes very little to the global uncertainty in the SIAF ($\pm$0.003 W m$^{-2}$ K$^{-1}$ -- the very narrow dark red distribution in the left Fig. 9A). %---------------------------------------------------------------------------------------- % ACKNOWLEDGEMENTS %---------------------------------------------------------------------------------------- \bibliographystyle{ametsoc2014} \bibliography{dhoe_bib_full} %%%%%%%%%%%%%%%%%%%%%%%% TABLES %%% \begin{table}[t] \caption{SIAF values (in W m$^{-2}$ K$^{-1}$) for the (top) Arctic and (bottom) Southern Ocean derived from (left) Observations and model simulations of (middle) 4XCO$_{2}$ and (right) historical simulations. Each value shows the central estimate and 2$\sigma$ range across the bootstrapping Monte-Carlo simulations for the observations and inter-model spread for the models. The top row in each hemisphere shows the full calculation using the model specific RS and IS. The second row shows the impact of inter-model differences in IS as calculated using the model specific IS and the observed RS. The third row shows the impact of inter-model differences in RS as calculated using the model specific RS and the observed IS.}\label{t1} \begin{center} \begin{tabular}{c|c|c|c} \hline\hline \multicolumn{4}{c}{\textbf{Arctic}} \\ \hline\hline & Observations & 4XCO$_{2}$ & Historical\\ \hline \hline Full Calculation & 0.16 $\pm$ 0.04 & 0.13 $\pm$ 0.09 & 0.12 $\pm$ 0.13 \\ \hline \begin{tabular}{@{}c@{}}IS contribution \\ RS$_{OBS}$ $\times$ IS \end{tabular}& & 0.13 $\pm$ 0.08 & 0.12 $\pm$ 0.12 \\ \hline \begin{tabular}{@{}c@{}}RS contribution \\ RS $\times$ IS$_{OBS}$ \end{tabular}& & 0.16 $\pm$ 0.04 & 0.16 $\pm$ 0.04 \\ \hline\hline \multicolumn{4}{c}{\textbf{Southern Ocean}} \\ \hline\hline & Observations & 4XCO$_{2}$ & Historical\\ \hline \hline Full Calculation & & 0.08 $\pm$ 0.13 & \\ \hline \begin{tabular}{@{}c@{}}IS contribution \\ RS$_{OBS}$ $\times$ IS \end{tabular}& & 0.07 $\pm$ 0.11 & \\ \hline \begin{tabular}{@{}c@{}}RS contribution \\ RS $\times$ IS$_{OBS}$ \end{tabular}& & 0.08 $\pm$ 0.04 & \\ \hline\hline \end{tabular} \end{center} \end{table} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \begin{table}[t] \caption{\TO{Global albedo feedback (GAF)} in CMIP5 climate models calculated using the methodology of this study -- with a model specific RS from the isotropic model-- compared to that calculated using RS from the GFDL surface albedo kernel for all models. The CMIP5 ensemble mean and 2 $\sigma$ are shown for each hemisphere and divided into ocean and full domains.}\label{t2} \begin{center} \begin{tabular}{c|c|c} \hline\hline \multicolumn{3}{c}{\textbf{Northern Hemisphere}} \\ \hline\hline & Ocean Domain & Total \\ \hline \hline This Study & 0.15 $\pm$ .10 & 0.27 $\pm$ 0.18 \\ \hline GFDL RS kernel & 0.11 $\pm$ 0.06 & 0.20 $\pm$ 0.11 \\ \hline\hline \multicolumn{3}{c}{\textbf{Southern Hemisphere}} \\ \hline\hline & Ocean Domain & Total \\ \hline \hline This Study & 0.09 $\pm$ 0.16 & 0.10 $\pm$ 0.17 \\ \hline GFDL RS kernel & 0.08 $\pm$ 0.14 & 0.09 $\pm$ 0.15 \\ \hline\hline \multicolumn{3}{c}{\textbf{Global}} \\ \hline\hline & Ocean Domain & Total \\ \hline \hline This Study & 0.24 $\pm$ 0.15 & 0.37 $\pm$ 0.19 \\ \hline GFDL RS kernel & 0.19 $\pm$ 0.11 & 0.29 $\pm$ 0.13 \\ \hline\hline \hline \end{tabular} \end{center} \end{table} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% FIGURES %%%%%%%%%%%% \begin{sidewaysfigure}\label{arc_kernel} \centering \includegraphics[width=50pc]{F1.pdf} \caption{Arctic summertime (MJJA) surface albedo radiative sensitivity (RS) calculated from radiative kernels (upper panels) and estimated from the climatological radiative fields using the idealized isotropic radiation model (lower panels) in the same models. The squared spatial correlation coefficient between the kernel isotropic methods in the same model are provided in the middle and the Arctic domain averaged values are shown in the upper right of each panel.Note that there is no RS calculation from the GFDL model because mean state fields from this simulation were not saved.} Observational estimates from CERES EBAF satellite data and the isotropic model are shown to the right. \end{sidewaysfigure} \begin{sidewaysfigure}\label{so_kernel} \centering \includegraphics[width=50pc]{F2.pdf} \caption{As in Fig. 1 but for the NDJF RS in the Southern Ocean. Domain averaged surface albedo feedbacks exclude the Antarctic continent. Note that the figures are ordered by domain average RS over the Southern Ocean and this order differs from Fig. 1.} \end{sidewaysfigure} \begin{sidewaysfigure}\label{cartoon} \centering \includegraphics[width=50pc]{F3.pdf} \caption{Schematic of the single layer isotropic model modified from T07.} \end{sidewaysfigure} \begin{sidewaysfigure}\label{arc_kernel_cs} \centering \includegraphics[width=50pc]{F4.pdf} \caption{Comparison of Arctic Summertime (MJJA) full sky and clear sky surface albedo kernels.} \end{sidewaysfigure} \begin{sidewaysfigure}\label{arc_kernel_param} \centering \includegraphics[width=40pc]{F5.pdf} \caption{Comparison of the atmospheric opacity parameters that result from the application of the isotropic model to the Arctic summertime (MJJA) mean state radiative fields in the different climate models and observations. (Top Row) All-sky RS repeated from Fig. 1. (Second Row) Cloud reflectivity defined as the isotropic reflectivity applied to the all-sky radiative fields minus that defined from the clear-sky fields with the latter shown in the third row. (Bottom row) All-sky absorptivity. The (full Arctic) domain average is shown in the upper right of each panel. The four models for which kernels are available are shown to the left columns and the observational calculation from CERES data is shown to the right.} \end{sidewaysfigure} \begin{sidewaysfigure}\label{arctic_iso} \centering \includegraphics[width=50pc]{F6.pdf} \caption{Arctic summertime (MJJA) radiative sensitivity estimated using the isotropic model and the climatological radiation fields for CMIP5 historical simulations. Models are ordered as in reading a book (left to right then down) according to the domain average albedo feedback. Asterisks denote the models for which radiative kernel calculations are available that have been repeated from Fig. 1. The dark purple line shows the sea ice edge designated by the MJJA 50 $\%$ sea ice concentration contour. The full domain spatial average is shown in the upper left corner of each panel in black, the Arctic ocean average is shown in the lower right corner in blue and, the spatial average over the sea ice is shown in the lower left corner in purple. Observational estimates from CERES satellite data are shown in the bottom right panel.} \end{sidewaysfigure} \begin{sidewaysfigure}\label{so_iso} \centering \includegraphics[width=50pc]{F7.pdf} \caption{As in Fig. 6 but for the NDJF Southern Ocean RS. Domain averaged surface albedo feedbacks exclude the Antarctic continent. The dark purple line shows the sea ice edge designated by the NDJF 50 $\%$ sea ice concentration contour. Note that models are ordered by Southern Ocean domain averaged RS and this order differs from that in Fig. 6. Asterisks denote the models for which radiative kernel calculations are available that have been repeated from Fig. 2.} \end{sidewaysfigure} \begin{figure}\label{obs_maps} \centering \includegraphics[width=24pc]{F8.pdf} \caption{Spatial maps of observational estimates of summertime (MJJA) radiative sensitivity(RS, top), ice sensitivity (IS, middle) and the radiative impact of surface albedo change ( $RI_{TOA,\alpha}$, bottom). The RS is calculated from the isotropic shortwave model applied to the CERES data. The IS is calculated from Observational Best Estimate (OBE) surface albedo change between 1982 and 2016 divided by the global mean surface temperature change. The left panels show the central estimates of each quantity and the right panels show the uncertainty (2 standard deviations, $\sigma$) calculated from a Monte Carlo bootstrapping re-sampling with replacement as described in the Appendix.} \end{figure} \begin{figure}\label{feedback_spread} \centering \includegraphics[width=37pc]{F9.pdf} \caption{Estimates of global (and annual) SIAF from climate models and observations using the radiative sensitivity (RS) from the isotropic model applied to the climatology and the change in surface albedo under external forcing normalized by the global mean temperature change. (Upper left) Arctic sea ice changes over the historical (2007 to 2016 minus 1982 to 1991 averages). The black bars show the CMIP5 model distribution using the climate model specific radiative sensitivity and ice changes, the blue bars show the distribution using the model specific sea ice changes and observational RS and the red bars show the distribution using the observational sea ice change and model specific radiative sensitivity. Solid vertical lines show the model mean of each distribution. The dashed vertical line shows the observational estimate. The overlayed dark and thinner distribution shows the histogram of observational estimates of ice albedo feedback calculated from a Monte Carlo re-sampling of subsets of the ice albedo data and radiative data; the black distribution shows the impact of uncertainties in the observational RS and IS combined, the blue distribution shows the impact of the IS uncertainty only and the red shows the impact of the RS uncertainty only. (Lower Left) As in the above panel except using the modeled changes in the 4XCO$_{2}$ simulations. (Lower Right) Distribution of surface albedo feedback in the Southern Ocean diagnosed from 4XCO$_{2}$ normalized sea ice changes. Because the observational estimate of sea ice changes over the historical simulation is not statistically significant, the red distribution is calculated from the model specific radiative sensitivity and the model mean normalized sea ice change.} \end{figure} \begin{figure}\label{soden_iaf} \centering \includegraphics[width=37pc]{F10.pdf} \caption{Comparison of ice albedo feedback calculated from CMIP5 4XCO$_{2}$ using (ordinate) the method of \citet{soden_kernel} with RS in all models set to the GFDL surface albedo kernel versus (abscissa) the method introduced here with RS calculated from the model specific climatological radiative fluxes via the isotropic model. The blue markers show the contribution of the ocean domain only and the red markers show the full domain. All values shown are the contribution to the global mean. Dots show individual models and filled squares show the ensemble average with bars showing one standard deviation of the mean. The upper left shows the NH, the upper right the SH and the bottom the global mean. The black line is the 1:1 line.} \end{figure} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% appendix figure \begin{figure}\label{arc_kernel_donohoe} \centering \includegraphics[width=37pc]{AP1.pdf} \appendcaption{A1}{Arctic summertime (MJJA) surface albedo radiative sensitivity (RS) calculated from radiative kernels (upper panels) and estimated from the climatological radiative fields using the idealized isotropic radiation model of T07 (middle panels) and \citep{donohoe_battisti_11a} in the same models (bottom panels). The squared spatial correlation coefficient between the kernel isotropic methods in the same model are provided in the middle and the Arctic domain averaged values are shown in the upper right of each panel.} \end{figure} \begin{figure}\label{ker_iso_scatter} \centering \includegraphics[width=32pc]{AP2.pdf} \appendcaption{A2}{(Top panel) Scatter plot of MJJA radiative sensitivity calculated by (ordinate) radiative kernels and (abscissa) the isotropic model from the mean state in the same climate model. All four climate models and Arctic gridpoints considered collectively. (Bottom panel) Scatter plot of MJJA radiative sensitivity calculated from radiative kernels in one model versus the radiative sensitivity calculated from radiative kernels in a different model (selected at random). The dashed black line shows the 1:1 line.} \end{figure} \begin{figure}\label{ker_iso_scatter_SO} \centering \includegraphics[width=32pc]{AP3.pdf} \appendcaption{A3}{(Top panel) Scatter plot of NDJFM radiative sensitivity calculated by (ordinate) radiative kernels and (abscissa) the isotropic model from the mean state in the same climate model. All four climate models and Southern Ocean gridpoints considered collectively. (Bottom panel) Scatter plot of NDJF radiative sensitivity calculated from radiative kernels in one model versus the radiative sensitivity calculated from radiative kernels in a different model (selected at random). The dashed black line shows the 1:1 line.} \end{figure} \begin{figure}\label{comp_delta_alb} \centering \includegraphics[width=32pc]{AP4.pdf} \appendcaption{A4}{Comparison of the (MJJA) surface albedo changes (1982-2016) calculated from the NSIDC sea ice concentration data (left) and the APP-x surface albedo data (right).} \end{figure} \end{document}