The large spread of model equilibrium climate sensitivity (ECS) is mainly caused by the differences in the simulated marine boundary layer cloud (MBLC) radiative feedback. We examine the variations of MBLC fraction in response to the changes of sea surface temperature (SST) at seasonal and centennial time scales for 27 climate models that participated in the Coupled Model Intercomparison Project phase 3 and phase 5. We find that the intermodel spread in the seasonal variation of MBLC fraction with SST is strongly correlated with the intermodel spread in the centennial MBLC fraction change per degree of SST warming and that both are well correlated with ECS. Seven models that are consistent with the observed seasonal variation of MBLC fraction with SST at a rate −1.28 ± 0.56%/K all have ECS higher than the multimodel mean of 3.3 K yielding an ensemble‐mean ECS of 3.9 K and a standard deviation of 0.45 K.

1 Introduction Large intermodel disagreement in equilibrium climate sensitivity (ECS), the global‐mean surface temperature change in response to the doubling of CO 2 , has persisted for decades [Arrhenius, 1896; Callendar, 1938; Charney et al., 1979; Meehl et al., 2007; Stocker et al., 2013]. Although cloud feedback with a primary contribution from marine boundary layer clouds (MBLC) has been identified as a leading factor that causes the large spread in model ECS [Cess et al., 1990; Bony and Dufresne, 2005; Andrews et al., 2012; Zelinka et al., 2012], constraining and improving the simulations of cloud feedback in general circulation models (GCMs) remain a challenge. Recently, a number of studies demonstrated that there exist quantities in present‐day climate simulations linked to cloud feedback in future climate change. Constraining these quantities using observations thus provides practical means to constrain ECS [Fasullo and Trenberth, 2012; Klein et al., 2013; Sherwood et al., 2014; Su et al., 2014; Tian, 2015]. All these emergent constraints are based on the multiyear mean climate states simulated in the models, including mean relative humidity, cloud fraction, circulation strength, and precipitation bias. An unanswered question is why climatological mean states bear any relations to future cloud changes in response to surface warming. Intuitively, one would expect the physical processes that govern cloud changes on long time scale, such as centennial time scale of 100 years or more, to be equally active on shorter time scales, given that the response time of atmospheric dynamic and radiative processes related to clouds is on the order of hours to days. Thus, there have been studies that examined short‐term cloud feedback on the interannual and decadal time scales [Desseler, 2010; Clement et al., 2009] to gain insights into cloud feedback for long‐term climate change. Based on the similarities in tropic environmental changes between seasonal variation and global warming, Fasullo and Trenberth [2012] used relative humidity to constrain the model ECS. Tsushima et al. [2013] found that the models that best capture the seasonal variations of cloud regimes tend to have higher climate sensitivity. The relationship between the seasonal and long‐term variations of clouds with sea surface temperature (SST), however, has not been established despite the seasonal variation of clouds that has been well observed [Klein and Hartmann, 1993], and a strong correlation of snow albedo feedback on the seasonal time scale and in the climate change context was reported for models [Hall and Qu, 2006]. This study aims to fill the gap to establish a linkage between seasonal cycle and cloud changes at centennial time scale, with a particular focus on the MBLC fraction, one of the most important parameters that determine the shortwave cloud radiative feedback. Once the linkage is established, the observations of seasonal variation of MBLC are used to discriminate the low cloud shortwave feedback and thus models' ECS. Unlike previous studies that used climatological mean states to constrain ECS, this study tackles the seasonal variabilities in current climate and their relevance to long‐term climate change.

2 Data and Methods We quantify the MBLC fraction changes in response to SST variations at both seasonal and centennial time scales using outputs from 16 models that participated in the World Climate Research Programme's (WCRP's) Coupled Model Intercomparison Project phase 5 (CMIP5) and 11 models from the Coupled Model Intercomparison Project phase 3 (CMIP3) archive. These models (listed in Table 1) are chosen based on the availability of their ECSs from the literature [Gettelman et al., 2012; Andrews et al., 2012; Sherwood et al., 2014; Su et al., 2014]. The model data are downloaded from the Earth System Grid Federation (ESGF) (http://esgf.org and ftp://ftp‐esg.ucllnl.org). This study uses the following variables: “surface temperature,” “cloud area fraction at each vertical layer,” and “vertical component of velocity.” For CMIP5 models, we generate present‐day climate from the historical run data (except for CNRM‐cm5, whose seasonal response is derived using Atmospheric Model Intercomparison Project (AMIP) run because no 3‐D cloud fraction is available for its historical run) covering 1980–2004 and future climate from the Representative Concentration Pathway 4.5 (RCP4.5) run, the Representative Concentration Pathway (RCP) for the greenhouse gas radiative forcing to reach 4.5 W/m2 at the end of 21st century, covering 2074–2098 [Taylor et al., 2012]. For the CMIP3 models, we used the 20th Century experiment for current climate covering period 1970–1999 and 1%/year CO 2 increase experiment to doubling for climate change scenario (http://www‐pcmdi.llnl.gov/ipcc/standard_output.html#Experiments). Table 1. Horizontal Resolutions, Vertical Levels, ECS, and Number of Levels Below 700 hPa for the 16 CMIP5 and 11 CMIP3 Models Used in This Study CMIP Phase Modeling Center Acronym−Model Name Resolution (Longitude × Latitude × Levels) ECS (K) Number of Levels Below 700 hPa CMIP3 CCCMA‐cgcm3.1 3.75° × 3.6802°, L31 3.4 15 CMIP3 GISS‐e‐r 5° × 3°, L17 2.7 5 CMIP3 IAP‐fgoals1 2.8125° × 6.1315°, L26 2.3 6 CMIP3 INM‐cm3 5° × 4°, L21 2.1 7 CMIP3 IPSL‐cm4 3.75° × 2.5352°, L19 4.4 6 CMIP3 MIRCO‐medres 2.8125° × 2.7673°, L20 4.0 6 CMIP3 MPI‐cham5 1.875° × 1.849°, L31 3.4 9 CMIP3 MRI‐cgcm2 2.8125° × 2.7673°, L17 3.2 4 CMIP3 NCAR‐ccsm 1.4062° × 1.389°, L26 2.7 6 CMIP3 NCAR‐pcm1 2.8125 ° × 2.7673, L18 2.1 6 CMIP3 UKMO‐hadcm3 3.75° × 2.5°, L19 3.3 6 CMIP5 CCCMA−canesm2 2.8125° × 2.79°, L35 3.69 15 CMIP5 CNRM−cm5a 1.40625° × 1.4°, L31 3.25 9 CMIP5 CSIRO‐QCCCE−mk3.6 1.875° × 1.86°, L18 4.08 7 CMIP5 GFDL−cm3 2.5° × 2°, L48 3.97 12 CMIP5 GISS‐e2‐h 2.5° × 2°, L29 2.30 9 CMIP5 GISS−e2‐r 2.5° × 2°, L29 2.11 9 CMIP5 IAP‐fgoals‐g2 2.8125° × 6.1315°, L26 3.45 6 CMIP5 INM−cm4 2° × 1.5°, L21 2.08 7 CMIP5 IPSL−cm5a 3.75° × 1.8947°, L39 4.13 9 CMIP5 MIROC−miroc5 1.40625° × 1.4°, L40 2.72 6 CMIP5 MIROC−miroc‐esm 2.8125° × 2.79°, L80 4.67 14 CMIP5 MOHC−hadgem2‐es 1.875° × 1.25°, L38 4.59 12 CMIP5 MPI−esm‐lr 1.875° × 1.87°, L47 3.63 9 CMIP5 MRI−cgcm3 1.125° × 1.12°, L48 2.60 10 CMIP5 NCAR−cam5 1.25° × 0.9424°, L30 4.10 9 CMIP5 NCC−noresm1‐m 2.5° × 1.8947°, L26 2.80 6 The observational data include the cloud fraction from the merged CloudSat radar and Cloud‐Aerosol Lidar and Infrared Pathfinder Satellite Observations (CALIPSO) data (RL‐GEOPROF‐LIDAR, version R04) described by Mace et al. [2009] and the monthly mean surface temperature measured by Advanced Microwave Scanning Radiometer–Earth Observing System (AMSR‐E) (http://nsidc.org/data/amsrel1a.html) obtained from the NASA Observations for Model Intercomparison Projects (Obs4MIPs) (https://www.earthsystemcog.org/projects/obs4mips/) hosted on the ESGF. We also used the European Centre for Medium‐Range Weather Forecasts (ECMWF) ERA‐Interim global reanalysis vertical pressure velocity [Dee et al., 2011]. The monthly mean could fraction data are produced using the CloudSat/CALIPSO daily cloud fraction data from July 2006 to December 2010 regridded onto a 144 (longitude) × 91 (latitude) × 42 (altitude) grid. The uncertainties of CloudSat/CALIPSO cloud fraction data are estimated to be 5% with a bias of −5 to −7% [Mace et al., 2009; Su et al., 2013; Dolinar et al., 2014], and AMSR‐E SST data have a systematic bias of ~0.03 K and a standard deviation of ~0.4 K [Gentemann et al., 2010]. Additional biases due to undersampled diurnal cycle for A‐Train satellites in the monthly mean cloud fraction are estimated to be no more than ~2% over the tropical oceans [Jiang et al., 2012]. These data have been used for model‐observation comparisons [Su et al., 2013; Jiang et al., 2012; Li et al., 2012]. MBLC, topped by the low‐level trade inversion, occurs most frequently in large‐scale subsidence regions over subtropical oceans [Norris, 1998]. Following Bony and Dufresne [2005], we use the monthly mean midtropospheric (500 hPa) vertical pressure velocity (ω 500 ) to determine large‐scale subsidence regions. The MBLC fraction is computed as the total cloud area fraction below 700 hPa over oceans under a random‐overlap assumption [Manabe and Strickler, 1964; Ramanathan et al., 1983; Stephens, 1984]. The number of levels below 700 hPa is listed in Table 1. We have also performed the analysis using the maximum‐overlap assumption [Stephens, 1984] and found that the results are insensitive to the cloud‐overlap assumptions. of a physical quantity X over N years in the following way (1) To average out interannual variations, we compute a monthly climatologyof a physical quantityoveryears in the following waywhere (lon, lat) is the longitude and latitude at the center of a grid box and mon labels the month in a year. X over these regions is given by (2) This study examines the variations at both seasonal and centennial times scales of averaged MBLC fraction over the oceanic subsidence regions within 20°N−40°N and 20°S−40°S. The tropical band between 20°S and 20°N is excluded for its weak seasonal cycle. The spatial average ofover these regions is given by To separate MBLC from low‐level clouds over convective regions, we further limit the spatial averages to the regions of large‐scale subsidence, indicated by ω 500 that changes from month to month. Thus, in equation 2, only grid points with positive ω 500 are included in the spatial averages.

3 MBLC Fraction Variation at the Seasonal and Centennial Time Scales With the methods in the previous section, we computed the seasonal variations and centennial changes of MBLC fraction with changing SST averaged over oceanic subsidence regions within 20°N−40°N and 20°S−40°S for the 27 models. We also computed corresponding averages of observational MBLC fraction and SST using the CloudSat/CALIPSO cloud fraction, AMSR‐E SST, and the ERA‐Interim ω 500 from July 2006 through December 2010. Figure 1 displays the spatial average of MBLC cloud fractions for northern (Figure 1a) and southern (Figure 1b) subtropics, respectively, as functions of the underlying averaged SST in each calendar month for current (blue squares) and future (red circles) climate simulated by CMIP5 models. The first 14 subplots are for 14 CMIP5 models (CNRM‐cm5 and INM‐cm4 are not shown due to not having valid 3‐D cloud fraction data for RCP4.5), and the last subplot displays observational data for comparison. The CMIP3 results are similar but not shown for compactness. Figure 1 Open in figure viewer PowerPoint Seasonal variations of MBLC fraction averaged over oceanic subsidence regions within (a) 20°N–40°N and (b) 20°S–40°S as function of averaged underlying SST for both current and future climate. The first 14 subplots are for 14 GCMs from CMIP5, and the last subplot shows the seasonal variations of observed MBLC derived from 4 year CloudSat/CALIPSO MBLC cloud fraction as a function of averaged SST measured by AMSR‐E. Blue squares and red circles represent present and future climate, respectively. Models produce diverse simulations of MBLC: the overall range is 20%–65% (30%–75%) for northern (southern) hemisphere compared with the observation's 35%–50% (45%–55%). It is clear that MBLC is not uniquely determined by SST; other factors like seasonal variations of circulation, water vapor amount, lower tropospheric stability, or estimated inversion strength (EIS) also affect the amount of MBLC [Qu et al., 2014]. This is especially true for northern hemisphere, where a hysteresis behavior shows in the MBLC‐SST relations for all the models and the observations. For southern hemisphere, models CCCMA‐canesm2, GFDL‐cm3, IAP‐fgoal‐g2, MIROC‐esm, MIROC‐miroc5, MOHC‐hadgem2‐es, and NCAR‐cam5 and the observations all show that the average of MBLC over subsidence oceans in 20°S−40°S is mainly determined by the corresponding SST. For the rest of models, the MBLC fractions still have significant dependency on physical quantities other than SST. The different seasonal MBLC‐SST relations among the models and observations demonstrate that current GCMs do not simulate MBLC consistently. Despite the differences and complexities, in general, the averaged MBLC fraction bears an approximately linear relationship with SST on seasonal time scale in all models and the observations. We will focus on this overall correlation and will not examine the cause and effect of the interaction between MBLC and SST. Majority of the models agree with the observation, showing a decreasing MBLC fraction with increasing SST (negative slope) in both hemispheres. Interestingly, the seasonal cycle of MBLC in the future climate is an evident shift from the current seasonal cycle; i.e., the future seasonal variation with SST is at the same rate as the current climate except a shift in the mean temperature associated with global warming. For most of the models, the future climate shifts from the current climate along the overall trend of the seasonal variation, i.e., the overall seasonal variation is aligned with the shift of the MBLC under global warming. α season between and over the 12 months via (3) We define the seasonal MBLC fraction response to SST change as the regression slopebetweenandover the 12 months via from current to future climate, The centennial scale MBLC fraction response to surface warming for a calendar month (mon) of a year is computed as the ratio of the change offrom current to future climate,to the change of α century is the average of the MBLC fraction responses for all 12 months of a year, (4) The mean centennial scale MBLC fraction response to surface warmingis the average of the MBLC fraction responses for all 12 months of a year, The seasonal MBLC fraction response α season and centennial responses α century are computed for 16 CMIP5 and 11 CMIP3 models, separately for both hemispheres. Figure 2 displays the relationships between the seasonal and centennial MBLC fraction response to SST in northern (Figure 2a) and southern (Figure 2b) hemispheres. The intermodel variations of the response sensitivities α season and α century are well correlated with correlation of 0.6 (0.73) for northern (southern) hemisphere, above the 99.9% statistical significance threshold ~0.59 for 25 models (for future climate, the cloud fraction of INM‐cm4 is unphysically large and CNRM‐cm5 does not output 3‐D cloud fraction). Because the models are not completely independent as noted by Knutti et al. [2013], we adopt a bootstrap method [Efron, 1979] to compute empirical probability distributions of the correlations between α season and α century for both hemispheres using 40,000 samples generated by randomly drawing models from the 25 models. Figures 2c and 2d display the estimated empirical probability density functions of the correlation, which show that correlation values are mostly between 0.4 and 0.9 (0.6 and 0.9) for northern (southern) hemisphere. There is no negative correlation in the 40,000 samples giving a very high confidence of positive correlation. The double‐peak for northern hemisphere is due to model 3, which behaves like an outlier; had we excluded it, the correlation would have been 0.72. Two vertical dashed green lines mark the correlation values 0.57 (0.43) and 0.49 (0.31) at which the cumulative probability from the left is respectively 0.01 and 0.001 for the northern (southern) hemisphere. The identity function y = x is shown in Figures 2a and 2b, indicating that α season and α century are close to each other in magnitude, especially for southern hemisphere. Figure 2 Open in figure viewer PowerPoint Sensitivities of averaged MBLC response to SST at centennial time scale (α century ) versus seasonal scale (α season ) for the oceanic subsidence region in (a) 20°N−40°N and (b) 20°S−40°S, respectively, for the 25 models (CNRM‐cm5 and INM‐cm4 are not shown due to not having valid 3‐D cloud fraction data for future climate). The correlations and the y = x line are shown. (c and d) Empirical probability density functions of the correlations between α century and α season in northern and southern hemispheres estimated using a bootstrap method. Green vertical lines mark the correlations at which the cumulative probabilities from the left side are 0.001 and 0.01, respectively. Blue vertical lines mark the correlations for using all the 25 models as shown in Figures 2a and 2b. The linkage between the seasonal and centennial changes of MBLC fraction with SST is not surprising because the typical response time of MBLC to external forcing is much smaller than a month [Garratt, 1992]. Hence, seasonal sensitivity α season can serve as a proxy for the long‐term response α century , which is intimately related to the low cloud feedback and thus the model climate sensitivity [e.g., Bony and Dufresne, 2005; Clement et al., 2009; Zelinka et al., 2012; Fasullo and Trenberth, 2012; Sherwood et al., 2014; Su et al., 2014]. It is also interesting that the sensitivities α season and α century for northern hemisphere are similar to those for southern hemisphere; the values for both hemispheres are well correlated across different models with confidence level above 99.9% (see supporting information for more details).

4 Correlation With ECS and Observational Constraints Previous studies [Bony and Dufresne, 2005; Andrews et al., 2012; Zelinka et al., 2012] showed that the boundary layer could feedback explained most of the spread in ECS. Therefore, the MBLC fraction response to SST, critical to shortwave cloud feedback, is expected to be correlated with model ECS. Figure 3a confirms this: the intermodel variation of α century (average of both hemispheres) is well correlated with ECS with a correlation of −0.69. The greater reduction of MBLC fraction leads to a higher ECS. Similarly, α season is also well correlated with model ECS. The intermodel spread in α season explains about 41% variance of intermodel spread in ECS with a correlation coefficient of −0.64 (Figure 3b). Statistically, the probability for α century (α season ) and ECS being negatively correlated is greater than 99.9% assuming that models are independent. Figure 3c (Figure 3d) displays the probability density function of the correlation between α century (α season ) and ECS estimated from the bootstrap method with randomly drawing of 40,000 samples showing that correlation −0.68 (−0.64) is near the center of the distribution with 99% confidence to have correlation < −0.47 (−0.37) and with 99.9% confidence to have correlation < −0.37 (−0.27) as marked by the green lines. For this correlation, MBLC response at seasonal scale α season provides a metric of shortwave cloud radiative feedback and opens a new way to constrain model climate sensitivity using temporal variabilities in current climate, different from previous studies that focused on current climate means [Fasullo and Trenberth, 2012; Sherwood et al., 2014; Su et al., 2014; Tian, 2015]. Figure 3 Open in figure viewer PowerPoint Relations between model ECSs and MBLC response sensitivities to SST change at (a) centennial scale α century and (b) seasonal scale α season . (c and d) Empirical probability density functions for the α century ‐ECS correlation and α season –ECS correlation, respectively, as estimated using the bootstrap method. The gray area in Figure 3b marks the observed MBLC response sensitivity using CloudSat/CALIPSO cloud fraction data and AMSR‐E SST data with the uncertainties (±3σ). Green vertical lines mark the correlations at which the cumulative probabilities from the right side are 0.001 and 0.01, respectively. Blue vertical lines mark the correlations for using all the available models as shown in Figures 3a and 3b. We compute the observed α season using CloudSat/CALIPSO cloud fraction and AMSR‐E SST. The uncertainty of α season is estimated as the year‐to‐year standard deviation of α season for the 4 years from 2006 to 2010, which include two EL Niños and two La Niñas [Su and Jiang, 2013]. The larger value out of the northern and southern hemispheres is used for the uncertainty estimate. Shown in Figure 3b, the gray area represents ±3σ (standard error) (0.56%/K) around the observed α season of −1.28 K. Model seasonal responses α season has a range from −2.1 to 0.73%/K with 18 models capturing the decreasing MBLC with SST. Seven models, namely, CCCMA‐cgcm3.1, UKMO‐hadcm3, CCCMA‐canesm2, CSIRO‐mk3.6, IPSL‐cm5a, NCAR‐cam5, and MOHC‐hadgem2‐es, have seasonal MBLC response to SST consistent with the observations. Assuming that these models produce a realistic low cloud radiative feedback, the “best estimate” of ECS has a mean of 3.9 K with a standard deviation of 0.45 K. This best estimate of ECS is larger than the multimodel‐mean of 3.3 K and has a much smaller standard deviation than 0.81 K from the 27 models, consistent with previous studies [Fasullo and Trenberth, 2012; Sherwood et al., 2014; Tsushima et al., 2013; Su et al., 2014; Tian, 2015]. Similar results hold if we use the AMIP experiments to estimate seasonal response α season . While the result of a higher ECS is consistent with previous studies, this study is unique in that (1) it uses the temporal variabilities of current climate instead of the mean state; (2) it uses the seasonal MBLC fraction response to SST, which is directly related to the low cloud feedback, to constraint the ECS, and hence provides a clearer physical meaning than other proxies. Because accurate observational record of MBLC fraction is rather short, the robust seasonal variation of MBLC fraction with SST becomes an attractive metric for model evaluation and constraint of ECS.

5 Conclusions Using model simulations from 27 climate models, we estimate the MBLC fraction response to surface warming over oceanic subsidence regions within 20°N–40°N and 20°S–40°S at both seasonal and centennial (climate change) time scales. We find that MBLC fraction response to SST warming at seasonal and centennial scales is strongly correlated as in the case of snow albedo feedback [Hall and Qu, 2006]. Therefore, for the MBLC, seasonal response can be a proxy for the response under global warming in these models. Because the model climate sensitivity differences are primarily caused by the differences in simulating the shortwave radiative feedback from the MBLC, the intermodel spread in the seasonal MBLC response is strongly correlated with that in ECS. Using CloudSat/CALIPSO cloud fraction and AMSR‐E SST, we find the models that simulate a seasonal MBLC response to SST consistent with the observation yield an ensemble‐mean ECS of 3.9 K, about 0.6 K higher than the multimodel‐mean ECS from all the 27 models analyzed. Our results support the previous studies that favor the high ECS and underscore the importance of MBLC variations to climate feedback.

Acknowledgments We thank Mark Zelinka and Joel Norris for their useful discussions and the two anonymous reviewers for their constructive suggestions to improve the manuscript. This work is conducted at JPL/Caltech under a contract with NASA, supported by the ROSES NDOA, MAP, and NEWS programs. The model data are downloaded from the ESGF (http://esgf.org and ftp://ftp‐esg.ucllnl.org). We acknowledge the modeling groups, the Program for Climate Model Diagnosis and Intercomparison and the WGCM, for their roles in making available the WCRP CMIP3/CMIP5 multimodel data set. Support of this data set is provided by the Office of Science, U.S. DOE. The merged CloudSat radar and CALIPSO lidar data (RL‐GEOPROF‐LIDAR, version R04) are downloaded from http://www.cloudsat.cira.colostate.edu/data‐products/, and the monthly mean SST measured by AMSR‐E is obtained from the Obs4MIPs project (https://www.earthsystemcog.org/projects/obs4mips/) hosted on the ESFG. The ERA‐Interim global reanalysis vertical velocity data are downloaded from the ECMWF data server via http://www.ecmwf.int/en/research/climate‐reanalysis/era‐interim.

Supporting Information Filename Description grl53577-sup-0001-supinfo.pdfPDF document, 67.6 KB Figure S1 Please note: The publisher is not responsible for the content or functionality of any supporting information supplied by the authors. Any queries (other than missing content) should be directed to the corresponding author for the article.