Shortened Duration of Global Warming Slowdowns with Elevated Greenhouse Gas Emissions

Continuous emissions of anthropogenic greenhouse gases (GHGs) and aerosols in the last 160 years have resulted in an increasing trend of global mean surface temperatures (GMSTs). Due to interactions with natural variability, rates of the combined anthropogenically and naturally induced warming trends are characterized by significant slowdowns and speedups on decadal timescales. Here, by analyzing observed and model-simulated data, we investigate how the duration of these episodes will change with different strengths of GHG and aerosol forcing. We found that the duration of warming slowdowns can be more than 30 yr with a slower rate of anthropogenic emissions but would shorten to about 5 yr with a higher one. This duration reduction depends on both the magnitude of the climate response to anthropogenic forcing and the strength of the internal variability. Moreover, the warming slowdowns can still occur even towards the end of this century under high emissions scenarios but with significantly shortened duration.


Introduction
Although concentrations of atmospheric greenhouse gases (GHGs) and aerosols have been increasing continuously, the rising of the global mean surface temperature (GMST) is not monotonic but shows periodically significant slowdowns and speedups, e.g., the warming slowdown (sometimes referred to as the warming hiatus) from 1998 to 2013 (Easterling and Wehner, 2009;Schmidt et al., 2014;Fyfe et al., 2016) and the speedup since 2014 (Rahmstorf et al., 2017) in association with phase transitions of the Interdecadal Pacific Oscillation [IPO, or Pacific Decadal Variability (PDV); Hu and Fedorov, 2017;Su et al., 2017]. The studies suggest that these decadalto-multidecadal timescale GMST variations mainly result from combined effects of external forcing and internally generated climate variability (Meehl et al., 2009(Meehl et al., , 2011(Meehl et al., , 2013(Meehl et al., , 2016aHansen et al., 2011;Solomon et al., 2011;Held, 2013;Kosaka and Xie, 2013;England et al., 2014;Santer et al., 2014;Watanabe et al., 2014;Dai et al., 2015;Nieves et al., 2015;Steinman et al., 2015;Smith et al., 2016;Chen and Tung, 2018;Haustein et al., 2019;Wu et al., 2019). The natural internal variability includes enhanced ocean heat uptake due to changes of the internal climate processes associated in part with Atlantic Multidecadal Variability [AMV, often associated with the Atlantic Multidecadal Oscillation (AMO)] and PDV (often associated with the IPO; Meehl et al., 2009Meehl et al., , 2011Meehl et al., , 2013Meehl et al., , 2016aDai et al., 2015;Nieves et al., 2015;Steinman et al., 2015;Chen and Tung, 2018;Haustein et al., 2019;Wu et al., 2019). PDV is characterized, for example, in its negative phase by surface cooling over the eastern Pacific (Kosaka and Xie, 2013) associated with intensified trade winds Watanabe et al., 2014) and a slowdown in the rate of GMST increase. External forcing includes changes in natural forcings (solar and volcano; Hansen et al., 2011;Solomon et al., 2011;Santer et al., 2014) and anthropogenic GHGs and aerosols (Smith et al., 2016;Wu et al., 2019). Alternative explanations suggested that either external forcing alone may be able to explain most multidecadal variations in GMST and AMV is primarily controlled by external forcing (Haustein et al., 2019), or potential biases in observed data could also have contributed to this variability (Karl et al., 2015). Nevertheless, the rate of warming in the early 21st century is much less than that of the late 20th century based on multiple observed datasets (Chen and Tung, 2018). Accompanying with global warming slowdown in the early 21st century, opposite seasonal changes are observed in recent decades, such as more frequent severe winter cold-air outbreaks and heavy snowfalls over northern Eurasia and North America in winter (Cohen et al., 2012(Cohen et al., , 2018, but more frequent and intensified summer heatwaves across most regions in the northern midlatitudes (Schär et al., 2004;Barriopedro et al., 2011). Here, we will not further explore the underlying physical mechanisms driving these periodic warming slowdowns/speedups under a rising influence of anthropogenic forcing. Instead, we focus on how the occurrence and duration of the warming slowdowns will change as climate change intensifies in the 21st century.
According to the IPCC Fifth Assessment Report (IP-CC AR5), GMST will continue to rise with an amplitude ranging from 0.3-1.7°C for low emissions scenario to 2.6-4.8°C for high emissions scenario by the end of this century (relative to the base period of 1986-2005; Collins et al., 2013). Although the multimodel ensemble mean time evolution of GMST from models participating the Coupled Model Intercomparison Project phase 5 (CMIP5) shows a monotonic rise, there is still variability on decadal-to-multidecadal timescales if individual simulations are used (Fyfe et al., 2013;Marotzke and Forster, 2015). The chances for the occurrence of significant speedups/slowdowns of warming in the future would depend on the strength of the anthropogenic forcing in the future projected climate scenarios (Maher et al., 2014) and the phases of natural internal climate variabilities. A previous study suggests a nearly zero chance for a decade-long warming slowdown to occur in the high-emissions scenario by the end of the century (Maher et al., 2014); however, this will not exclude the potential occurrence of a warming slowdown shorter than a decade. Because these speedup/slowdown warming periods could generate significant impacts on regional and global climates (Cohen et al., 2014), it is important for us to understand how the characteristics of these phenomena would change in the future climate. Therefore, the interplay between natural internal climate variability and external forcings (e.g., GHGs and aerosols) may not only be important to the current and past climates, but also important to the future projected climate.
Here, we will analyze the surface temperature changes from observations and model simulations and explore the changes in the duration of these variabilities featuring with special focus on warming slowdowns. The model simulations include historical and projected periods during 1850-2100, and the idealized 1% CO 2 runs from models participating in both CMIP5 (Taylor et al., 2012) and CMIP6 (Eyring et al., 2016). Although there are several major future climate scenarios, only the highest GHG emissions scenarios will be used in the main text: Representative Concentration Pathway 8.5 (RCP8.5; Riahi et al., 2011) from CMIP5 and Shared Socioeconomic Pathway 5 (SSP5-8.5;O'Neill et al., 2016) from CMIP6. SSP5-8.5 is an updated scenario from RCP8.5 with an elevated CO 2 concentration of a few hundred ppmv by the end of the 21st century (O'Neill et al., 2016) due to different assumptions used to derive these two scenarios, especially the industrial CO 2 emissions (such as the CO 2 emissions is about 7 Gt yr −1 higher in SSP5-8.5 than RCP8.5 after about 2070; Meinshausen et al., 2019).

Gridded observational and model datasets and annual mean estimates
We primarily used the Hadley Centre-Climate Research Unit combined land surface air temperature and sea surface temperature (HadCRUT4; Morice et al., 2012) version 4.6.0.0. HadCRUT4 is a 5°lon × 5°lat gridded dataset based on the observed global historical monthly surface temperature anomalies (relative to a 1961-1990 reference period) from 1850 to 2018, neither interpolated nor variance adjusted). As pointed out previously (Wu et al., 2019), the spatial and temporal coverage of the surface temperature observations is not uni-form, especially before the 1920s, and in regions poleward of 60°. As before Wu et al. (2019), the annual mean data is derived from the grid points where data are available at least eight months per year. Regions poleward of 60° are excluded from our calculation of the GMST due to the significant lack of observed samples. Three additional datasets of global monthly surface temperature from 1881 to 2018 are also used. One is the global 2° × 2° grid NASA Goddard Institute for Space Studies (GISS) Surface Temperature Analysis dataset (GISTEMP; Hansen et al., 2010), which covers most regions of the globe by sampling at the station level for land, and using ship-based and satellite-based measurements for ocean. The global 1° × 1° grid Berkeley Earth Surface Temperatures (BEST) are merged land-ocean dataset that combines the land analysis with an interpolated version of HadSST3 from the USA's NOAA and NASA, and the UK's CRU/Hadley Centre (Rohde et al., 2013). The third dataset is NOAA's Merged Land Ocean Surface Temperature Analysis (NOAAGlobalTemp), which combines land surface air temperatures primarily from the Global Historical Climatology Network-Monthly (GHCN-M) version 3 with SSTs of the ERSSTv3b analysis into a comprehensive global surface temperature dataset spanning 1880 to the present at monthly resolution, on a 5° × 5° grid (Vose et al., 2012).
Simulations of 16 models from CMIP5 and 18 models from CMIP6 were used (Table 1). Idealized one percent CO 2 increase per year compound (1pct_CO 2 hereafter) and historical experiments during 1850-2005 and RCP8.5 experiments during 2006-2100 from the CMIP5 models were used. As for the CMIP6 models, their simulations include 1pct_CO 2 and historical experiments during 1850-2014 and SSP5-8.5 scenario experiments during 2015-2100.

Equivalent CO 2 (CO 2ea ) data derived
Radiative forcing is an important feature of GHGs and aerosols to influence on the balance of incoming and outgoing energy in the earth-atmosphere system. In order to consider the integrated effects of anthropogenic forcing due to all well-mixed GHGs and anthropogenic aerosols, the CO 2ea concept is employed. The total effective radiative forcing (ERF) induced by GHGs (e.g., CO 2 , CH 4 , NOx, and CFCs) and aerosols from 1881 to 2011 based on CMIP5 in Table AII.1.2 (Collins et al., 2013) are converted to CO 2ea data by using Eq. (1).
where CO 2orig is 278 ppmv in the year 1750 and the ERF is relative to a zero baseline in 1750 and represents additional forcing compared to pre-industrial levels. The CO 2ea data during 1881-2011 are directly derived from Eq. (1). It is noted that there exists a highly linear correlation relationship between CO 2 and CO 2ea during the period of 1991-2011 ( Fig. 1).
This relationship is used to estimate CO 2ea data by using the given atmospheric CO 2 concentrations from 2012 to 2100 in this work.

Internal variability of GMST on different timescales
The internal variability of GMST on different timescales is estimated by using the time series of 1881-2018 annual GMST residuals from which the CO 2ea -induced GMST are substracted. The two-standard deviation value of the given N-yr smoothed residual surface temperature is regarded as the internal variability from its lowest to highest values or from its highest to lowest values on Nyr timescales. If it is greater than surface temperature  warming due to the increased anthropogenic forcing, the internal variability induced by, for example, AMV/PDV could overcome the GHG-induced warming and could contribute to a global warming slowdown.

GMST related to GHGs and aerosols variations in the 21st century
The global warming mainly comes from increase of GHGs including CO 2 and other well-mixed gases (e.g., CH 4 , NO x , and CFCs). However, changes in anthropogenic aerosols can also modulate GMST over the 20th century, especially in the period of 1940-1975 when GMST was nearly unchanged (Easterling and Wehner, 2009). In order to estimate the collective effects of the anthropogenic forcing (GHGs and aerosols) on GMST, we adopt a concept of CO 2ea . Here, CO 2ea includes both GHGs and anthropogenic aerosols, slightly different from a previous study (Wu et al., 2019). To derive CO 2ea , the accumulated ERF from both GHGs (including CO 2 and other well-mixing GHGs) and anthropogenic aerosols is first converted to that of CO 2 based on the AR5 in Table AII.1.2 (Collins et al., 2013) from 1881 to 2010. Then the sum of the total radiative forcing is converted to CO 2ea concentration time series. The CO 2ea data from 2011 to 2100 are derived from the newly available CO 2 observations (2011-2014) and scenarios of RCP8.5 and SSP5-8.5 (2015-2100) by using the approximately linear relationship between CO 2ea and CO 2 during 1991-2011 from AR5-published data (see Fig. 1). As shown in Fig. 2, CO 2ea increases from 282 ppmv in 1880 to 421 ppmv in 2018 with the corresponding observed CO 2 concentration increasing from 291 ppmv in 1880 to 407 ppmv in 2018. However, the projected CO 2ea increases from 431 ppmv in 2021 to 1052 ppmv in 2100 for RCP8.5 (CMIP5), and from 432 to 1292 ppmv for SSP5-8.5 (CMIP6) due to different assumptions used to derive these scenarios (Meinshausen et al., 2019). In terms of the rate changes, the rising of CO 2ea changes from about 0.33 ppmv yr −1 before the 1950s to about  2.31 ppmv yr −1 after the 1970s, and is projected to be 5−7 ppmv yr −1 in the mid-21st century, and over 10 ppmv yr −1 for RCP8.5 and about 16 ppmv yr −1 for SSP5-8.5 by the end of the 21st century (Table 2). How will these changes in anthropogenic forcings (GHGs and aerosols) and natural internal decadal variability affect the GMST? To explore these, it becomes important to separate the contributions from anthropogenic forcing due to GHGs and aerosols and those from natural internal climate variability to the GMST or regional surface air temperature. Different methodologies have been used in previous studies (Dai et al., 2015;Meehl et al., 2016a;Wu et al., 2019). A commonly used method (Dai et al., 2015;Meehl et al., 2016a) is based on climate model simulations in which the externally-forced GMST changes from observations are represented by using the multimodel ensemble mean of 20th century model simulations, this method is assuming that internal climate variability would be smoothed out by this multimodel ensemble due to the random phases of the natural internal climate variability. However, the externallyforced GMST changes deduced this way include not only effects of GHGs and anthropogenic aerosols but also effects of natural forcings (solar and volcanoes). In our previous work (Wu et al., 2019), we found a more physically-based quasi-linear relationship between transient GMST and well-mixed GHG changes for both observations and model simulations. In this study, by considering the collective radiative effects of GHGs and anthropogenic aerosols, a similar quasi-linear relationship is found to exist between CO 2ea (GHGs plus aerosols) and GMST in observations and model simulations (Fig. 3). Figure 3b shows the scatter plot of GMST versus the corresponding concentrations of the atmospheric CO 2ea in each year from 1880 to 2018 based on HadCRUT data.
Other observed datasets are also tested (Fig. 4) and the results are similar to HadCRUT (Fig. 3). The fitted linear line in Fig. 3b is calculated by using the scatter points of 9-yr running mean of Surface Air Temperature (SAT) versus CO 2ea . It is highly correlated (0.91) to the yearly data of 9-yr running mean of GMST (Fig. 3b) and is statistically significant at the 99% confidence level based on the Student's t-test. This significant linear relationship between GMST and CO 2ea also exists in the 20th century historical simulations of CMIP5 and CMIP6 as shown in Figs. 3c, d from 1850 to 2014. To extend the CMIP5 historical experiment to 2014, the experiment using RCP4.5 is used. The use of RCP4.5 does not affect this extension since the difference of the GHG concentration among different RCPs is very small during 2006-2014. To be consistent with observations, GMST for the model simulations is also defined as an area-weighted mean between 60°S and 60°N. Additional analysis by including the surface temperature of both poles from model simulations does not alter our results. In Figs. 3c, d, the red points show the CMIP5 and CMIP6 multimodel ensemble means of the yearly SAT from 1881 to 2014 versus the corresponding yearly CO 2ea , respectively. Although the transient GMST changes on decadal or longer timescales for both observations and model simulations can be scaled linearly by changes of CO 2ea concentration, the resulting rate of GMST change per 100 ppmv CO 2ea is closer to observed for CMIP6 model simulations than CMIP5 although they both fall into the uncertainty range of the observations, with CMIP5 on the high end (Fig. 3) along with NOAA dataset (Fig. 4). This may imply that the simulated 20th century climate could be better in CMIP6 than in CMIP5 and the underlying physical processes leading to this better simulation is needed to be studied further, but beyond the scope of this study.
As indicated by our previous study (Wu et al., 2019), this linear relationship between GMST and CO 2ea is valid almost perfectly for a CO 2ea change of around/less than 150 ppmv, but a logarithmic relationship would work better with a high rate of CO 2ea change. As shown in Fig.  5, a logarithmic fit between CO 2ea and GMST exists in both multimodel ensemble means (CMIP5 and CMIP6 models) for simulations including the idealized experiments with 1pct_CO 2 , and the historical experiments plus RCP8.5/SSP5-8.5 projections (1850-2100). The GMST variations in the 1pct_CO 2 experiments show a much larger spread and higher sensitivity to GHGs in CMIP6 models than CMIP5 models (Fig. 5), which is associated with the larger spread of the equilibrium climate sensitivity in CMIP6 models due to a significant increase of the equilibrium climate sensitivity at the high end from some CMIP6 models relative to the corresponding CMIP5 models (Stevens et al., 2017;Meehl et al., 2020;Zelinka et al., 2020). For example, there are 14 out of the 40 CMIP6 models whose equilibrium climate sensitivity is over 4.5°C, but only 2 in all the CMIP5 models. A study indicates that this higher equilibrium climate sensitivity in CMIP6 models is mainly from changes of low cloud feedback in extra tropical regions (Zelinka et al., 2020). Interestingly, the multimodel ensemble mean change of GMST in the historical simulations is 0.73°C per 100ppmv CO 2ea increase for CMIP6 models, which is nearly identical to that derived from observations (0.71°C per 100-ppmv CO 2ea increase), but 0.93°C per 100-ppmv CO 2ea increase for CMIP5 models which is closer to that derived from GISS dataset. This controversy between historical simulations and idealized 1pct_CO 2 simulations may reflect the changes of the anthropogenic aerosol optical depth used in CMIP6 (Stevens et al., 2017) and may also be related to the much stronger cooling effect in association with the stronger indirect effects of aerosols in some CMIP6 models. These indirect effects could be tuned to make the simulated GMST close to ob-  Fig. 3. Concentration of the atmospheric CO 2ea , and observed and simulated global mean surface temperature (GMST) anomalies. GMST is defined as the area-weighted mean surface temperature between 60°S and 60°N. (a) Annual-mean (thin black line) and its 9-yr smoothed (thick black line) time series of GMST anomalies relative to the 1881-1910 mean derived from HadCRU data (methods), and atmospheric CO 2 (red line) and CO 2ea (blue line) concentrations. (b) Red and blue points denote the annual mean from 1880 to 2018 and its 9-yr smoothed GMST change with CO 2ea , respectively. (c, d) CMIP5 and CMIP6 historical simulations from 1881 to 2014, and the reference climate is the mean from 1881 to 1910. Solid black lines in (b-d) denote the linear trend of the 9-yr smoothed GMST; red dots in (c, d) represent the ensemble mean of yearly GMST, and the shading represents the GMST spread of the 16 CMIP5 and 12 CMIP6 models. MME represents the multimodel ensemble mean.
servations (Golaz et al., 2013;Mauritsen and Roeckner, 2020). The logarithmic relationship between CO 2ea and GMST in observations can be derived by using the Had-CRUT GMST during 1951-2018 (Fig. 6, also see Fig. 7 for other observed datasets) since this is a period with rapid increase of the global GHG and aerosol emissions (Fig. 2). In comparison with the estimated linear and logarithmic fits between GMST and CO 2ea in observations, it shows a similar correlation between the linear fit and Annual mean 9-yr running average Fig. 4. As in Fig. 3, but for concentration of the atmospheric CO 2ea , and observed GMST anomalies from BEST, GISTEMP, and NOAAGlobal-Temp datasets. All the GMST are defined as the area-weighted mean surface temperature between 60°S and 60°N.

APRIL 2021
actual GMST (0.914) to that of logarithmic fit (0.910) for this time period (Fig. 6a). However, the extrapolation of this linear fit to 2100 will result in a significant overestimation of the GMST change in comparison to CMIP6 models (6.1 vs. 4.8°C), but an underestimation from the logarithmic fit (3.4 vs. 4.8°C) if SSP5-8.5 scenario is realized. If lower scenarios are realized, the logarithmic fit derived from observations could do a much better job to project/predict future GMST. Moreover, the logarithmic relationship between CO 2ea and GMST shown in Figs. 5, 6 suggests that with an increasing rate of CO 2ea change, the GMST change per unit CO 2ea will gradually reduce. As demonstrated in the previous study (Washington et al., 2000), a slower rate of CO 2ea changes gives more time for the climate system to equilibrate with the external forcing (CO 2ea ) and a faster rate of CO 2ea changes leads to a significantly lagged response since there is less time for the climate system to respond to forcing changes. The underlying physical mechanisms are the ocean heat uptake and changes of the ocean circulation in response to the changes in CO 2ea forcing, which is associated with the vast heat capacity of the ocean in comparison to that of the land and where the heat absorbed by the ocean is deposited. Table 2 lists the mean rate of CO 2ea changes and the GHG-induced GMST changes in six periods from 1880 to 2100. The results clearly demonstrate that the rate of GMST increase gradually slows down with an enhanced CO 2ea incre- ment. For example, with an increase of CO 2ea of 0.33 ppmv yr −1 before the 1960s, the GMST change is about 0.76 ± 0.02°C per 100-ppmv CO 2ea increase. This reduces to 0.59 ± 0.01°C in recent decades with an increase rate of CO 2ea changes to about 2.3 ppmv yr −1 . This GMST change rate reduces further to around 0.40 ± 0.01°C in the mid-21st century with a CO 2ea increase of about 7 ppmv yr −1 , and 0.20 ± 0.004°C by the end of the 21st century with a CO 2ea increase of about 16 ppmv yr −1 .
Although GMST change rate declines with an enhanced increment of CO 2ea for a unit CO 2ea change, GMST changes for a given fixed time period increase as total CO 2ea rises. As shown in Table 2, the increasing rate of GMST was about 0.025 ± 0.001°C decade −1 before the 1960s, rising to 0.137 ± 0.003°C decade −1 in recent decades, and is expected to further rise to 0.280-0.359°C decade −1 in the mid-to-late 21st century (2061)(2062)(2063)(2064)(2065)(2066)(2067)(2068)(2069)(2070)(2071)(2072)(2073)(2074)(2075)(2076)(2077)(2078)(2079)(2080). It is interesting that this decadal GMST increase rate reduces by the end of the 21st century to 0.252-0.321°C decade −1 , which is related to the flattened rate of CO 2ea increment in this time period as shown in Fig. 1 (dotted  lines). This suggests that the changes of GMST depend not only the amount of CO 2ea , but also the rate of CO 2ea amount changes. Overall, this result indicates that with a slower rate of CO 2ea change, the increase of GMST per unit CO 2ea change is larger than that with a faster rate of  Fig. 6. Observed GMST and CO 2ea , and CO 2ea -induced GMST changes and internal variability. (a) Relationships between observed GMST and CO 2ea . Black dots denote the changes of annual mean HadCRUT GMST anomalies (60°S-60°N, relative to 1881-1910 mean) with CO 2ea from 1881 to 2018, the red line shows the GMST changes fitted with logarithmic CO 2ea , the shading area denotes the fitted uncertainty, and the green dash line shows the GMST changes linearly with CO 2ea . (b) Increase of GMST with different rates of CO 2ea increment. The black line denotes two-standard deviations of the internal-variability induced GMST variations on timescales from annual to 43 years, and 5 shaded color belts represent anthropogenic forcing induced GMST variations with 5 different mean CO 2ea increment rates for the 5 periods (1880-1950, 1981-2018, 2021-2061, 2061-2080, and 2081-2100). The first three periods are based on observations and the last three periods are from the SSP5-8.5 scenario.
APRIL 2021 CO 2ea changes since a lower rate CO 2ea change means that it will take a longer time to make a unit CO 2ea change. However, the total GMST change for a given time period is larger for a faster rate of the CO 2ea change than for a slower one.

Duration of global warming slowdown periods in the 21st century
To assess the influence of the internal variability on GMST, the CO 2ea -induced GMST is removed from the raw annual mean observed data (HadCRUTv4 data from 1880 to 2018) based on the method outlined in the previous section. A slowdown in the warming rate occurs only when the GMST downward trend induced by natural internal variability is larger than the anthropogenic forcinginduced GMST upward trend. Here, we use the twostandard deviation value of the residual GMST after the removal of the anthropogenic forcing-induced GMST to represent the intensity of natural internal variability-induced GMST variations, which represents a 95% confidence level based on the Student's t-test. To derive the changes of natural internal variability on different timescales, the standard deviation of the residual GMST is calculated based on different running means of the residual GMST time series. Such as, to derive the natural internal variability on the 5-yr timescale, we first apply a 5-yr running mean to the annual mean residual GMST time series and then calculate the standard deviation of this new time series to represent the natural internal variability-induced GMST variation on a 5-yr timescale. Here, running-mean operation is based on the consideration to more sample numbers of time series on a given timescale and to avoid the uncertainty due to periods artificially chosen. As shown in Fig. 6, we applied this operation 22 times (with a smooth window ranging from 1 to 43 yr) to get the natural internal variability-induced GMST changes on timescales of 1, 3, 5, …, and 43 yr. Here, the timescales up to 43 yr are to taken into account the duration of three global warming slowdown periods (1896-1910, 1941-1975, and 1998-2013) occurred since 1891 in observational record (Meehl et al., 2016b). The internally-generated decadal timescale GMST variability mainly comes from natural internal sea surface temperature oscillations such as PDV (Meehl et al., 2009) with a typical period of 20-30 yr and AMV (Trenberth and Shea, 2006) with an estimated period of 60-80 yr. As shown in Fig. 6b, the amplitude of the internal variability-induced GMST changes decreases as the timescale of these internal variability lengthens, implying that the GMST changes induced by lower frequency variability are smaller than those induced by high frequency variability.
To get the duration changes of the warming slowdown periods, we calculate the GMST changes induced by GHGs and aerosols over periods of 1-43 yr. When the GMST changes induced by CO 2ea are larger than twostandard deviations of the natural internal variability-induced GMST for a given time period (varying from 1 to 43 yr, see methods), and this is defined as the possible maximum duration of the slowdown period. Since the rate of CO 2ea change is not constant during the historical and future projected periods (Fig. 2), CO 2ea increase rates  Fig. 6, but for the observed GMST and CO 2ea , CO 2ea -induced GMST changes, and internal variability from BEST, GISTEMP, and NOAAGlobalTemp datasets.
for five periods of 1880-1950, 1981-2018, 2021-2061, 2061-2080, and 2081-2100 in Fig. 6b (Table 2) are chosen to show how these different CO 2ea rates would affect the occurrence and duration of slowdown periods. The changes of GMST with CO 2ea increase for the given different ranges of CO 2ea variations in the five periods (Table 2) can be estimated from the time series of GMST versa CO 2ea (in Fig. 6a) that are derived from the logarithmic relationship between CO 2ea and GMST. Then, the CO 2ea -induced GMST variation with time for the given different CO 2ea increments per year in the five periods (Table 2) can be derived. Five color belts in Fig. 6b denote the GMST increase rate with time for those different conditions of CO 2ea increase per year, and the width of color belts represents the uncertainty caused by GMST changes per unit CO 2ea . These five rates of GMST increase are correspondent to 0.33 ppmv yr −1 during 1880-1950, about 2.31 ppmv yr −1 during 1980-2018, about 6.92 ppmv yr −1 during 2021-2060 from SSP5-8.5, about 13.40 ppmv yr −1 during 2061-2080 from SSP5-8.5, and about 16.28 ppmv yr −1 during 2081-2100 from SSP5-8.5, respectively ( Fig. 2 and Table 2).
As shown in Fig. 6b, with a rate of CO 2ea change at 0.33 ppmv yr −1 , it will take up to 31 yr in order for the GMST change induced by GHGs and aerosols to be larger than the two-standard deviation of the natural internal variability-induced GMST change. This suggests that if a slowdown occurs under this CO 2ea increase rate, potentially the slowdown period could last as long as 31 yr. This assessment is consistent with the observed slowdown periods that occurred during the 1940s-1970s when the observed GMST was nearly flat (no trend) or with a slight downward trend. With a rate of CO 2ea increase of about 2.31 ppmv yr −1 , the duration of a possible slowdown reduces to about 12 yr, which is nearly identical to the recent warming slowdown (1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013). This duration reduces further to 5-6 yr in the mid-21st century and by the end of the 21st century based on SSP5-8.5 and RCP8.5 scenarios. Our results indicate that higher rates of GHGs and aerosols increase will reduce the duration of future slowdown periods. Moreover, slowdown periods could still occur even towards the end of this century with the highest emissions scenarios, but would be much shorter than recent slowdown periods. In other words, periodic slower rates of change in GMST are inevitable given the influence of internal variability on the GMST changes for both the past and future climates. The only difference is that the duration of these events would not be a decade or longer, but much shorter than a decade with future higher rates of change of GHGs. Moreover, these shorter slowdown warmings are likely to occur less frequently (Maher et al., 2014).
Here, our analyses mainly focus on the occurrence of slowdown warmings for the high emissions scenarios of RCP8.5 and SSP5-8.5. There are other relative low emissions scenarios of RCP2.6, RCP4.5, RCP6.0, and SSP126, SSP245, SSP370 in the 21st century. In those conditions of future scenarios, the anthropogenically induced global warming rates are estimated lower than those for RCP8.5 and SSP5-8.5 scenarios, and the duration of slowdown warmings for other relative low emissions scenarios would be longer than 5-6 yr in the mid and end of the 21st century.

Summary
In summary, we evaluated the potential changes of warming slowdown durations for the climate of the past century and a half as well as the climate from now to the end of the 21st century. We found that the duration of warming slowdowns would decrease as anthropogenic forcing increases. This reduction depends on the magnitude of the CO 2ea changes for different emissions scenarios and the strength of the internal variability. Under the high emissions scenarios, warming slowdowns could still occur even towards the end of the 21st century, though with shorter durations. The internally generated climate variability such as PDV and AMV can cause the occurrence of warming slowdowns, but it is noticeable that they can also result in the warming speedup when the warming slowdowns end. It is worth pointing out that the global mean aerosol effect on ERF was considered when we built the CO 2ea time series. However, the regional heterogeneous distribution of the anthropogenic aerosols is neglected due to the difficulties in deriving the aerosol optical depth from the observed and projected time periods. Although the heterogeneous anthropogenic aerosols will lead to distinct duration changes on regional scales, the effect of this on the estimations for the GMST is expected to be small (Wu et al., 2019). To better evaluate the regional impact of both anthropogenic and volcanic aerosols on the warming slowdowns on regional scales, further studies are needed.