Global mean thermosteric sea level projections by 2100 in CMIP6 climate models

Most of the excess energy stored in the climate system is taken up by the oceans leading to thermal expansion and sea level rise. Future sea level projections allow decision-makers to assess coastal risk, develop climate resilient communities and plan vital infrastructure in low-elevation coastal zones. Confidence in these projections depends on the ability of climate models to simulate the various components of future sea level rise. In this study we estimate the contribution from thermal expansion to sea level rise using the simulations of global mean thermosteric sea level (GMTSL) from 15 available models in the Coupled Model Intercomparison Project Phase 6 (CMIP6). We calculate a GMTSL rise of 18.8 cm [12.8–23.6 cm, 90% range] and 26.8 cm [18.6–34.6 cm, 90% range] for the period 2081–2100, relative to 1995–2014 for SSP245 and SSP585 scenarios respectively. In a comparison with a 20 model ensemble from Coupled Model Intercomparison Project Phase 5 (CMIP5), the CMIP6 ensemble mean of future GMTSL (2014–2100) is higher for both scenarios and shows a larger variance. By contrast, for the period 1901–1990, GMTSL from CMIP6 has half the variance of that from CMIP5. Over the period 1940–2005, the rate of CMIP6 ensemble mean of GMTSL rise is 0.2 ± 0.1 mm yr−1, which is less than half of the observed rate (0.5 ± 0.02 mm yr−1). At a multi-decadal timescale, there is an offset of ∼10 cm per century between observed/modelled thermosteric sea level over the historical period and modelled thermosteric sea level over this century for the same rate of change of global temperature. We further discuss the difference in GMTSL sensitivity to the changes in global surface temperature over the historical and future periods.


Introduction
About 93% of the excess energy stored in the climate system due to anthropogenic greenhouse gas emissions has been absorbed by the oceans, leading to thermal expansion and sea level rise (Cheng et al 2017, Oppenheimer et al 2019. For more than 600 million people living in low-elevation coastal areas future sea level rise is one of the main damaging aspects of climate change (Oppenheimer et al 2019). Crucial decisions about adaptation to sea level rise for populous coastal megacities and communities in low lying small islands, where a considerable fraction of global economic activity and critical infrastructure exists, will be made based on future sea level projections (Vousdoukas et al 2018, Jevrejeva et al 2019, Oppenheimer et al 2019, Abadie et al 2020. In a warming climate, global sea level is expected to rise and future sea level projections are made using the conventional method by simulating contributions from individual sea level components, such as thermal expansion, melting ice from glaciers and ice sheets, changes in land water storage and summing them up (Church et al 2013, Oppenheimer et al 2019. The robustness of and confidence in projections of future sea level change depends on the ability of climate models to reproduce the components of sea level rise over the 20th century and simulate future changes across a range of emission scenarios (Church et al 2013, Oppenheimer et al 2019. New results for all sea level components under the most recent climate scenarios will be simulated by, or derived from results in the new Coupled Model Intercomparison Project Phase 6 (CMIP6) (Eyring et al 2016), including future contribution from the ice sheets in Ice Sheet Model Intercomparison Project for CMIP6 (Nowicki et al 2016).
In this study we analyse global mean thermosteric sea level (GMTSL), one of the components of sea level rise, resulting from thermal expansion due to heat uptake of the ocean. Developing a better understanding of heat uptake and heat transport by the ocean are crucial to simulate the climate system response to the radiative forcing in climate models (Church et al 2013, Oppenheimer et al 2019, Zanna et al 2019, Gregory et al 2020. Furthermore, increasing our understanding of large uncertainty in the components of future sea level rise, including the GMTSL, is crucial for effective communication of sea level projections to policy makers, coastal engineers and the general public (Jevrejeva et al 2019, Oppenheimer et al 2019. Our main objective is to estimate the thermosteric contribution of future global sea level rise from a new generation of climate models and new scenarios, SSP245 and SSP585 (Eyring et al 2016, O'Neill et al 2016. We evaluate the CMIP6 simulations of GMTSL against available observations (Ishii and Kimoto 2009, Levitus et al 2012, Cheng et al 2019 to assess the performance of the CMIP6 models. In addition, we compare GMTSL projections to the Coupled Model Intercomparison Project Phase 5 (CMIP5) outputs previously used for sea level projections (Taylor et al 2012, Yin 2012, Church et al 2013, Kopp et al 2014, Jackson and Jevrejeva 2016, Slangen et al 2017. Climate models considered in this study are atmosphere-ocean coupled general circulation models (AOGCMs) and Earth system models (ESMs) from the CMIP6 archive (Eyring et al 2016).
We expect CMIP6 GMTSL projections to be different to those of CMIP5 due to a new generation of climate models (e.g. Eyring et al 2016), as well as a new set of scenarios of concentrations, emissions, and land use (O'Neill et al 2016). We discuss the differences between GMTSL projections, suggesting several explanations for large uncertainties in simulation of GMTSL. Our study demonstrates the need to further improve climate model performance given the importance of oceans in the storage and redistribution of heat defining future climate change (Church et al 2013, Winton et al 2013, Melet and Meyssignac 2015

Global mean thermosteric simulations in CMIP6
We analyse GMTSL simulations from 15 climate models that are directly available as 'zostoga' in the CMIP6 database (https://esgf-node.llnl.gov/ search/cmip6/; Eyring et al 2016). Zostoga in CMIP6 models represents a part of the global mean sea level change due to changes in ocean density arising from changes in temperature (Griffies et al 2016. Table SI1 (available online at stacks.iop.org/ERL/16/014028/mmedia) displays the list of all 15 CMIP6 models currently available and used in this study. We use 'historical' simulations for the period 1850-2014 (Eyring et al 2016) and scenarios SSP245 and SSP585 simulations for the period 2015-2100. The new scenarios are based on a matrix that uses the shared socioeconomic pathways (SSPs) and forcing levels of the Representative Concentration Pathways (RCPs) discussed in O'Neill et al (2016). SSP245 and SSP585 correspond to RCP4.5 and RCP8.5 in the CMIP5 framework respectively. Only the first ensemble member 'r1i1p1f1' of the 15 CMIP6 models have been used for both historical and future scenarios. In cases where models provide 'f2' instead of 'f1' , that particular forcing was used. We apply a drift correction calculated from the full length pre-industrial control run corresponding to each climate model, for GMTSL simulations following the approach by Sen Gupta et al (2013) (see SI).

GMTSL simulations in CMIP5
The zostoga 'historical' simulations and simulations with RCP4.5 and RCP8.5 future scenarios from 20 CMIP5 models are used in this study (see table SI2 for model information). The 'historical' simulation spans 1850/1860-2005, while the RCPs span 2006-2100. As in the case of CMIP6, only the first member ensemble 'r1i1p1' was used for both historical and RCP simulations. Models were drift corrected following the same method as that for CMIP6 (see SI).

GMTSL from observations
We compare historical CMIP6 and CMIP5 simulations with updated observational GMTSL datasets from Levitus et al (2012), Ishii and Kimoto (2009) and Cheng et al (2019)

GMTSL simulations
Using 15 CMIP6 model simulations we compute a multi-model ensemble mean (MEM) for GMTSL rise of 18.8 cm with the 5-95th percentile of 12.8-23.6 cm for SSP245, and a MEM of 26.8 cm with the 5-95th percentile of 18.6-34.6 cm for SSP585 for the period 2081-2100 relative to 1995-2014 (figure 1(a), table SI1). Figure 1 We analyse the difference between the CMIP6 and CMIP5 MEMs of GMTSL in 'historical' simulations and for the future projections with SSP245/RCP4.5 and SSP585/RCP8.5 scenarios respectively (figure SI1). While there is a broad agreement between the CMIP6 and CMIP5 MEMs over the projection period, a clear discrepancy occurs between the two MEMs during the historical period. The considerable disagreement between the model simulations is also reflected in the magnitude of 5-95th percentile ranges in the CMIP6 and CMIP5 models. By 2100 CMIP6 model simulations show larger spread compared to the CMIP5, though the two overlap considerably (figure 1). Over the 1901-1990, the CMIP6 the 5-95th percentile range is nearly half that of CMIP5 (figure SI1).
Although a portion of their uncertainties overlap the CMIP5 ensemble contains a greater number of models (20) than CMIP6 (15). We test the hypothesis that the trend estimated for the MEM from CMIP5 is statistically distinct from that of CMIP6 by randomly drawing a subset from the 20 models in CMIP5 and estimating each MEM. We perform a bootstrap analysis by randomly sampling all combinations of 15 out of 20 CMIP5 models. This is then used to estimate the ensemble mean and compare with the ensemble mean of CMIP6 models. All possible combinations (>15 000) of 15 CMIP5 models were used. Figure 2 shows trend histograms of the MEM GMTSL from these 15 randomly selected CMIP5 models over three time periods: 1901-1990, 1957-2005 and 2015-2100 (RCP8.5/SSP585 projections used here). The ensemble mean CMIP6 (red) and observational GMTSL (black) trends are also displayed as vertical lines. Figures 2(a) and (b) shows that over the historical time period the rate of the CMIP6 MEM is distinctly smaller than the CMIP5 MEM, irrespective of which 15 models are selected from CMIP5. We point out here that the range of CMIP5 MEMs of GMTSL most of the time overestimates the observational trend (figure 2(b); figure SI2, in which the uncertainty for the trend in observations are included). This result contrasts with the projection period (figure 2(c)) where the rate of MEM GMTSL from CMIP6 is larger than all combinations of CMIP5 MEMs. The CMIP6 MEM shows a higher rate of GMTSL rise than the CMIP5 MEM over the future projections for both the scenarios (table 1). In the case of the historical simulations (period 1901-1990), the CMIP6 MEM rate is lower than that of the CMIP5 ensemble. The uncertainties of the rate of the CMIP6 and CMIP5 MEM, as shown in the table 1, are evaluated using Monte Carlo methodology. Around 10 000 random time series are produced by sampling a normal distribution centred on the ensemble mean time series. The uncertainty of the trend is then defined as two times the standard deviation of the normal distribution of the trends of the sampled time series (corresponding to the 95% confidence interval).
While these results reveal that both MEMs are distinct (figure 2), their variances overlap considerably (figure SI1). Visually (figure SI1) it appears that the variance of CMIP6 is smaller than of CMIP5 in the historical period, but larger than CMIP5 for both future scenarios. We evaluate the statistical significance of the CMIP6 variance relative to the variance in CMIP5 by again considering the difference in the number of the available models for CMIP6 and CMIP5 (figure SI3). For future GMTSL projections the CMIP6 variance is larger compared to the variance in CMIP5 and there is no intersection with the range of possible variances estimated from CMIP5 model subsets. We calculate the relative uncertainty (the variance divided by the median versus time) of CMIP6 and CMIP5 GMTSL for historical and projected periods (figure 3). For future projections, the relative uncertainties in both CMIP6 and CMIP5 models are almost constant over time, and independent of scenario for each CMIP respectively. Unlike the projection period, CMIP6 and CMIP5 show relative uncertainties that are time varying and non-linear over the historical period.

Comparison with observational GMTSL
We compare GMTSL ensemble mean from CMIP6 models with the ensemble mean of three observational datasets (Ishii, 0-2000  We compare the rate of GMTSL rise from CMIP6 MEM of 0.2 ± 0.1 mm yr −1 with observational steric sea level rise rate of 0.5 ± 0.02 mm yr −1 (Cheng et al 2017) over the period 1940-2005 (table 2), further indicating that CMIP6 models are underestimating observed changes in global steric sea level. There is a caveat in using dataset from Cheng et al (2017), as it is global steric sea level estimated from temperature and salinity related density changes; however, in CMIP6 only GMTSL is available and used in our study. Nonetheless, for the global scale the differences between global mean steric and GMTSL are almost negligible (Gregory and Lowe 2000, WCRP    Rate from CMIP6 MEM over the same period is 0.2 ± 0.1 mm yr −1 , which is also twice smaller compared to Storto et al (2019). These comparisons with reconstructed and reanalysis global mean thermosteric data over longer time show that the CMIP6 MEM is underestimating the observed GMTSL rate.

GMTSL change as a function of global surface temperature
The delayed response time of the ocean to the forcing means that thermosteric sea level rise in 2100 does not immediately reflect the forcing occurring at that time. Instead the entire pathway since the forcing change was introduced is important (e.g. Bouttes et al 2013).
Here we consider the relationship between the average rate of GMTSL rise and global surface temperature increase representing the entire preceding century. We estimate the rate of each model's GMTSL and global surface temperature for CMIP6 and CMIP5 models and their MEMs over the period  and 2015-2100 (figure 5) and introduce a GMTSL sensitivity as a function of the changes in the GMTSL (as rate in cm per century) to changes in global surface temperature ( • C per century). To place these results in a wider context we calculate sensitivities of GMTSL as a function of global surface temperature for the 20th century, early deglacial period, and a long-term climate stabilisation case. For the 20th century, we use global temperature reconstructions Had-CRUT4 (Morice et al 2012) and GISTEMP (Lenssen et al 2019, GISTEMP Team 2020 and GMTSL observational datasets (Ishii andKimoto 2009, Cheng et al 2017). For paleo conditions in the early part of the last deglaciation we have used a global temperature reconstruction (Marcott et al 2011) and global steric sea level, derived from North Atlantic hosing experiments under this climatology (Fluckiger et al 2006). We also estimate GMTSL sensitivity using 1000 year idealised climate model simulations that followed the SRES A1B scenario (close to RCP6.0 scenario) stabilising from 2100 to 3000 (Meehl et al 2007). In the latter two cases, while centennial rates of global surface temperature change are comparable to today, the alternate physical mechanism (Fluckiger et al 2006) and long-term ocean inertia illustrate contrasting GMTSL responses yet both differ from presentday rates or CMIP projections.
In figure 5 GMTSL sensitivities calculated for CMIP6 (crosses) are generally higher than those for CMIP5 (circles) for all future scenarios. As expected, there is a quasi-linear increase in GMTSL sensitivity to the changes in global surface temperature, which themselves are reflection of changes in both temperature and thermosteric sea level to radiative forcing (e.g. Bouttes et al 2013). However, we observe an offset of ∼10 cm per century between historical (since 1901, modelled: purple, observed: black symbols) and projected thermosteric sensitivities for the same rate of change of temperature (e.g. for RCP2.6). Despite this offset, CMIP6 and CMIP5 historical sensitivities show significant overlap though CMIP5 MEM GMTSL sensitivity is higher. When compared to the GMTSL sensitivities simulated by CMIP models with estimates from other sources, the observationally driven thermosteric sea level sensitivities are larger by a few cm per century. The initial response of climate models (Meehl A in figure 5) to SRES A1B aligns well with those in CMIP6 and CMIP5. As the forcing is held constant after 2100 (Meehl et al 2007), the climate model's temperature equilibrates quickly resulting in a very small rate of change while the inertia of the ocean heat uptake results in a thermosteric sensitivity with the same magnitude as that of the CMIP5 and CMIP6 future scenarios (Meehl B in figure 5).

Discussion
Under the SSP245 and SSP585 scenarios, the CMIP6 models display larger spread in future thermosteric sea level rise by 2100 compared to CMIP5. As discussed in O'Neill et al (2016), we expect CMIP6 climate projections to be different from those for CMIP5 due to a new generation of climate models (Eyring et al 2016), as well as a new set of scenarios of concentrations, emissions, and land use. One possible explanation for the larger spread in CMIP6 could be that the range of the projected global surface warming is wider in the SSP245 and SSP585 compared to the similar scenario used in CMIP5 (Zelinka et al 2020), as changes in GMTSL are related to global surface warming (figure 5). Global surface temperature projections from the selected 15 CMIP6 models, used in this study, show higher temperature changes by 2100 compared to the CMIP5 simulations (figure SI5), supporting the results for 27 CMIP6 in Zelinka et al (2020). In addition, for the models used in this study, the CMIP6 mean equilibrium climate sensitivity (ECS) is 3.8 ± 1.2 K, while it is 3.3 ± 0.6 K in CMIP5 (tables SI1 and SI2, based on Meehl et al 2020). The transient climate response (TCR) for the models used in the study shows a mean of 2.1 ± 0.5 K in CMIP6 while it is 1.8 ± 0.3 K for CMIP5 (tables SI1 and SI2, based on Meehl et al 2020). Thus, the MEM of CMIP6 TCR and ECS are both larger and with a larger spread than for CMIP5.
The response of GMTSL to radiative forcing differs from that of global surface temperature (e.g. Bouttes et al 2013). GMTSL depends upon the emissions pathway, in fact Bouttes et al (2013) and Melet and Meyssignac (2015) showed that while forcing is positive, GMTSL rise is linearly dependent upon the  (1901-1990, magenta) and future (2015-2100) scenarios RCP2.6 (yellow), SSP245/RCP4.5 (blue) and SSP585/RCP8.5 (red). Observational sensitivities (black triangles) using temperature/thermosteric sea level reconstructions (HadCRUT4/ Cheng et al 2017 andGISTEMP/Ishii andKimoto 2009). Initial (2000-2100, Meehl A, light grey) and long run (2100-3000, Meehl B, dark grey) sensitivity of climate model simulations under SRES A1B from Meehl et al (2007). Sensitivities estimated during early deglacial conditions (light blue) from global temperature reconstruction (Marcott et al 2011) and North Atlantic hosing experiments (Fluckiger et al 2006). integral of the forcing. Global surface temperature on the other hand depends roughly linearly upon the forcing (Gregory and Forster 2008). These relations led us to consider the effect of global surface temperature change upon GMTSL over multi-decadal time scales, and whether such analysis helps us to understand both the underestimation of historic GMTSL by both CMIPs and the greater variance of CMIP6 compared to CMIP5.
While this study has primarily considered the sensitivity of GMTSL to global surface warming, ocean heat uptake (and thereby thermosteric sea level) is also influenced by regional changes in ocean dynamics such as the North Atlantic Deep Water formation (Knutti and Stocker 2000) and weakening of AMOC circulation (Cheng et al 2019), ocean mixing (Watanabe et al 2020) and subduction processes (Griffies and Greatbatch 2012) though such an investigation into their influence is beyond the scope of this paper. While heat redistribution and transport in the ocean play an important part in its warming (Winton et al 2013, Zanna et al 2019, Cheng et al 2019, other reasons for greater sensitivity and ensemble spread include model resolution, numerical schemes and parameterizations (Kirtman et al 2012, Eyring et al 2016. Using the directly calculated zostoga parameter from CMIP6 may introduce a bias for the historical period in some models due to inclusion/exclusion of marginal seas in calculation of the GMTSL, as it has been suggested in the studies using selected CMIP5 models (e.g. CMIP5 models GFDL-CM3, MIROC-ESM, and GISS-E2-R in Slangen et al 2017; not used in this study). The difference up to 0.1 mm yr −1 in GMTSL trends for some individual model simulations could be explained by omitting volcanic forcing in pre-industrial runs, but including it in historical simulations. This might lead to differences during historical simulations and projections, as discussed in Melet and Meyssignac (2015) and Slangen et al (2017).
Estimated MEM rates of GMTSL in CMIP6 are smaller than those for CMIP5 from 1901 to 1990 and the difference in rates is statistically significant, which contrasts with the CMIP6 rate for the projections (2015-2100) being higher than in CMIP5. If we check individual models, it is hard to explain why some models do not show any GMTSL rise over the 20th century (e.g. NorESM2, figure SI6 and table SI1) while their future projections are very close to the ensemble mean. If we consider the link between the model simulations of GMTSL to the ECS (not shown here), it is only for the high emission scenario (SSP585) that we can demonstrate an obvious connection between the highest estimate for projected GMTSL and a high ECS in the models. However, for the historical period ) there is no noticeable relationship between the GMTSL simulations by CMIP6 models and ECS. This could be explained by the differences in the CMIP6 model response to the radiative forcing for the historic period (e.g.  and forcing for the future projections by 2100, supported by the differences in relative uncertainties for historic and for future periods (figure 3).
Forcing in future scenarios (e.g. SSP245 and SSP585) is dominated by the increase in CO 2 and other greenhouse gas concentrations (O'Neill et al 2016). For the 'historical' simulations, in addition to the anthropogenic CO 2 , natural forcing occurs as well including the variability of solar irradiance, 'volcanic' aerosol injected into the stratosphere by explosive volcanic eruptions, and anthropogenic aerosols (O'Neill et al 2016;Meehl et al 2007).
The climate response to volcanic forcing differs between atmosphere and ocean. For a few years following an eruption, volcanic aerosols cause a net negative radiative forcing of the climate system by reflecting sunlight (shortwave radiation) though this is partly offset by absorption of outgoing longwave radiation by the volcanic aerosol (Forster and Taylor 2006). In the ocean, the response to a single (Pinatubo-like) eruption comprises two primary time scales: one fast (1-2 years) and one slow (decadal). Over the fast time scale, the ocean sequesters cooling anomalies induced by the eruption into its depth, enhancing the damping rate of sea surface temperature relative to that which would be expected if the atmosphere acted alone (Gupta and Marshall 2018). However, because of its large effective heat capacity, the ocean memory of successive eruptions enables an accumulation larger than any single event (e.g. Hansen et al 2002, Grinsted et al 2007, Stenchikov et al 2009. This complex response by the ocean to the forcing during the historical period could provide some explanation for relatively lower sensitivity of GMTSL rise to the changes in global surface temperature ( figure 5).
Further studies that focus on developing an understanding of large uncertainties in heat uptake by the ocean in climate models (AOGCMs and ESMs) within CMIP6 would greatly contribute to the improvement of future sea level projections given the important role played by oceans in the storage and redistribution of heat and to modulate future changes in climate (Melet and Meyssignac 2015, Gregory et al 2020, Todd et al 2020.

Conclusion
We estimate a GMTSL rise of 18.8 cm (12.8-23.6 cm, 5-95th percentile) and up to 26.8 cm (18.6-34.6 cm, 5-95th percentile) by the end of the 21st century with SSP245 and SSP585 scenarios, respectively. Our study provides evidence that despite huge improvements in current AOGCMs and ESMs (O'Neill et al 2016, Eyring et al 2016 CMIP6 models show larger disagreement in projecting the GMTSL by 2100 compared to CMIP5 simulations, which will by extension create larger uncertainties in total, global and regional sea level projections. Over the historical time (since 1900) the comparisons with observed and reconstructed GMTSL show that overall CMIP6 models and CMIP6 MEM are underestimating the observed GMTSL.
Given the importance of GMTSL to sea level projections, an improved understanding of the causes of model spread and the role of natural variability and external forcing therein needs to be a high priority in the climate and sea level rise modelling community. As more outputs from CMIP6 models become available, more detailed analysis of GMTSL will be possible. On the other hand, evaluating model performance (by comparison with instrumental records) will remain limited by the lack of long-term historical observations of temperature and salinity in the ocean.

Data availability statement
The data that support the findings of this study are available upon reasonable request from the authors.