1 Introduction

Although a large inter-model uncertainty remains, discussions on future carbon emission pathways are assisted by the observation that the transient climate response to cumulative carbon emissions (TCRE) is largely unaffected by time and pathways (e.g., Collins et al. 2013). For example, in searching for pathways to maintain a given temperature target, one can first estimate the total allowable carbon emissions and then consider how to distribute the emissions for each year.

According to Rogelj et al. (2019), the remaining carbon budget (Blim) can be calculated from the TCRE using Eq. (1):

$$ {B}_{\mathrm{lim}}=\left({T}_{\lim -}{T}_{\mathrm{hist}-}{T}_{\mathrm{nonCO}2}-{T}_{\mathrm{ZEC}}\right)/\mathrm{TCRE}-{E}_{\mathrm{esfb}} $$
(1)

where Tlim, Thist , TnonCO2, and TZEC are the temperature change limit, the historical human-induced warming to-date, the non-CO2 contributions to a future temperature increase, and the zero-emissions commitment (ZEC), respectively, and Eesfb is an adjustment term for sources of unrepresented Earth system feedback.

In the Special Report on Global Warming of 1.5 °C by the Intergovernmental Panel on Climate Change (IPCC) (SR1.5; Rogelj et al. 2018), experts calculated the remaining carbon budget from Eq. (1) using TCRE values that were estimated in the IPCC AR5 (i.e., 0.8–2.5 K/GtC; Collins et al. 2013). Although this estimate was the result of expert judgment based on various existing studies using stabilizing, overshooting, and increasing scenarios, the TCRE range of 0.8–2.4 K/GtC estimated by Gillet et al. (2013) using Earth systems models (ESMs) seems remarkably influential. In addition, Gillet et al. (2013) used a scenario with a monotonic increase of 1% per annum (ppa) in the atmospheric CO2 concentration (pCO2). This corresponding to the threshold exceedance budgets (TEBs), or total carbon emissions of the underlying scenarios up to the time of exceeding the threshold, as defined by Rogelj et al. (2016). In other words, only transient warming is considered. The remaining carbon budget is more closely related to another type of budget, which was termed the threshold avoidance budget (TAB) by Rogelj et al. (2016) for the non-overshooting case, and which reflects the total carbon emissions of the underlying scenarios up to the time of peak warming (i.e., slow warming is also considered). For the average value, however, the results of the Zero Emissions Commitment Model Inter-comparison Project (ZEC-MIP; MacDougall et al. 2020) indicated that warming after the cessation of CO2 emissions (i.e., the ZEC) can be approximately estimated as zero (TZEC = 0). In that case, the TAB can be considered to be equal to the TEB.

Tachiiri et al. (2015) demonstrated that the uncertainty of the TCRE increases (i.e., large TCREs get larger and small TCREs get smaller) when the pCO2 is stabilized, although the average (or median) TCRE remains unchanged. Hence, the TABs can be smaller than the TEB if we need to achieve the target with a high probability, whereas the TABs and TEBs are expected to be similar when considering approximately the 50th percentile.

Although the initial main focus of researchers was the quasi-constancy of the TCRE, the scenario dependence of the TCRE has been previously reported (Hajima et al. 2012; Krasting et al. 2014). More recently, MacDougall (2017) concluded that the TCRE is only nearly constant (i.e., < 5% change) when changes in the pCO2 are 0.3–1.2% p.a. at 400 ppm (or 0.5–2.5% p.a. at > 1000 ppm). For overshooting scenarios, Zickfeld et al. (2016) indicated some increase in the TCRE and Tachiiri et al. (2019) quantified a 10–20% increase in the TCRE under a decline-after-a-peak pCO2 (i.e., declining part of an overshooting) scenario. By contrast, MIROC-ES2L (Hajima et al. 2020), which is the new model in the MIROC-ESM series, has shown a non-increase in the TCRE in the decline part of an overshooting scenario (Hajima, personal communication). Thus, whether the TCRE increases in the declining part of overshooting scenarios (i.e., the existence of a scenario dependence) depends on the model. While the quasi-constancy of the TCRE is obviously interesting and helpful, TCRE changes, which have not been the focus thus far, could be important when estimating the remaining carbon budget. Therefore, investigating the magnitude of scenario dependence, and identifying those cases in which the scenario dependence is enhanced, can be important for estimating the remaining carbon budget for given temperature targets.

The present study aims to (1) confirm the magnitude of the scenario dependence of the TCRE and (2) investigate in which cases the TCRE has a high (or low) scenario dependence. To achieve these aims, a single model ensemble perturbing physical and biogeochemical parameters using an Earth system model of intermediate complexity (EMIC) is conducted to investigate the magnitude of the scenario dependence of the TCRE and its correlations with each parameter. In addition, the mechanism of scenario dependence and its impact on the estimation of the remaining carbon budget are discussed.

2 Methods

The present study uses the MIROC-lite/JUMP-LCM (Tachiiri et al. 2010, 2013) model, which contains a fully coupled core general ocean circulation model, with a biogeochemical nutrition–phytoplankton–zooplankton–detritus (NPZD) component, and an energy moisture balance atmosphere model (Oka et al. 2011). In addition, MIROC-lite uses a simple, one-layer-bucket land model to close the water budget by considering precipitation (snow and rain), snowmelt, evaporation, and runoff. However, the spatial distribution of climatic variables is insufficiently represented by this model—especially precipitation, which is difficult to represent with a two-dimensional atmosphere model. This was compensated for when calculating the land carbon uptake by using the loosely coupled Sim-CYCLE model, a process-based land ecosystem model that includes vegetation and soil (Ito and Oikawa 2002) in combination with once-yearly data interactions with the global mean surface air temperature (GMSAT). First, the annual mean GMSAT for the year in question was calculated from the output of MIROC-lite. In addition, a dataset was prepared from the output of the 1-ppa scenario using the medium-resolution version of MIROC3.2 (K-1 model developers 2004). The annual mean GMSAT data for the year in question were then examined, and the most similar datum to that of the MIROC-lite output was called from the prepared MIROC3.2 dataset. This information was then input to Sim-CYCLE. The variables called to run the vegetation model were cloud cover, precipitation, specific humidity, shortwave radiation, air temperature, and soil temperature at depths of 10 cm and 200 cm. This combined model has been employed in international model inter-comparison projects such as EMIC-AR5 (Weaver et al. 2012; Eby et al. 2013; Zickfeld et al. 2013) and ZEC-MIP (MacDougall et al. 2020).

For the model, 512 parameter sets were developed by using the Latin hypercube method to perturb the 11 parameters listed in Table 1 (i.e., parameter sets were generated with equal distribution in the 11-dimension parameter space). These parameters include the Gent–McWilliams thickness parameter that measures the eddy-induced velocity of the ocean (Gent and McWilliams 1990). In addition, a freshwater flux adjustment is necessary in order to realistically represent the meridional overturning circulation. For measuring the marine carbon cycle response, changing the photosynthesis parameters proved ineffective at varying the carbon cycle sensitivity to pCO2 and temperature. Therefore, an atmosphere–ocean gas exchange coefficient was implemented as a perturbation parameter.

Table 1 Nature and ranges of the parameters perturbed in the present study

The ranges of the parameters perturbed were selected on the basis of the original ranges given in Table B1 of Tachiiri et al. (2013) rather than the final range in Table 1 of the same work. Whereas Tachiiri et al. (2013) considered the representability of the Coupled Climate Carbon Cycle Model Inter-comparison Project (C4MIP) models (Friedlingstein et al. 2006), the ranges of parameter perturbation in the present work are based on previous studies. In particular, the range of horizontal diffusivity is determined on the basis of Collins et al. (2007) and Roach et al. (2018), and is significantly narrowed. The reasons for the varying ranges of the parameters are given in more detail in the Supplementary Information (SI_201015_1126.docx). In addition, while a temporally unchanged but spatially varying wind speed (namely, the average value from the ERA5 reanalysis for 1979–2020 (Copernicus Climate Change Service 2017)) is used to calculate marine ocean uptake, the coefficient of aerosol forcing is not applicable here since the scenarios only consider the pCO2 change.

As indicated in Fig. 1 and Table 2, the following concentration pathways were applied to the study scenarios:

  1. (a)

    I-1%: 1 ppa increase to 2 × CO2 (dark blue line)

  2. (b)

    I-0.5%: 0.5 ppa increase to 2 × CO2 (green line)

  3. (c)

    I-0.25%: 0.25 ppa to 2 × CO2 (red line)

  4. (d)

    I-2%: 2 ppa to 2 × CO2 (light blue line)

  5. (e)

    I-4%: 4 ppa for 40 years, up to around 4.5 × CO2 (purple line)

  6. (f)

    D-1%: 1 ppa decrease from 2 × CO2 after (a) (black line)

Fig. 1
figure 1

Plots of pCO2 levels against time (years) for the various study scenarios: I-1% (dark blue line); I-0.5% (green line); I-0.25% (red line); I-2% (light blue line); I-4% (Purple line); D-1% (black line)

Table 2 Study scenarios

In addition, a biogeochemically coupled run was performed for the I-1% scenario and a zero-emission experiment was initiated from the end of the I-1% scenario.

Before the experiments, as spin-up, all members of the atmosphere-ocean and land ecosystem parts were run for 3000 and 2000 years, respectively. The TCRE for each scenario was the average for years 21–280 for the I-0.25% scenario, 21–140 for the I-0.5% scenario, 21–70 for the I-1% scenario, 21–40 for the I-2% scenario, 15–24 for the I-4% scenario, and 1–50 for the D-1% scenario. Years 1–20 and the final 20 years of the D-1% scenario were removed to avoid the instability caused by very small cumulative carbon emissions. The average of the 10-year period for the I-4% scenario was centered on year 19, where the CO2 concentration level was doubled. Note that the TCRE for the D-1% scenario was calculated from the temperature anomaly and cumulative carbon emissions at the start of the precedent I-1% scenario so a combination of I-1% and D-1% could provide an example of an overshooting scenario.

To investigate the correlation between the CO2–carbon cycle and climate–carbon cycle feedback to the TCRE and its scenario dependence, a biogeochemically coupled run of the I-1% scenario was performed in which ecosystems were subject to an increasing pCO2 but the CO2 did not function as a greenhouse gas. The TCRE can be decomposed to Eq. (2):

$$ \alpha /\left(1+\beta +\alpha \gamma \right) $$
(2)

where α[K/GtC], β [GtC/GtC], and γ [GtC/K] are the linear transient climate sensitivity to pCO2 increase, the CO2–carbon cycle feedback, and the climate–carbon cycle feedback, respectively (Friedlingstein et al. 2006; Gregory et al. 2009). Hence, these are assumed to directly influence the change in the TCRE. Here, in calculating β and γ, the effects of temperature change on the biogeochemically coupled run is confirmed to be small, and β and γ are simply calculated using Eqs. (3) and (4) (Arora et al. 2020):

$$ \upbeta =\frac{\Delta {C}_b}{\Delta CA} $$
(3)
$$ \upgamma =\frac{\Delta {C}_f-\Delta {C}_b}{\Delta {T}_f} $$
(4)

where ΔCb and ΔCf are the total change in land and ocean carbon storage in the biogeochemically and fully coupled runs, respectively; ΔCA is the change in the atmospheric carbon, and ΔTf is the temperature change in the fully coupled run. Here, values for the 70th year were used.

To examine the relationship between the scenario dependence of the TCRE and the ZEC, a ZEC experiment was performed. In contrast to ZEC-MIP (MacDougall et al. 2020), which starts the zero-emission period when the cumulative emissions reach their given levels (i.e., 750, 1000, and 2000 GtC), the start of the zero-emission period in the present study is the 70th year of the I-1% scenario. Moreover, the ZEC values up to the temperature peak (the same as Eq. (1)), the 25th (the average of the 20th–29th), and the 50th (the average of the 45th–54th) years are used for the analysis, whereas ZEC-MIP uses the averages of 10 years centering on the 25th, 50th, and 90th.

For three parameters, namely the equilibrium climate sensitivity (ECS), the maximum photosynthetic rate, and the specific leaf area, values near the edges of the ranges are given smaller weights, as shown in Figure S1. Thus, for the ECS, the weight linearly increases from 0 to 1 as the ECS increases from 1 to 2; the weight then remains at 1 throughout the 2–4.5 K range before linearly decreasing to 0 beyond that range (Figure S1a). The coefficient of maximum photosynthetic rates of between 2 and 3 have smaller weights (Figure S1b), and specific leaf areas of between 1.5 and 2.5 have smaller weights (Figure S1c). The total weight is calculated as the product of the weights for these three parameters. For the members which develop negative cumulative carbon emissions during the experiments (as happens for two members in the D-1% scenario), the end of the period for calculating the average TCRE is moved to 20 years prior to the year in which the TCRE becomes negative.

3 Results and discussion

3.1 Scenario dependence of the TCRE

The distribution in the TCRE change for each scenario relative to that of the I-1% scenario is presented in Fig. 2. Here, many members under the I-0.5% scenario (green line) have similar (slightly larger) TCRE values than those under the I-1% scenario; the average ratio of the TCRE to that of I-1% is 1.05, the peak is in the 1.0–1.1 bin, and 50% of the members are within the range of ± 5% difference from those of the I-1% scenario. Under the I-2% scenario (light blue line), the distribution is shifted toward small TCREs with a peak in the 0.9–1.0 bin; the average is 0.91, and 12% are within the range of ± 5% difference from those under the I-1% scenario. The other peaks occur at 0.7–0.8, 1.0–1.1, and 1.1–1.2 for the I-4% (purple line), I-0.25% (red line), and D-1% (black line) scenarios, respectively, with corresponding averages of 0.79, 1.09, and 1.12. The corresponding fractions of members within the range of ± 5% difference from the I-1% scenario are 0.5%, 29%, and 21%, respectively. These results clearly indicate that the TCRE tends to increase under scenarios with a slowly increasing pCO2 and decrease under scenarios with a rapidly increasing pCO2. In the case of D-1%, with a 1-ppa decease from the non-equilibrium 2 × CO2 level, the TCRE is even larger than that of the I-0.25% scenario. The large similarity between I-1% and I-0.5% is consistent with the results of MacDougall (2017), who concluded that the TCRE is only close to constant (i.e., < 5% change) at a pCO2 of 0.3–1.2% p.a. at 400 ppm.

Fig. 2
figure 2

Weighted probability distributions of the relative change in the TCRE from that of the I-1% scenario. The numbers in parentheses indicate the average among all members and the percentages of the members within a ± 5% difference from the value under the I-1% scenario: I-0.5% (green line); I-0.25% (red line); I-2% (light blue line); I-4% (purple line); D-1% (black line; starting from a non-equilibrium state)

3.2 Correlations of the parameters to the scenario dependence of the TCRE

The coefficients of correlation between each parameter and the TCRE for the I-1% scenario are presented in Table 3. Also shown are the correlations between each parameter and ratio of the TCRE value of each alternative scenario to that of the I-1% scenario. Unsurprisingly, the correlations of many parameters display opposite signs depending on whether the specific scenario has a larger or smaller rate of increase in the pCO2 relative to that of the I-1% scenario. This is a direct result of using the ratio of the specific scenario TCRE to that of the I-1% scenario as the metric.

Table 3 Correlation between each parameter, and TCRE in the I-1% scenario and the relative TCRE change from the I-1% scenario

An examination of Table 3 indicates that the ECS parameter is especially strongly correlated to the TCRE value of the I-1% scenario and to the ratio of the TCRE values of all other scenarios to that of the I-1% scenario. Among the other parameters, relatively large correlations with the scenario dependence of the TCRE are observed for the plant respiration parameter, the coefficient for marine carbon uptake, the maximum photosynthetic rate, and the specific leaf area, particularly for rapid pCO2 increase scenarios. Among the physical parameters, the vertical diffusivity of the ocean displays a significant correlation. For the D-1% scenario, the maximum photosynthetic rate and plant respiration parameter display relatively large correlations with the scenario dependence of the TCRE.

For comprehensive sensitivities, the CO2–carbon cycle feedback (β) is fairly strongly correlated with the TCRE under the I-1% scenario and is relatively unchanged under all other scenarios except for the D-1% scenario. The climate–carbon cycle feedback (γ) displays a weak correlation to the TCRE and its scenario dependence, although the unweighted case shows a statistically significant correlation in slow pCO2 increase scenarios (not shown). For fast pCO2 increase scenarios, the correlation of the ECS is slightly reduced and β instead displays a slightly larger correlation relative to that in slow scenarios. When the TCRE is large under the I-1% scenario, the TCRE is seen to increase under slower scenarios (D-1%, I-0.25%, and I-0.5%) and to decrease under faster scenarios (I-2% and I-4%), with very large absolute coefficients of correlation (0.81–0.94 and 0.84–0.95 for Pearson’s and Spearman’s coefficients, respectively).

The correlation between the ECS and the TCRE under the I-1% scenario, and between the ECS and the ratio of the TCRE value of each alternative scenario to that of the I-1% scenario, are indicated by the scatter plots in Fig. 3. The TCRE is seen to be almost independent of scenarios that are slower than the I-1% scenario for models with an ECS of 2–3 K (i.e., the ratio of the TCRE to that of the I-1% scenario was is close to 1) and has a value of approximately 1.5–1.6 K/1000 GtC (see also Figure S2a, which shows that the deviation from the I-1% scenario is smallest for an ECS of 2–3 K under slow scenarios). For the I-2% and I-4% scenarios, a smaller ECS results in a lesser scenario dependence. The ECS correlation is seen to decrease under fast scenarios, and this is attributed to the relatively small variation in the temperature anomaly among the members in rapid scenarios, particularly when the ECS is high (Fig. 3e, f).

Fig. 3
figure 3

The correlation between the ECS and a the TCRE under the I-1% scenario, and bf the relative change in the TCRE from that of the I-1% for the other scenarios: b D-1% (starting from a non-equilibrium state), c I-0.25%, d I-0.5%, e I-2%, and f I-4%. The red, blue, and black circles indicate members with weights of > 0.5, 0.1–0.5, and < 0.1, respectively. Note that each sub-figure has a different y-axis

The ECS dependences of the temperature and cumulative carbon emissions for the various scenarios are presented in Figure S3a-d. The large correlation between the ECS and the TCRE change results from the combined influence of temperature and cumulative carbon emissions. High ECS models have previously been shown to have low transient climate responses (TCRs) to the ECS ratio (Millar et al. 2015, Tsutsui 2017, 2020); i.e., these models are prone to longer responses. This tendency is confirmed in the present study (Figure S3f). Thus, members of high ECS respond more slowly than members of low ECS, generating a smaller temperature rise for a given level of emission (smaller TCRE) in fast pCO2 increase scenarios than in slow pCO2 increase scenarios. This higher realized fraction (i.e., ratio to the value when the system is equilibrated for a given pCO2 condition) for warming in slow scenarios becomes even more remarkable in the declining part of overshooting scenarios (Fig S3c); in overshooting, the system is warmed up when it experiences a high pCO2 condition, causing high realized warming fraction in the following period. Actually, a similar effect occurs for CO2 too (Fig S3d), but as the temperature is more sensitive, the TCRE rise in the declining part of the overshooting scenarios becomes the highest among the scenario considered here. This is regarded as the mechanism by which the ECS correlates with the scenario dependence of the TCRE. Indeed, Tachiiri et al. (2015) also showed a positive relationship between the ECS and TCRE, and this proposed mechanism probably explains their observation that stabilization of the pCO2 resulted in enhanced TCRE increases during the pCO2-increasing period for high-TCRE models and TCRE decreases for low-TCRE models. In brief, the results presented herein can explain why the use of a high ECS model (MIROC-ESM, 4.7K; Andrews et al. 2012) revealed an increase in the TCREs under overshooting scenarios (Tachiiri et al. 2019), whereas the relatively low climate sensitivity of the MIROC-ES2L model (ECS = 2.7 K; MacDougall et al. 2020) did not give that outcome.

For each scenario, the deviations (root-mean-square errors) of the TCREs from that of the I-1% scenario are related to the ECS, CO2–carbon cycle feedback (β), and climate–carbon cycle feedback (γ) in Figure S2, while the correlations between the TCRE and the β and γ are indicated by the scatter plots in Figures S4 and S5, respectively. Thus, a larger β is seen to result in a relatively small scenario dependence for the TCRE, as does a large (small negative) γ value (Figure S2). Here, it is worth noting that the members with γ<−300 PgC/K have small ECS, thus resulting in small weights and a small temperature rise.

Care should be exercised on the limitations of EMICs or differences in the behaviors between EMICs and full ESMs. Tokarska et al. (2016) found that the TCRE was stable for the Representative Concentration Pathway 8.5 scenario, but tended to decrease for EMICs, when the cumulative carbon emissions were large. In addition, MacDougall et al. (2020) identified a positive relationship between the ZEC and the ECS in ESMs, but this tendency disappeared when EMICs were added. Moreover, although the correlation between the ECS and ZEC was small, the same study showed a clear positive relationship between these values when different ECS versions of the same EMICs were compared, thus demonstrating that the use of a single-model ensemble (with varying parameters) can provide different results from those of a multi-model study.

3.3 Implications for carbon budget estimation

The values of the weighted 1st–99th TCRE percentiles for each scenario are plotted in Fig. 4. Thus, for the ensemble members of the lowest TCRE values, the scenario dependence of the TCRE is small. Among the 33rd, 50th, and 67th percentiles, which present the remaining carbon budget in SR1.5 (Rogelj et al. 2018), the 33rd percentile is relatively stable except for a slight decrease in the I-4% scenario. For the 50th and 67th percentiles, the TCRE is slightly increased under I-0.25% (by 8% and 13%, respectively), while 12% and 20% increases are observed under the D-1% scenario. For higher percentiles, the scenario dependence of the TCRE is significant, with a tendency toward low values under rapidly increasing pCO2 scenarios, high values under slowly increasing pCO2 scenarios, and even higher values under the decline-after-a-peak (D-1%) scenario. In this experiment, the 50th percentile in the dataset corresponds to the 1.6 K/1000GtC for the I-1% scenario. The 95% ranges (2.5–97.5 percentiles) are 0.9–1.6, 1.0–2.0, 1.0–2.3, 1.0–2.6, 0.9–3.0, and 0.9–3.3 K/TtC for the I-4%, I-2%, I-1%, I-0.5%, I-0.25%, and D-1% scenarios, respectively. Although the sample size is small, a linear regression is attempted in Figure S6 in order to estimate the TCRE using the ECS and the rate of increase in the pCO2. This enabled the derivation of Eq. (5):

$$ \mathrm{TCRE}\left[\mathrm{K}/\mathrm{TtC}\right]=\left(-0.091\ \ln\ (I)+0.26\right)\ast \mathrm{ECS}\left[\mathrm{K}\right]+\max \left(0,\min \left(0.77,0.22\ast \ln (I)+0.76\right)\right) $$
(5)
Fig. 4
figure 4

Weighted values of the 1st–99th percentiles of the TCRE under each scenario. D-1% scenario starts from a non-equilibrium state

where I is the rate of increase in pCO2 [ppa].

The remaining carbon budget calculated from the TCRE for at the 33rd, 50th, and 67th percentile of the I-1%, I-0.5%, and I-0.25% scenarios, along with the ratio of the remaining carbon budget to that of the I-1% at the same percentile, are presented in Table 4. The two remaining carbon budget values quoted for each scenario are those obtained before (left) and after (right) subtracting the estimated contribution of unpresented Earth system feedback (100 GtCO2 = 27.3 GtC; Rogelj et al. 2018). The I-2% and I-4% scenarios are not included because they are too fast for practical carbon budget estimation for temperature targets such as 1.5 °C and 2 °C. Meanwhile, the D-1% scenario is excluded because Eq. (1) is not applicable to it (Rogelj et al. 2019). Here, the method of Tachiiri et al. (2019) was used to estimate the remaining carbon budget for the 1.5 °C target, i.e., using TlimThist = 0.53 °C in Eq. (1), assuming that the non-CO2 greenhouse gas contribution was 30%, and that the 2011–2017 emissions were 79.4 PgC. The TZEC was expected to be scenario dependent but since only one scenario was assessed in this study, its value was assumed to be zero, following MacDougall et al. (2020). In this study, for the I-1% scenario, the ZEC is found to be 0.07, 0.11, and 0.15 °C for the 33rd, 50th, and 67th percentile, respectively. Due to the scenario dependence of the TCRE, the remaining carbon budget could be reduced by 18% and 22% before and after subtracting the estimated contribution of unpresented Earth system feedback for the most extreme case (i.e., the 67th percentile I-0.25% scenario).

Table 4 Remaining carbon budget (RCB) for the 1.5 °C target based on the TCREs for the I-1%, I-0.5%, and I-0.25% scenarios

The ratio of the TCREs under the D-1% and I-1% scenarios is used as an indicator of the scenario dependence of the TCRE, and its relationship to the ZEC at the peak temperature is plotted in Fig. 5. Similar plots at 25 and 50 years after starting the zero-emission period are presented in Figure S7. In each case, a positive relationship between the ZEC and the scenario dependence of the TCRE can be observed. This is expected because the warming in the zero-emission period is an indicator of scenario dependence in terms of the temperature and the TCRE, as shown in Figure S8. Note that the ratio of the TCRE under the D-1% and I-1% scenarios is used in this analysis because the extreme scenario of D-1% generates a large spread in TCRE.

Fig. 5
figure 5

The relationship between the ZEC and the scenario dependence of the TCRE. The ZEC is presented at the peak temperature after the cessation of CO2 emissions. Ratio of TCREs for D-1% (starting from a non-equilibrium state) to that of I-1% scenarios is used as an indicator of the scenario dependence. The red, blue, and black circles indicate members with weights of > 0.5, 0.1–0.5, and < 0.1, respectively

Although the ZEC was considered in Eq. (1), the TCRE was assumed constant (having no scenario dependence). Nevertheless, as shown in Fig. 5, a non-zero ZEC is associated with the scenario dependence of the TCRE. Thus, for cases with no scenario dependence for the TCRE, and no ZEC, the ZEC term can be omitted from Eq. (1) and carbon budget estimation becomes relatively trivial. However, when a ZEC and scenario dependence of the TCRE occur, inclusion of the ZEC term may be insufficient. It is necessary to consider what specific scenario was used to estimate the TCRE and, if the scenario is significantly different, the TCRE value must be adjusted. This is particularly important when the ECS is large.

The pathways for maintaining the 2 °C target were previously analyzed with a realistic magnitude of overshoot by Tokarska et al. (2019) to indicate identical cumulative CO2 emissions for both overshooting and non-overshooting pathways. They concluded that no corrections are needed for ambitious mitigation scenarios with the low levels of overshoot presented in SR1.5. The present study confirms their output (which is available at http://data.iac.ethz.ch/Tokarska_et_al_2019_Overshoot_UVicESCM/, and was accessed on May 8, 2020) and found that the TCRE is slightly larger in the negative emission period, moderate in the slow emission period, and small in the rapid emission period. More importantly, however, the differences in the TCREs between the scenarios and periods were small (< 5%), which is presumed to lead to almost identical cumulative carbon emissions for both pathways. For the experiments described herein, the relationship between the ECS and the scenario dependence of the cumulative carbon emission for the 1.5 °C target shown in Fig. 6 indicates that this outcome only occurs when the ECS is low. By contrast, when the ECS is high, an overshooting scenario could generate a significant reduction in the carbon budget. Thus, more research is needed to determine whether a correction is required in carbon budget estimation methods for overshooting pathways.

Fig. 6
figure 6

Relationship between the ECS and the scenario dependence of the cumulative carbon emission for the 1.5°C target: D-1% to I-1% (black), D-1% to I-0.5% (red), and D-1% to I-0.25% (blue). Members with a low ECS (approximately < 2.5 K) do not reach the 1.5 °C warming level from the initial state

4 Conclusions

Considering the importance of the TCRE in estimating the remaining carbon budget, the magnitude of scenario dependence and its relation to physical and biogeochemical parameters were investigated via ensemble experiments using an EMIC. The results indicated that, on average, the TCRE decreases in scenarios with a rapid rate of CO2 increase, and increases in scenarios with a slow rate of CO2 increase. Moreover, the declining part of an overshooting scenario was shown to generate an even larger TCRE than the slow CO2 increase scenario.

Among the parameters considered, the ECS displayed by far the strongest correlation with the scenario dependence of the TCRE. This is attributed to the longer response time in the high ECS model and can explain why the use of a high ECS model (MIROC-ESM, 4.7 K; Andrews et al. 2012) revealed an increase in the TCRE under overshooting scenarios (Tachiiri et al. 2019), whereas this was not observed when using the MIROC-ES2L (ECS = 2.7 K). Here, when the ECS was high, the ECS correlation was reduced for fast scenarios, and the correlations of the biogeochemical parameters were strengthened compared to those in slow scenarios.

In the slower scenarios, the TCRE increase was found to be more significant when the ECS was high. As such, the remaining carbon budget for the 1.5 °C target could be significantly reduced in the slow pCO2 scenarios. For example, the 67th percentile values for the I-0.25% scenario before and after subtracting the estimated contribution of the unpresented Earth system feedback were reduced by 18% or 22% relative to the case using the value for the I-1% scenario, respectively, thus indicating that the scenario dependence of the TCRE could cause non-negligible changes in the remaining carbon budget. For an overshooting scenario (i.e., I-1% + D-1%), the cumulative carbon emission to 1.5 °C warming was found to be reduced compared to the slow increasing (e.g., I-0.25%) scenarios for high ECS cases. In the present study, the scenario dependence of the realized warming fraction demonstrated that the ZEC is also likely to be scenario dependent, and could therefore be included in calculating the carbon budget. However, the incorporation of the ZEC term into the carbon budget estimation remains as future work because the specifics are beyond the scope of the present study.

The quasi-constant TCRE concept is both scientifically interesting and practically helpful, and there are many other uncertainties in carbon budget estimation (e.g., GMST observation). However, an investigation of the subtle differences in the TCRE is meaningful for both attaining a deeper understanding of the climate–carbon cycle system and for more practical carbon budget estimation.

Notably, the present study was conducted using a single EMIC. The difference in the behaviors between ESMs and EMICs is sometimes non-negligible and the results of a single-model ensemble can differ from those of a multi-model ensemble. The methods used to simplify the EMIC and select the range of parameter variation could also affect the results. However, considering the significance of the impact, and the tendency for Coupled Model Inter-comparison Project Phase6 (CMIP6) models to have larger ECS values than previous versions (Zelinka et al. 2020; Meehl et al. 2020), an investigation into the scenario dependency of the TCRE in state-of-the-art ESMs could be very important future work.