The social cost of carbon dioxide under climate-economy feedbacks and temperature variability

Our setup recognizes that deterministic assessments of the SCCO₂ carry only very limited information. PAGE-ICE uses Monte Carlo sampling of over 150 parameter distributions (Yumashev et al 2019) to provide distributions of the results. All results presented use 50 000 Monte Carlo draws (and 100 000 for PAGE09, using @RISK within Excel), with draws taken from the same superset to be able to compare SCCO₂ distributions across models. The PAGE-ICE model has been translated into the Mimi modeling framework, using the same validation process as for Mimi-PAGE (Moore et al 2018). Model code and documentation are available from the GitHub repository, https://github.com/openmodels/MimiPAGE2020.jl.

To estimate the marginal damage of an additional tonne of CO₂, PAGE-ICE is run twice, with one run following the exogenously specified emission pathway and the second run adding a CO₂ pulse. The SCCO₂ is then calculated as the difference in global equity-weighted damages between those two runs divided by the pulse size, discounted to the base year (2015). Equity weighting of damages follows the approach by Anthoff et al (2009a) using a mean (minimum, maximum) elasticity of marginal utility of consumption of 1.17 (0.1–2.0), and equity-weighted damages are discounted using a pure time preference rate of 1.03% (0.5%, 2.0%). For all our results, we rely on a 75Gt pulse size in the first time period of PAGE-ICE (mid-2017–2025), representing an annual pulse size of 10 Gt CO₂. In this setup, we found that the choice of pulse size can have an effect on the SCCO₂ estimates, and we explore these sensitivities in appendix A.2, alongside a general sensitivity analysis of PAGE-ICE's model parameters in appendix A.3.

We provide results for a selection of scenarios across climate outcomes and socioeconomic developments, based on the Tier 1 scenarios of ScenarioMIP (O'Neill et al 2016). SSP1-1.9 and SSP1-2.6 are generally well aligned with the Paris Agreement (IPCC 2018, Rogelj et al 2018), while SSP5-8.5 features very high radiative forcing and rapid GDP growth. The 'middle-of-the-road' socioeconomic pathway SSP2 is combined with the emission pathway RCP4.5 scenario for the central values presented in this analysis. Since reference SSP scenarios are provided until 2100, we extend these by making the assumption that regional GDP per capita and population growth rates in the different model regions converge toward the global mean. We implement this by defining a region's post-2100 growth rate as $g_{r, t} = (1 - DR - CR) \cdot g_{r, t-1} + CR \cdot \bar{g}_{t-1}$, where CR and DR are the respective universal rates of convergence and decay, and $\bar{g}$ is the global mean growth rate. The subscripts $r,t$ represent region and time period, respectively. We fit the convergence and decay rates based on SSP growth rates up to 2100 using Stan, a Bayesian MCMC system (appendix A.5).

We implement the persistence parameter following Estrada et al (2015) into the growth system of Burke et al (2015) such that: $\mathrm{GDP}_{r,t} = \mathrm{GDP}_{r,t-1} \cdot (1 + g_{r,t} - \rho \cdot \gamma_{r,t-1})$, where g is the growth rate, γ represents the contemporary economic damages in % of GDP returned by the market damage function and ρ specifies the share of economic damages that persist and thus alter the growth trajectory in the long run. Note that this approach nests the extreme assumptions of zero persistence usually made in CB-IAMs and the assumption of full persistence from the empirical literature (Burke et al 2015, Ricke et al 2018) for ρ = 0 and ρ = 1, respectively. For consistency with the underlying RCP-SSP scenarios, we also rescale greenhouse gas emissions proportionally to the change in GDP, such that emission intensities of economic output remain unchanged. To estimate the distribution of ρ, we calculate the ratio of the long-run marginal impact of temperature on economic growth to the immediate marginal impact using the historical panel data by Burke et al (2015) (appendix B.3). For estimating the long-run marginal impact, we follow the literature in estimating regression models featuring temperature lags (Dell et al 2012, Burke et al 2015), and we additionally apply a more recent approach using low-pass filtering by Bastien-Olvera and Moore (2021) for robustness checks (appendix B.4).

For our aggregated analysis, we do not include any sectoral distinction for the persistence parameter ρ since PAGE-ICE models economic damages on the aggregated macro level using the damage function by Burke et al (2015). By taking one global value distribution for the ρ parameter, we assume that the persistence of damages is similar in each region modeled. Reduced significance for non-global estimates of partial persistence hinders the use of region-specific parameters. For all central results, damage persistence also remains constant over time following from the limited evidence for successful adaptation to date (Burke et al 2015, Burke and Emerick 2016). We explore these simplifications in detail in appendix B.5, where we consider a case that sees economically persistent damages only in lower-income regions, following suggestive evidence for regional heterogeneity, for instance due to higher climate vulnerability and reduced adaptive capacity (Byers et al 2018, Andrijevic et al 2020). Furthermore, we discuss the impacts of possible future adaptation reducing the persistence of temperature impacts on GDP (see also appendix B.6).

To allow for modeling annual temperature anomalies, the inputs to the temperature and GDP modules of PAGE are annualized by a combination of exponential interpolation for population, sea-level rise, and abatement costs, and linear interpolation for GDP growth rates. Annualized temperatures are subsequently used for calculating the damages on an annual basis for the non-market, market, and discontinuity modules. Sea-level rise and consequent damages are modelled to depend on the global climatic mean temperatures modeled in PAGE-ICE, rather than annual on mean temperatures.

There is no strong evidence for the increase or decrease in global annual temperature variability with global mean temperature increases (Huntingford et al 2013, Sippel et al 2015). Therefore, we model the magnitude of interannual temperature variability to be constant over time, based on the analysis of regional means in a recent spatially complete dataset (Ilyas et al 2017) that is suitable for the analysis of temperature variability because its underlying statistical simulation approach allows for a more adequate approximation of the expected local variability and its uncertainty (Beguería et al 2016). We model global (1a ) and regional (1b ) temperature variability as follows:

$$\begin{align} &T_{g,t+1} \sim \mathcal{N}(\alpha_g + \beta_g \Theta_{gt} + \gamma_g T_t,\,\sigma_{g}^{2}) \end{align} \tag{ 1a }$$

$$\begin{align} &T_{r,t+1} \sim \mathcal{N}(\alpha_r + \beta_r T_{rt} + \gamma_g T_{g,t+1},\,\sigma_{r}^{2}) \end{align} \tag{ 1b }$$

here, T is the realized annual temperature, Θ is the mean temperature realized in the PAGE climate module prior to variability, σ is the standard deviation of internal variability, and the subscripts $g,r,t$ represent global, regional, and time (in years). In the absence of auto-regressive feedback, $\alpha_g = \alpha_r = \gamma_g = \gamma_r = 0$, $\beta_g = 1$, and $\beta_r = AF_r$, the regional amplification factor from PAGE-ICE. The temporal standard deviations σ are derived from annual spatially aggregated regional mean temperatures. These regional temperatures are derived from median gridded observations of 10 000 statistical ensemble members (Ilyas et al 2017), which are area-weighted and linearly detrended for a 30 year climatic period. Further information including uncertainties in these observations are found in appendix C.1. With auto-regression, we fit all parameters to historical temperatures using least-squares regression, and applying a LOESS of global temperatures as $\Theta_{gt}$ (appendix C.2). As a robustness check, we include a simple alternative implementation of variability (appendix C.3).

Assuming the SSP2-4.5 scenario, the mean value of the SCCO₂ in 2020 (expressed in $2015 USD throughout) from the previous version of the model, PAGE09, is $158 (figure 2). We note this is greater than the SCCO₂ mean value of $106 previously reported for PAGE09 (Hope 2013) due to the use of the SSP2-4.5 scenario instead of SRES A1B (22% increase) and the change in monetary units (23% increase from inflation). We break down the remaining changes from PAGE09 to PAGE-ICE into a series of key steps (figure 2). (a) Updating the CO₂ model, climate sensitivity, and other climate parameters to conform to IPCC AR5, increases the mean SCCO₂ to $217. (b) Applying a fat-tailed (Gamma) distribution to sea-level rise increases mean SCCO₂ to $228. (c) The permafrost feedback, represented in PAGE-ICE through both CO₂ and CH₄ cycles, increases it to $244. (d) Modeling non-linear transitions in the sea ice and land snow albedo feedback slightly reduces the SCCO₂ to $239. (e) We reduce the size, threshold, and lag of discontinuity impacts, since Arctic feedbacks and catastrophic sea level rise are now explicitly modeled, reducing the mean SCCO₂ slightly to $213. (f) PAGE-ICE has a longer simulation period (until 2300), which further increases the mean SCCO₂ to $239. (g) The baseline income levels, used for calibrating vulnerability, are updated to those observed in 2015, increasing mean SCCO₂ to $245. (h) Adjusting the vulnerability of regions in PAGE-ICE to market and non-market damages to be common and equal to the vulnerability of the European Union region increases the mean SCCO₂ to $281. Finally, (i) by changing damages to empirical estimates (Burke et al 2015) imposed as non-persistent changes to GDP, the final mean SCCO₂ for PAGE-ICE is calculated as $307 per tonne.

Zoom In Zoom Out Reset image size

The uncertainty range is considerable, with an interquartile range of $147–$349, and a 5%–95% range of $82–$831, based on 50 000 Monte Carlo simulations. Nevertheless, for the SSP2-4.5 scenario a clear increase in the SCCO₂ is shown, with PAGE-ICE having a mean SCCO₂ that is double the value for PAGE09. Strikingly, this increase can be attributed almost entirely to increased social costs in the Global South. The mean SCCO₂ attributed to the Global North sees no notable change, mainly due to the regional empirical damage estimates in PAGE-ICE, which allows for benefits of global warming in cooler regions and comes with higher damages in warmer regions compared to the temperature damage function in PAGE09.

Based on panel regression analysis using historical temperature and GDP data (Burke et al 2015), we estimate that on average 50.1% (34.5%–69.0%, interquartile range) of GDP damages persist as impacts on growth (appendix B.3). Introducing this empirically derived level of damage persistence affects both the size of global GDP and its regional distribution by vast amounts, dominating the effect of any other modeling choice. This impact is a result of deteriorating growth trajectories, particularly in Africa, South Asia, and Latin America. In contrast, cooler regions such as Eastern Europe and Northern Asia ('Russia+') and the European Union (EU) experience additional growth. Through the growth rate feedback, these GDP impacts accumulate over time and thus the introduction of damage persistence leads to a considerable redistribution of income to the detriment of poorer regions, an effect which is exacerbated over time. To identify the implications of different levels of persistence on GDP outcomes, we present results both for various fixed degrees of persistence from 0% to 100% (figures 3(a) and (b)) and for our empirically derived persistence distribution (figures 3(c) and (d)).

Zoom In Zoom Out Reset image size

For this persistence distribution, mean global GDP in 2100 is 30% lower than the SSP2-4.5 growth path without persistence ¹⁰ . By 2200, the mean growth loss due to the persistence of climate impacts under SSP2-4.5 exceeds the GDP growth rate in the whole Global South, causing economic contraction. However, the share of runs of negative growth on a global scale as early as 2100 exceeds 8%. By 2200, the 'India+', 'Africa+' and 'Latin America' regions reach income levels associated with extreme poverty in 12%, 11% and 4% of Monte Carlo simulations, respectively. In contrast, mean economic output in the 'Russia+' region is 52% higher with the possibility of Russia+ dominating the global economy near the end of the model horizon for high persistence levels.

If a mere 10% of economic damages were to persist via reduced growth, we find a fifteenfold increase in the mean SCCO₂ (figure 3(e)). For higher levels of persistence, however, model estimates of the SCCO₂ are found to decrease. This effect is caused by the total climate damages in PAGE being capped to the exogenous statistical value of society ($6.13 × 10¹⁶ following Yumashev et al 2019). For high levels of damage persistence, the projected damages go beyond this predefined limit of the model, such that an additional tonne of CO₂ leads to no increase in overall damages and the SCCO₂ becomes zero, while the effects of climate change are actually catastrophic in these runs (Weitzman 2014).

Although this is the case in more than half of the Monte Carlo simulations under our empirical persistence distribution, we still estimate a mean SCCO₂ of [Line:91 Col:11763372 which exceeds the estimate without damage persistence by an order of magnitude. Performing a sensitivity analysis shows that this estimate is highly sensitive to equity weighting parameters and also varies across methodologies to estimate persistence (appendices B.2 and B.4). The substantial increase in the SCCO₂ is primarily driven by the adverse implications for the Global South. Therefore, the resulting SCCO₂ is very similar if only lower-income regions are assumed to suffer from persistent economic damages (appendix B.5).

Both current temperature-induced economic damages and the persistence of damages are derived from historical data which do not show an increase in resilience to temperature shocks over time (Burke et al 2015, Burke and Emerick 2016). As a result, our approach implicitly assumes that future adaptation will remain at the levels observed hitherto for economic damages. Yet, while there is no clear evidence for adaptation to the market impacts of climate change, it seems unlikely that countries would not invest considerable resources in reducing the lasting effects of temperature increases if actual damages in the (far) future were as big as simulated. Therefore, we further explore the effects of adaptation reducing the damage persistence by a constant annual rate, similar to Moore and Diaz (2015). As figure 4 shows, if current persistence is around 50% (as suggested by the mean of our empirical distribution), then the level of persistence would need to decrease by more than 2% per year via adaptation efforts to reduce the model's SCCO₂ to below $2500. As discussed above, the combination of high persistence and low to no adaptation leads to a substantial share of the Monte Carlo runs producing SCCO₂ estimates of zero, meaning that the methodological decision behind PAGE-ICE's damage cap increasingly dominates the model results in the lower-right corner of figure 4. Yet, even if current persistence is only 20% and decreases by an annual 0.5%, more than 5% of Monte Carlo reach PAGE-ICE's damage cap and the SCCO₂ exceeds $5000. For our empirical persistence distribution, we find that a reduction by at minimum 3% per year would be needed to see a SCCO₂ of less than $600 (appendix B.6). In other words, this would require lowering the persistence of temperature-related economic impacts by half within less than 25 years.

Zoom In Zoom Out Reset image size

We first implement stochastic annual regional temperatures to the base PAGE-ICE model, and then extend this to the version with persistent climate damages. Variability produces small changes in the means and medians of SCCO₂ distributions for all RCP-SSP scenario combinations (figures 5 and 6), with a mean increase of $21 (6.5%) for SSP2-4.5. Accounting for temperature anomalies in model runs leads to a widening of realized temperature ranges, leading to more extreme annual temperatures occurring more frequently, both in the near term and in the long run (figures 5(a) and (b)). This captures better the combination of model uncertainty and internal variability of the climate system for annual global mean temperatures. These wider ranges result in higher mean damages from the non-linear market and non-market damage functions and a greater likelihood of discontinuity damages.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Adding variability also widens the SCCO₂ uncertainty ranges. Most of the additional spread of the SCCO₂ distributions is explained by the interaction of temperature variability with the convex empirical temperature-damage function used for market damages. The introduction of this symmetric temperature variation widens tails at both ends of the SCCO₂ distributions, with the possibility of SCCO₂ values being negative due to regionally beneficial temperature draws in the short term. Introducing the temperature anomalies to an annual version of PAGE with constant damage persistence does not result in an increase in the mean SCCO₂ because of the damage cap.

For the SSP1-1.9 scenario, all model specifications without growth effects produce lower values of the mean SCCO₂ compared to SSP1-2.6 and SSP2-4.5. This scenario features the same socioeconomic projections as SSP1-2.6, but with stronger emissions abatement. PAGE with annual variability shows most clearly how more stringent climate mitigation leads to reduced tail risks. For all scenarios, except SSP1-1.9, including the possibility of warm and cold years and periods disproportionately pushes the 75th percentile further from the median than it affects the 25th percentile, while the 5th and 95th percentile are affected more symmetrically. Thus, the consideration of temperature variability increases extremes on both ends, as well as correcting the most likely values upward (figure 6).

For all the PAGE-ICE-based model specifications, we find higher mean values for SSP2-4.5 than for SSP5-8.5. While the projected climate forcing and subsequent impacts are higher for the RCP8.5 scenario, SSP5 projects higher GDP projections than SSP2, leading to a stronger discounting of future damages, producing lower SCCO₂ values in PAGE-ICE. Moreover, because radiative forcing is logarithmic in concentration (Shine et al 1990), marginal impacts are reduced for higher emission trajectories. For SSP1-2.6, climate change-induced damages are much lower compared to SSP2-4.5, leading to significantly lower SCCO₂ values if no damage persistence is considered. The SCCO₂ values for model versions with growth effects for different scenarios are driven by model runs reaching the damage cap, explaining the lower values for pathways with higher emissions when damages are likely to be higher.

To estimate the marginal damage of an additional tonne of CO₂, PAGE is run twice, with one run based on the exogenous emission scenario and the second run adding a CO₂ pulse in the first model time step in 2020, which ranges from mid-2017 to 2025. This requires a decision regarding the time and the size of the said pulse. While a pulse of merely one tonne of CO₂ captures a truly marginal effect, larger pulse sizes are typically used (Otto et al 2013, Rose et al 2017). The results displayed in figure A2 show that the SCCO₂ is more stable for pulses that are not truly marginal. While small pulses are too variable and sensitive for the specific model and pulse size setting, larger pulse sizes do not accurately capture marginal damages and could see unexpected slight increases in the SCCO₂ as seen for PAGE-ICE around 100 Gt. This increase, which is still well within the confidence interval, might be attributed to the triggering of a discontinuity effect. To strike a balance, all our SCCO₂ estimates are based on an annual pulse size of 10 Gt of CO₂ in model year 2020, totalling a 75 Gt CO₂ pulse over the period modeled.

Zoom In Zoom Out Reset image size

To estimate population and income levels past 2100, we fit a model to the available SSP data and extrapolate it. Both population and income use the same model, defined in terms of growth rates. The model postulates that the changes in growth rates are explained by a rate of convergence and a rate of decay.

The model is as follows:

$$\begin{align*} \textrm{Growth}_{i,t} & = (1 - \beta - \delta) \textrm{Growth}_{i,t-1}\nonumber\\ & \quad + \delta \textrm{MeanGrowth}_{t-1} \end{align*}$$

where i indexes the region, t indexes years, and $\textrm{MeanGrowth}_{t-1} = \sum_i \frac{\mathrm{Population}_{i,0}}{\sum_i \mathrm{Population}_{i,0}} \textrm{Growth}_{i,t-1}$. Above, δ is the rate of convergence, and β is the decay rate.

As SSP data are not available for every year, fitting the exact expression above requires a model with dynamics. We use a two-step approach, fitting the model using Stan, a computational Bayes system. The first step uses the available data directly, fitting the following

$$\begin{align*} & \textrm{Growth}_{i,s} \sim \mathcal{N}(\textrm{Growth}_{i,s-1} (1 - \Delta t (\beta + \delta))\nonumber\\ & \quad + \textrm{MeanGrowth}_{s-1} \Delta t \delta, \sigma_i) \end{align*}$$

where s indexes time-steps, $\Delta t$ is the number of years between time-steps, and country i has uncertainty σ_i . We apply a prior that both β and δ are between 0 and 0.5.

Next, we fit the full model, using the results of the simplified model to improve the Bayesian model convergence. In this case, for a given MCMC draw of β and δ, we calculate the entire time series:

$$\begin{align*} & \widehat{\textrm{Growth}}_{i,t} \sim \mathcal{N}\left(\widehat{\textrm{Growth}}_{i,t-1} (1 - \beta - \delta)\right.\nonumber\\ & \quad \left.+ (\widehat{\textrm{Growth}}_{i,t-1} \cdot w_i) \delta, \sigma_i\right) \end{align*}$$

starting with $\widehat{\textrm{Growth}}_{i,0}$ as known from the SSP dataset.

The probability evaluation is over both the performance of the fit and the priors:

$$\begin{align*} \textrm{Growth}_{is} &\sim \mathcal{N}\left(\widehat{\textrm{Growth}}_{i,t(s)}, \sigma_i\right) \\ \beta &\sim \mathcal{N}\left(\mu_\beta, \sigma_\beta\right) \\ \delta &\sim \mathcal{N}\left(\mu_\delta, \sigma_\delta\right) \\ \log{\sigma_i} &\sim \mathcal{N}\left(\mu_{\sigma,i}, \sigma_{\sigma,i}\right) \end{align*}$$

where $\mu_\cdot$ is the mean estimate for the corresponding parameter, and $\sigma_\cdot$ is the standard deviation across its uncertainty. The prior for σ_i is defined as a log-normal, centered on the mean of the estimates of log σ_i .

The estimates for each SSP are shown below (table A3), with visualizations of time series data for SSP2 and SSP5 in figure A5.

Zoom In Zoom Out Reset image size

In addition to the results presented in the main text, we have run all different model versions for four separate scenarios (table A4).

The PAGE model is composed of several modules for increased functionality and increased flexibility for future development and specific applications, and is based on Mimi-PAGE (Moore et al 2018) and PAGE-ICE (Yumashev et al 2019). To introduce annual temperature variability, an annual temperature model is required, and hence all downstream and related upstream modules involved in the climate-economy feedback were annualized. This was done by exponentially interpolating between PAGE-ICE time steps for population, sea level, and abatement costs, and linearly interpolating GDP and emissions growth rates. These choices are made to stay as close to the original design as possible. An annual model comes with reduced modeling artefacts that plague models with multi-year time steps. In non-annual models, rather than having smooth damage distributions, one finds multimodal damage distributions for binary damage elements, for instance in discontinuity damages modules where damages can only be triggered in the chosen model years, rather than any year modeled. Our annual model reports SCCO₂ distributions that differ from the original model. This is due to changes in the discontinuity sector and the model design of PAGE09 where each time step year covers half of the time to the next analysis year.

In PAGE-ICE, economic growth is exogenously specified by:

$$\begin{equation} \mathrm{GDP}_{r,t} = \mathrm{GDP}_{r,t-1} \cdot (1 + g_{r,t}) \end{equation} \tag{ B1 }$$

where $g_{r,t}$ is the SSP-specific growth rate for region r at time step t. Previous work has introduced growth effects as a share of damages affecting capital stocks or total factor productivity, but as such factors are not explicitly modeled in PAGE, we implement a persistence parameter (Estrada et al 2015) into the parameterized growth system (Burke et al 2015) such that:

$$\begin{equation} \mathrm{GDP}_{r,t} = \mathrm{GDP}_{r,t-1} \cdot (1 + g_{r,t} - \rho \cdot \gamma_{r,t-1}) \end{equation} \tag{ B2 }$$

where ρ specifies the share of economic damages γ that affect the growth rate. If, for instance, market damages equal 2% of GDP, then for ρ = 0.5, economic growth decreases by one percentage point. This approach has three advantages. Firstly, it introduces only one additional parameter thus facilitating the specification of its Monte Carlo distribution. Nonetheless, it can easily be adjusted to new research about region-specific vulnerabilities or hypotheses about future developments by differentiating ρ over regions and time. Secondly, it nests both level effects and the growth system by Burke et al (2015) for ρ = 0 and ρ = 1, respectively, and hence covers a range of approaches applied in previous studies. Thirdly, it mirrors the share of damages inflicted upon productivity or capital stocks in previous papers (Moyer et al 2014, Dietz and Stern 2015, Moore and Diaz 2015) which facilitates comparisons. It is important to note, however, that the γ returned by the damage function in Burke et al (2015) is expressed in % of GDP per capita, not in % of GDP. By applying γ directly to the GDP growth rate and leaving population exogenous, we assume that the adverse GDP per capita impact of a temperature stems only from a reduced GDP and not from an increase in population, in line with indications that a temperature shock would result in a minor decrease in population (Carleton et al 2020).

Damage persistence can have catastrophic impacts on the welfare of regions, to the extent that some regions approach zero consumption. Even if this is not the case, counterfactual GDP losses, defined as the difference between counterfactual GDP levels in the absence of climate change and GDP levels actually realized, might still exceed a region's consumption. Both outcomes would lead to infinite equity-weighted damages in utility terms in PAGE. To see this, note that damages are equity-weighted based on a region's consumption level (Anthoff et al 2009a) relative to the EU in 2015 such that:

$$\begin{equation} \tilde{D}_{r,t} = \frac{C_\mathrm{EU,2015}^\eta}{(1-\eta)} \cdot \left[C_{r,t}^{1-\eta} - (C_{r,t} - D_{r,t})^{1-\eta} \right] \end{equation} \tag{ B3 }$$

where $\tilde{D}_{r,t}$ is equity-weighted damages per capita, $C_{r,t}$ and $D_{r,t}$ are unweighted consumption and damages per capita and η is the elasticity of marginal utility of consumption of the underlying iso-elastic welfare function. If $D_{r,t}$ approaches $C_{r,t}$, this expression yields infinitely large values for $\tilde{D}_{r,t}$ and becomes numerically infeasible once $D_{r,t}$ exceeds $C_{r,t}$.

We therefore introduce two boundaries. First, we assume that a region's per capita consumption cannot fall below the World Bank's current threshold for extreme poverty of $1.90 (2011 PP) per day, converted to $2015 (PAGE-ICE's monetary unit) using the World Development Indicators' GDP deflator time series. As a result, per capita consumption in the model cannot drop below $740.65 (per year) which equals about 17% of the 2015 per capita consumption for the poorest region in PAGE-ICE, India+. Furthermore, we assume that only damages which equal up to 99% of a region's per capita consumption are subject to equity-weighting, following Dietz (2011). Damages that exceed this threshold are still taken into consideration but receive no equity weights. Overall damage estimates are highly sensitive to the level of such a threshold (Dietz 2011), a finding corroborated by our sensitivity analysis (see figure B1). In general, total damages increase with increased levels of equity weighting. However, the SCCO₂ Monte Carlo mean actually decreases because for higher total damages, more Monte Carlo runs reach PAGE-ICE's upper bound for overall damages and thus produce SCCO₂ estimates of zero.

Zoom In Zoom Out Reset image size

The two boundaries explained above are likely to cause considerable discontinuities in results. For example, the marginal damages of CO₂ will abruptly decrease once a region reaches the subsistence level of per capita consumption because additional temperature rises cannot depress economic growth further. Similarly, the equity-weighted marginal damages of an additional tonne of CO₂ will be extremely high for a region with damages just below the 99% threshold but will fall abruptly once the threshold is reached. To smooth out these effects, we introduce these boundaries as logistic paths such that the consumption and equity-weighted damages converge against the threshold without fully reaching it. This convergence system is expressed by the following equation:

$$\begin{equation*} \tilde{y} = \begin{cases} y & \textrm{if}\ y \geqslant \epsilon\\ -\theta + \dfrac{2(\theta - \epsilon) \cdot \mathrm{exp}(b(y - \epsilon))}{1 + \mathrm{exp}(b(y - \epsilon))} &\textrm{if}\ y \lt \epsilon \end{cases} \end{equation*}$$

where y is the original variable subject to threshold θ and $\tilde{y}$ is the adjusted variable converging asymptotically against θ. $\epsilon$ is a value close to θ which triggers the converging process, and b is defined as $b = \dfrac{2}{\theta - \epsilon}$.

This system ensures both that $\tilde{y}$ converges against θ as y goes to infinity and that the derivative of $\tilde{y}$ with respect to y equals 1 at $y = \epsilon$, thus avoiding breaks. Figure B2 illustrates the convergence system for the equity-weighting threshold discussed above, which, in our default settings, has an upper bound at θ = 0.99 and convergence starting at a neighborhood value of  = 0.9. For the per capita consumption threshold, equals 1.5 times the subsistence level.

Zoom In Zoom Out Reset image size

Figure B3 shows the impact of the convergence system on the SCCO₂ results. In general, the introduction of the convergence system decreases the SCCO₂, but the effect is modest compared to other model parameter settings (see for instance figure B1).

Zoom In Zoom Out Reset image size

The previous literature provides two general approaches to determining a point estimate for the persistence of market damages. In some studies, persistence is set to a low, but somewhat arbitrary value to explore growth effects without making strong assumptions regarding the magnitude of persistence (Dietz 2011, Moyer et al 2014). The more empirical approach relies on fixed effects regressions of economic growth on various variables, including temperature (Dell et al 2012, Burke et al 2015). A generic equation for such a model would be

$$\begin{equation} g_{i,t} = \beta \mathrm{Temperature}_{i,t} + \gamma \cdot X_{i,t} + \alpha_{i} + \delta_{t} + \epsilon_{i,t} \end{equation} \tag{ B4 }$$

where $g_{i,t}$ denotes the growth rate of country i at time t and $X_{i,t}$ denotes a vector of additional independent variables with the vector of coefficients γ. If not all of the damages indicated by β persist over time and the economy eventually returns to its original growth trajectory, then a negative GDP impact at time t is followed by a compensating impact in the subsequent years. This can be conceptualized by introducing a lag into the model equation:

$$\begin{align} g_{i,t} & = \beta_{1} \mathrm{Temperature}_{i,t} + \beta_{2} \mathrm{Temperature}_{i,t-1}\nonumber\\ & \quad + \gamma \cdot X_{i,t} + \alpha_{i} + \delta_{t} + \epsilon_{i,t}. \end{align} \tag{ B5 }$$

If damages are not fully persistent over time, then signs of β₂ will be positive. If, on the other hand, damages only reach their full extent over time (e.g. because some damages to infrastructure is exacerbated over time), the sign of β₂ will be negative as well. Of course, there is no imminent reason to assume that the full effects of a temperature shock are realized after only one year. Therefore, researchers typically include further lags into the equation and estimate the long-run marginal impact of temperature on growth (Dell et al 2012, Moore and Diaz 2015), which can identify damage persistence ρ based on its ratio to the immediate marginal impact. Note that for a model with k lags, such a ratio is defined as follows:

$$\begin{equation} \rho = \frac{\sum\limits_{j = 1}^k \beta_{j}} {\beta_{1}}. \end{equation} \tag{ B6 }$$

In the case of pure level effects, all immediate impacts of temperature shocks are compensated for over time and the numerator of ρ = 0 becomes zero. If, however, after k years, the long-run marginal impact is only half the immediate impact of a temperature shock, we can conclude that only 50% of economic damages persist over time. Yet, unlike the simple model given in equation (B4), the regression model by Burke et al (2015) which informs PAGE's market damage function is quadratic in temperature:

$$\begin{align} g_{i,t} & = \beta_{1,0} \cdot \mathrm{Temperature}_{i,t} + \beta_{2,0} \cdot \mathrm{Temperature}_{i,t}^2\nonumber\\ & \quad + \gamma \cdot X_{i,t} + \alpha_{i} + \delta_{t} + \epsilon_{i,t}. \end{align} \tag{ B7 }$$

This means that the marginal impact of temperature depends on the initial temperature. Hence, there is no universal ratio of long-run to immediate impacts for all regions but for a model with k lags, ρ is defined as follows:

$$\begin{equation} \rho = \frac{\sum\limits_{j = 0}^k \beta_{1,j} + 2\beta_{2,j} \cdot \mathrm{Temperature}_{i,t} } {\beta_{1,0} + 2\beta_{2,0} \cdot \mathrm{Temperature}_{i,t}}. \end{equation} \tag{ B8 }$$

However, the variation of ρ can be understood as a numerical, rather than theoretical, issue. To see this, note that in the case of a quadratic regression model, both the long-run and the immediate marginal impacts are linear functions of temperature. The immediate marginal impact, which is the denominator of ρ, approaches zero around 13°C for the coefficients from Burke et al (2015), leading to a damage persistence of $\pm \infty$ in the proximity of this temperature (see figure B4). Therefore, there is little use in aggregating the ratio of marginal impacts for different temperatures, e.g. by using population weights. Asymptotically, however, the ratio becomes numerically stable, as illustrated by figure B4. For this reason, we estimate the damage persistence as the ratio of long-run to immediate marginal impacts as the temperature goes to (minus) infinity, as indicated by the lagged regression model by Burke et al (2015).

Zoom In Zoom Out Reset image size

The question remains as to how many temperature lags one should include in the regression model. Ideally, we would include as many lags as possible to avoid a bias in our ρ estimate due to omitted future rebounds or, conversely, intensification of damages. In practice, however, every lag added to the model comes with a considerable amount of variance, such that the ρ estimate loses its meaning if too many imprecisely estimated lag coefficients are included. For this reason, we rely on the specification with only one lag, for which the quadratic lag coefficient is significant at the 5%-level, while the linear lag coefficient is not (see table B1). As the signs of the first temperature lag are opposite to the immediate temperature impact, the resulting damage persistence equals only 52.8% of the immediate impact.

The numerator and denominator of ρ are linear combinations of asymptotically Gaussian random variables and thus asymptotically Gaussian themselves. As the ratio of two dependent, asymptotically Gaussian variables, ρ does not, however, have a Gaussian distribution. To derive a parameter distribution, we randomly sample the regression coefficients displayed in table B1 from a Multivariate Gaussian distribution based on their variance-covariance matrix and compute our persistence estimate for each of them. We then use this empirical distribution of the persistence parameter ($N = 1 \times 10^{06}$) for Monte Carlo simulations in PAGE.

While we choose a model with only one lag because of the aforementioned bias-variance trade-off, it is important to note that including a second lag would strongly reduce a point estimate to 27.8%. A model with four lags even suggests a persistence of only 10.3%. Such lower values might be explained by two effects. Either the economy takes several years after a temperature shock to rebound, implying that our persistence estimate has an upward bias, or the reduced persistence only mirrors the additional noise introduced by further lags. Confidence intervals in table B1 show that with each additional lag, the estimate's distribution shifts towards zero and widens considerably. Including a second lag would extend the 90% confidence interval for ρ to $-41\%$, implying that current climate damages could push the economy to a trajectory of considerably increased growth in the longrun, which seems at odds with the literature on climate change impacts. However, there is, as of yet, no definitive answer as to how many lags one should include. Therefore, we explore the robustness of our approach by adopting another methodology in the following section.

Table B1. Results for the baseline regression model by Burke et al (2015) for various numbers of temperature lags. Each regression model further contains precipitation, country, and year fixed effects as well as country-specific time trends. Following Burke et al (2015), the regression model containing k temperature lags also features k precipitation lags. Here we only report the coefficients for temperature variables that factor into our estimates of damage persistence.

	(1)	(2)	(3)	(4)	(5)	(6)
	GDPpc Growth	GDPpc Growth	GDPpc Growth	GDPpc Growth	GDPpc Growth	GDPpc Growth
Temperature	0.0127$^{***}$	0.0136$^{***}$	0.0106$^{**}$	0.00 949$^{**}$	0.00 933$^{*}$	0.00 920$^{*}$
	(3.36)	(3.64)	(3.00)	(2.65)	(2.49)	(2.50)
Temperature squared	−0.000 487$^{***}$ (−4.11)	−0.000 517$^{***}$ (−4.33)	−0.000 456$^{***}$ (−3.97)	−0.000 441$^{***}$ (−3.75)	−0.000 446$^{***}$ (−3.62)	−0.000 459$^{***}$ (−3.87)
L.Temperature		−0.00 674	−0.00 413	−0.00 549	−0.00 578	−0.00 459
		(−1.56)	(−1.10)	(−1.35)	(−1.43)	(−1.36)
L2.Temperature			−0.00 613	−0.00 638	−0.00 681	−0.00 698
			(−1.63)	(−1.76)	(−1.86)	(−1.95)
L3.Temperature				−0.00 143	−0.000 885	−0.000 993
				(−0.57)	(−0.31)	(−0.38)
L4.Temperature					−0.00 111	−0.00 214
					(−0.42)	(−0.82)
L5.Temperature						0.00 176
						(0.37)
L.Temperature squared		0.000 244$^{*}$ (2.01)	0.000 206 (1.81)	0.000 240 (1.95)	0.000 245 (1.97)	0.000 228$^{*}$ (2.03)
L2.Temperature squared			0.000 123 (1.14)	0.000 126 (1.22)	0.000 133 (1.28)	0.000 146 (1.43)
L3.Temperature squared				−0.0000 509 (−0.59)	−0.0000 779 (−0.78)	−0.0000 882 (−0.92)
L4.Temperature squared					0.000 101 (1.01)	0.000 143 (1.35)
L5.Temperature squared						−0.0000 658 (−0.43)
Resulting ρ	100%	52.81%	27.82%	28.36%	10.25%	20.81%
5th Monte Carlo percentile	—	0.55%	−73.40%	−101.57%	−143.72%	−94.44%
95th Monte Carlo percentile	—	91.34%	87.77%	106.52%	94.79%	93.10%
N	6584	6519	6398	6277	6155	6031
bic	−19 806.5	−19 677.9	−19 377.4	−19 125.8	−18 744.0	−18 400.4
ll	10 127.4	10 080.4	9942.9	9829.6	9655.6	9500.5

t statistics in parentheses. Standard errors are clustered at the country level. $^{*}$  p < 0.05, $^{**}$  p < 0.01, $^{***}$  p < 0.001.

Zoom In Zoom Out Reset image size

Given the uncertainties outlined in the previous section, we explore the results under a different methodological approach to estimate persistence following Bastien-Olvera and Moore (2021). Their approach, similar to Burke et al (2015), regresses GDP per capita growth on temperature and precipitation, but it also applies a Butterworth low-pass filter to the climate variables. By doing this, the authors argue that the estimated impact on economic growth is identified only via long-term temperature variation. This leads to near-zero estimates if the temperature-based damages to GDP are in fact non-persistent. Filtering out high frequencies in temperature and precipitation is statistically beneficial because it does not require additional lag coefficients to be estimated and is therefore less susceptible to noise and multicollinearity if longer time horizons are considered. Therefore, Bastien-Olvera and Moore suggest comparing the estimated marginal impacts for different minimum periodicities to explore if growth impacts of temperature shocks intensify, remain constant, or converge over time.

To derive an estimate of persistence which we can compare with our results, we apply such a low-pass filter to the temperature and precipitation observations in the data by Burke et al (2015) for a minimum periodicity of 3, 5, and 10 years, respectively. Following Bastien-Olvera and Moore (2021), we discard countries with less than twice the minimum periodicity ($2 \times 10 = 20$) of subsequent temperature/precipitation observations for the filtering to work. This reduces the overall number of country-year records in the data only slightly, from 6584 to 6535. Figure B6 illustrates the effect of the low-pass filter for selected countries that are closely aligned with the model regions of PAGE-ICE. We then regress GDP per capita growth on the low-pass filtered climate variables using the main specification by Burke et al (2015) and, similar to the lag-based procedure, calculate the implied persistence for the 3 year, 5, year and 10 year filter by dividing the asymptotic marginal impact of the respective regression by the marginal impact of the model without filtering. Note that as this involves comparing the estimates returned by two different regression models (based on either unfiltered or filtered climate variables) the empirical persistence distribution cannot be derived from the variance-covariance matrix of a single regression. To overcome this, we draw with replacement from country clusters to create 5000 bootstrap samples following Burke et al (2015) and estimate the implied persistence for each bootstrap sample. The results of the regression models for different filter settings and the resulting persistence distribution are displayed in table B2, where the resulting persistence ρ for each filter period is calculated as the asymptotic marginal impact of temperature in the respective column, divided by the asymptotic marginal impact according to the coefficients in column (1).

Zoom In Zoom Out Reset image size

Table B2. Regressing GDP per capita growth in the Burke et al (2015) data on low-pass filtered climate variables.

	(1)	(2)	(3)	(4)
	Unfiltered	3 year filter	5 year filter	10 year filter
Temperature	0.012 933$^{***}$	0.012 211$^{**}$	0.009 065	0.009 761
	(3.41)	(2.79)	(1.71)	(0.90)
Temperature2	−0.000 490$^{***}$	−0.000 436$^{**}$	−0.000 299	−0.000 258
	(−4.11)	(−3.24)	(−1.88)	(−0.95)
Resulting ρ	100.00%	89.00%	61.09%	52.72%
5th percentile	—	69.33%	9.36%	−60.29%
95th percentile	—	102.40%	91.38%	118.25%
N	6535	6535	6535	6535
BIC	−19 634.1	−19 632.5	−19 620.8	−19 615.5
ll	10 045.5	10 040.3	10 034.4	10 031.8

t statistics in parentheses. Standard errors are clustered at the country level. Percentiles are estimated via 5000 cluster bootstrap samples. $^{*}$  p < 0.05, $^{**}$  p < 0.01, $^{***}$  p < 0.001.Author's calculations based on the data provided by Burke et al (2015).

Table B2 shows that similar to the results under the lag-based approach used in this manuscript, the implied persistence ρ decreases and confidence intervals (or technically, the interval between bootstrap percentiles) widen. For the 3 year filter, the point estimate for ρ is 89.00% which reduces to 52.72% for the 10 year filter. Still, for the 5 year filter, the coefficients for linear and squared temperature remain statistically significant at the 10% level (since the respective t-statistics exceed 1.7). With 9.36%, the 5th bootstrap percentile is well above zero and only 3.6% of the 5000 bootstrap samples yield a non-positive estimate for persistence. For the 10 year filter, temperature coefficients are not significantly different from zero and the 90% bootstrap interval ranges from $-60.29\%$ to 118.25%. The point estimate happens to be relatively close to the 52.81% based on the first temperature lag we have used in model runs (see previous section). This exploration indicates that for the data used in this study, the newly proposed method (considering longer time horizons for persistence via low-pass filtering climate variables) does not produce results that show substantially different damage persistence behaviour from the insights in our main approach (assessing persistence via lag coefficients).

We also explore how an empirical persistence distribution based on low-pass filtering would affect our findings with respect to the SCCO₂. To this end, we use the bootstrap distribution of ρ by comparing the marginal impacts of unfiltered and 5 year filtered temperature, for which temperature coefficients remain statistically significant at the 10% level. The low-pass filtering approach yields higher persistence levels and hence higher damages. However, as figure B7 illustrates, the Monte Carlo mean SCCO₂ decreases from $3372 under the lag-based approach to $3025 for the distribution using low-pass filtering. The primary reason is that the latter on average yields a higher persistence and hence more Monte Carlo runs in which the maximum amount of damages in the model is reached and the SCCO₂ becomes zero. More specifically, the share of Monte Carlo runs with a SCCO₂ of zero increases from 62% under the lag-based approach to 69% based on low-pass filtering.

Zoom In Zoom Out Reset image size

In our main model, we implement persistence as a global parameter. We recognize that this is a very strong assumption and does not represent heterogeneity across climate zones or income levels. Research on adaptive capacity research has identified many potential barriers to adaptation, which differ strongly depending on the context (Eisenack et al 2014). Thus, one might expect differences in the persistence parameter across countries (and thus across PAGE-ICE regions) due to, for instance, income differences or other social and geographic contexts.

To explore the empirical evidence for this, we repeat the persistence estimation procedure for two groups of countries: rich and poor countries. We follow Burke et al (2015) by assigning each country in the data to one of two subsamples based on whether its 1980 GDP per capita, with conversion based on purchasing power parity (PPP), was above or below the sample median. We then repeat our main lag-based approach to estimate persistence separately for each of these income-based subsamples and display the results in tables B3 and B4.

Since subsampling approximately halves the sample size, overall statistical power is diminished. Standard errors increase considerably for all numbers of lags under consideration. In particular, the coefficients on the first temperature lag are no longer statistically significant for either of the two subsamples. For both income groups, the signs of the 1st lag is opposite to the sign of the immediate temperature impact, suggesting less than full persistence. The point estimates for the persistence parameter ρ in the poorer and richer subsample are 42.48% and 61.97%, respectively. When two or more lags are used, the persistence estimate for the below-median sample becomes negative while the estimate for the richer subsample remains positive ,but much closer to zero. However, Monte Carlo experiments show that estimates for ρ vary considerably and even when using one lag only, the 5th Monte Carlo percentiles ranges below zero for the lower- and higher-income subsamples (at $-146.71\%$ and $-95.71\%$, respectively). Regression tables B3 and B4 thus remain largely inconclusive with respect to whether spatially heterogeneity should be implemented in the ρ parameter. Therefore, we conclude that while there might very well be regional differences in the actual persistence of climate damages, we do not find clear evidence in this data set for applying heterogeneous persistence under the lag-based approach in our model.

To explore the possibilities of regional heterogeneity in persistence further, we also apply the low-pass filtering approach outlined above separately to the income-based subsample of the Burke et al (2015) data and display the results in tables B5 and B6. Note that we remove countries with less than 20 temperature-precipitation observations for the low-pass filtering to work which reduces the initial sample size slightly. Similar to the lag-based approach, confidence intervals for the two subsamples are considerably wider compared to the full sample in table B2, reducing t statistics for the individual temperature coefficients and increasing the uncertainty about the persistence level ρ. However, we see different behavior in filters with larger periodicity. Point estimates for the higher-income countries decrease for a larger filter periodicity (from 69.76% for the 3 year filter to $-2.91\%$ for the 10 year filter in table B6) and even for the 3 year filter, the bootstrap interval is almost symmetric around zero, indicating low evidence for persistence. In contrast, the point estimates for the lower-income subsample range between 73.29% and 112.45% for different filter periods (see table B5). Notably, the poorer subset's 95% bootstrap interval does not include zero for the 3 year filter and even for the 5 year filter, 88.8% of the 5000 bootstrap samples yield a persistence estimate above zero. Therefore, we conclude that unlike the lag-based approach, estimating persistence via low-pass filters provides some suggestive evidence that poorer regions might be more vulnerable to persistent growth impacts of temperature shocks.

To illustrate the implications of regional differentiation for the SCCO₂, we carry out an additional robustness check by applying our main persistence distribution (derived via temperature lag coefficients) only to the four PAGE-ICE model regions with the lowest initial GDP per capita in the 2015 base year, namely China+, India+, Africa+ and Latin America. Simultaneously we thus assume for this run that the EU, Russia+, USA and Other OECD regions do not face alterations from their scenario-based long-term growth trajectory, though they will still suffer from contemporary market damages. Figure B8 illustrates the relatively small effect of this change, in part because Global North regions can see opposing effects of temperature increases. Due to generally slightly lower damages, the overall mean SCCO₂ would be slightly higher because fewer Monte Carlo runs reach the maximum amount of overall damages in the model. The share of runs with an SCCO₂ of zero decreases from 62.3% to 61.7%. Overall, our findings are primarily driven by persistence in regions with lower GDP per capita and do not change strongly if richer regions were less vulnerable to damage persistence (or much quicker to adapt).

Table B3. Results for the baseline regression model by Burke et al (2015) for countries with GDPpc below the median and various numbers of temperature lags. Each regression model further contains precipitation, country and year fixed effects as well as country-specific time trends. Following Burke et al (2015), the regression model containing k temperature lags also features k precipitation lags. Here we only report the coefficients for temperature variables that factor into our estimates of damage persistence.

	(1)	(2)	(3)	(4)	(5)	(6)
	GDPpc Growth	GDPpc Growth	GDPpc Growth	GDPpc Growth	GDPpc Growth	GDPpc Growth
Temperature	0.0236	0.0246	0.0096	0.0048	0.0029	−0.0040
	(1.34)	(1.48)	(0.60)	(0.32)	(0.19)	(−0.28)
Temperature squared	−0.0007 (−1.89)	−0.0007$^{*}$ (−2.08)	−0.0004 (−1.25)	−0.0003 (−0.99)	−0.0003 (−0.87)	−0.0002 (−0.51)
L.Temperature		−0.0150 (−0.88)	−0.0035 (−0.33)	−0.0049 (−0.46)	−0.0043 (−0.43)	0.0009 (0.10)
L2.Temperature			−0.0165 (−0.95)	−0.0180 (−1.01)	−0.0206 (−1.10)	−0.0240 (−1.32)
L3.Temperature				0.0061 (0.65)	0.0049 (0.44)	0.0046 (0.44)
L4.Temperature					0.0021 (0.17)	−0.0020 (−0.18)
L5.Temperature						0.0083 (1.18)
L.Temperature squared		0.0004 (1.18)	0.0002 (0.93)	0.0003 (1.09)	0.0002 (1.06)	0.0001 (0.69)
L2.Temperature squared			0.0003 (0.77)	0.0003 (0.82)	0.0004 (0.94)	0.0004 (1.16)
L3.Temperature squared				−0.0002 (−0.85)	−0.0002 (−0.81)	−0.0002 (−0.83)
L4.Temperature squared					0.0001 (0.31)	0.0002 (0.77)
L5.Temperature squared						−0.0003 (−1.59)
Resulting ρ	100.00%	42.48%	−16.29%	−15.89%	−69.83%	−94.29%
5th Monte Carlo percentile	—	−146.71%	−894.45%	−923.95%	−1246.10%	−1081.44%
95th Monte Carlo percentile	—	135.33%	905.03%	1022.92%	1345.78%	1332.94%
N	3429	3391	3332	3273	3214	3154
bic	−9681.9	−9564.9	−9420.3	−9264.0	−9086.2	−8911.6
ll	5048.5	5006.0	4945.4	4882.9	4805.6	4729.7

t statistics in parentheses. $^{*}$  p < 0.05, $^{**}$  p < 0.01, $^{***}$  p < 0.001.Author's calculations based on the data provided by Burke et al (2015).

Table B4. Results for the baseline regression model by Burke et al (2015) for countries with GDPpc above the median and various numbers of temperature lags. Each regression model further contains precipitation, country and year fixed effects as well as country-specific time trends. Following Burke et al (2015), the regression model containing k temperature lags also features k precipitation lags. Here we only report the coefficients for temperature variables that factor into our estimates of damage persistence.

	(1)	(2)	(3)	(4)	(5)	(6)
	GDPpc Growth	GDPpc Growth	GDPpc Growth	GDPpc Growth	GDPpc Growth	GDPpc Growth
Temperature	0.0085$^{*}$	0.0092$^{*}$	0.0084$^{*}$	0.0077	0.0074	0.0082$^{*}$
	(2.22)	(2.42)	(2.22)	(1.97)	(1.77)	(2.04)
Temperature squared	−0.0003$^{*}$ (−2.16)	−0.0004$^{*}$ (−2.24)	−0.0004$^{*}$ (−2.23)	−0.0004$^{*}$ (−2.18)	−0.0004$^{*}$ (−2.11)	−0.0004$^{*}$ (−2.43)
L.Temperature		−0.0044	−0.0032	−0.0045	−0.0053	−0.0051
		(−0.89)	(−0.66)	(−0.87)	(−0.98)	(−1.16)
L2.Temperature			−0.0067	−0.0065	−0.0062	−0.0069$^{*}$
			(−1.83)	(−1.85)	(−1.80)	(−2.06)
L3.Temperature				−0.0024	−0.0026	−0.0018
				(−0.82)	(−0.79)	(−0.58)
L4.Temperature					−0.0006	−0.0008
					(−0.23)	(−0.27)
L5.Temperature						−0.0015
						(−0.27)
L.Temperature squared		0.0001 (0.82)	0.0001 (0.69)	0.0002 (0.84)	0.0002 (0.94)	0.0002 (1.07)
L2.Temperature squared			0.0002 (1.69)	0.0002 (1.60)	0.0002 (1.40)	0.0002 (1.73)
L3.Temperature squared				−0.0001 (−0.38)	−0.0000 (−0.13)	−0.0001 (−0.44)
L4.Temperature squared					−0.0000 (−0.22)	0.0000 (0.00)
L5.Temperature squared						0.0001 (0.36)
Resulting ρ	100.00%	61.97%	1.21%	12.90%	18.08%	0.22%
5th Monte Carlo percentile	—	−95.71%	−323.86%	−391.63%	−370.90%	−229.34%
95th Monte Carlo percentile	—	147.38%	117.18%	148.72%	159.32%	109.67%
N	3155	3128	3066	3004	2941	2877
bic	−9923.1	−9894.5	−9721.7	−9606.1	−9396.5	−9228.0
ll	5167.0	5168.6	5093.7	5047.3	4953.8	4880.8

t statistics in parentheses. $^{*}$  p < 0.05, $^{**}$  p < 0.01, $^{***}$  p < 0.001.Author's calculations based on the data provided by Burke et al (2015).

Table B5. Regressing GDP per capita growth in the Burke et al (2015) data on low-pass filtered climate variables—below GDPpc median.

	(1)	(2)	(3)	(4)
	Unfiltered	3 year filter	5 year filter	10 year filter
Temperature	0.025 169	0.031 177	0.019 792	0.034 540
	(1.43)	(1.75)	(0.75)	(0.60)
Temperature2	−0.000 749	−0.000 842$^{*}$	−0.000 549	−0.000 806
	(−1.97)	(−2.19)	(−0.98)	(−0.67)
Resulting ρ	100.00%	112.45%	73.29%	107.71%
5th percentile	—	84.12%	−62.78%	−296.33%
95th percentile	—	164.32%	155.72%	340.61%
N	3391	3391	3391	3391
BIC	−9556.7	−9552.6	−9542.3	−9541.1
ll	4985.6	4983.6	4978.4	4977.8

t statistics in parentheses. Standard errors are clustered at the country level. Percentiles are estimated via 5000 cluster bootstrap samples. $^{*}$  p < 0.05, $^{**}$  p < 0.01, $^{***}$  p < 0.001.Author's calculations based on the data provided by Burke et al (2015).

Table B6. Regressing GDP per capita growth in the Burke et al (2015) data on low-pass filtered climate variables—above GDPpc median.

	(1)	(2)	(3)	(4)
	Unfiltered	3 year filter	5 year filter	10 year filter
Temperature	0.008 412$^{*}$	0.006 184	0.004 519	0.000 073
	(2.20)	(1.30)	(0.73)	(0.01)
Temperature2	−0.000 330$^{*}$	−0.000 230	−0.000 135	0.000 010
	(−2.13)	(−1.23)	(−0.58)	(0.03)
Resulting ρ	100.00%	69.76%	40.96%	−2.91%
5th percentile	—	−70.72%	−248.19%	−631.26%
95th percentile	—	97.68%	106.27%	117.74%
N	3144	3144	3144	3144
BIC	−9883.0	−9879.2	−9876.6	−9875.4
ll	5146.8	5145.0	5143.7	5143.1

t statistics in parentheses. Standard errors are clustered at the country level. Percentiles are estimated via 5000 cluster bootstrap samples. $^{*}$  p < 0.05, $^{**}$  p < 0.01, $^{***}$  p < 0.001.Author's calculations based on the data provided by Burke et al (2015).

Zoom In Zoom Out Reset image size

By using historical data to calibrate our estimate of damage persistence and extrapolating it into the distant future, we implicitly assume (in figure 3) that future persistence will equal the persistence levels observed over the past decades, consistent with Burke et al (2015). It is, however, possible that the increased economic damages due to climate change will trigger a rise in adaptation efforts which might not only reduce current damages, but also limit the vulnerability of capital stocks (both physical and human) to climate change. If this is the case, then our assumption of extrapolating previous persistence levels might turn out to be too pessimistic. In the absence of empirical evidence or mechanistic understanding of how temperature damages might persist in future scenarios, we explore the theoretical possibility that damage persistence (ρ) decays exponentially over time following a similar approach adopted by Moore and Diaz (2015) for the DICE CB-IAM. Figure B9 illustrates the impacts of varying decay rates on our SCCO₂ estimate under SSP2-4.5.

An annual decay of 3% or higher greatly reduces the estimated marginal damages of CO₂ emissions, especially later in the model horizon for the regions with the highest temperatures, when the highest climate change-induced economic damages are reached, but the degree of persistence has decayed to low levels. If damage persistence decreases by 4%–5% per year, the SCCO₂ equals $305–312 which is virtually identical to the SCCO₂ for the default PAGE-ICE model with no persistence.

A decay rate of 1% or 2% counterintuitively increases the SCCO₂ compared to our default with zero decay from $3372 to $6595 or $4094 respectively, because of interactions with the cap for overall damages in PAGE-ICE. With low decay rates, fewer Monte Carlo simulations reach the cap, which decreases the share of model runs with a SCCO₂ of zero and thus increases the Monte Carlo mean.

Zoom In Zoom Out Reset image size

The historical temperature dataset used is described in Ilyas et al (2017), which is a variant of HadCRUT4 (Morice et al 2012) and is available at https://oasishub.co/dataset/global-monthly-temperature-ensemble-1850-to-2016. A multi-resolution lattice kriging approach was used to estimate the uncertainties related to limited spatial coverage in HadCRUT4. 10 000 ensemble members together express probabilistic gridded temperatures.

The current body of literature does not conclusively answer whether interannual variability is increasing globally, but observation-based studies appear to indicate little change in global annual mean temperature variability (Alexander and Perkins 2013, Huntingford et al 2013). After Huntingford et al (2013), we compute interannual temperature variability in our dataset as the long-term average of the 11 year standard deviations after detrending the annual temperature anomaly data using a local 11 year running mean. Figure C1 shows the interannual variability for used period 1986–2015, compared to pre-1986 levels. We identify no clear overall positive or negative trend in interannual variabiility from the Ilyas et al (2017) dataset. Contrary to Huntingford et al (2013), we find strongly positive changes over the tropics, especially over Amazonia, whilst we see little increases in variability over Europe. Some of this difference probably arises from the longer records in our analysis, but the more rigorous treatment of uncertainties during interpolation used in this study may be another minor underlying difference (Beguería et al 2016). Nonetheless, even for a longer time series analysis, we do not find global or PAGE-ICE region-specific evidence that justifies increasing or decreasing temperature variability with rising global temperatures in simulations.

To find the level of interannual variability per region, we apply a land-sea mask to the coarse HadCRUT4 grid, followed by multiplying with the cosine of the latitude. Consequently, annual mean temperatures are calculated for the PAGE-ICE regions by spatial aggregation. A 30-year running window is used for linear detrending to remove the global warming signal, before the standard deviation is calculated over that 30 year window, following, for example, Huntingford et al (2013). We use the results for the period 1986–2015, which are found in table C1. The medians of the ensemble, which are provided separately, are used in the main model specification.

Zoom In Zoom Out Reset image size

Temperatures exhibit a considerable degree of autocorrelation across years. Table C2 displays regression models of temperature, showing autocorrelation levels of over 0.3. These are computed using the regional temperature data summarized in table C1. We use column (5) from table C2 as a basis for decomposing annual temperature variability into autocorrelation, smooth trends, and random variability. The smooth trend, which is represented in the regression model with locally estimated scatterplot smoothing (LOESS), is equivalent to the underlying climatic temperature computed in PAGE (PT_g ). Figure C2 shows several realizations of the autoregressive model, which roughly captures both the variability and temporal structure of observed temperatures.

Table C2. Models of annual temperature, accounting for trends and auto-correlation. Columns 1–3 represent different models of the smooth trends behind annual temperatures, using a simple trend (1), CO₂ concentrations (2), and a LOESS of observed GMST (3). Columns (4) and (5) show combinations.

	Dependent variable:
	GMSTt
	(1)	(2)	(3)	(4)	(5)
Year	0.006$^{***}$
	(0.0003)
GMST$_{t-1}$				0.402$^{***}$	0.315$^{***}$
				(0.071)	(0.074)
CO2 (ppm)		0.010$^{***}$		0.006$^{***}$
		(0.0004)		(0.001)
LOESS			1.007$^{***}$		0.693$^{***}$
			(0.033)		(0.080)
Constant	−10.822$^{***}$	−2.867$^{***}$	−0.009	−1.732$^{***}$	−0.005
	(0.628)	(0.112)	(0.014)	(0.226)	(0.014)
Observations	169	169	169	168	168
R2	0.652	0.825	0.847	0.853	0.862
Adjusted R2	0.650	0.824	0.846	0.851	0.860
Residual Std. Error	0.206	0.146	0.136	0.134	0.130

$^{*}$  p < 0.1; $^{**}$  p < 0.05; $^{***}$  p < 0.01.

Zoom In Zoom Out Reset image size

The amount of autoregressivity varies by PAGE region, from 0.012 (Other OECD) to 0.200 (India+) (see table C3). These levels of autoregression are low both because of the inherently greater variability of weather at a local scale and because the contemporaneous global temperature is controlled for. We use the results of column (5) in tables C2 and C3 in a two-step process to add autocorrelated temperatures to PAGE-ICE.

Table C3. Regression results of the autoregressive model by region.

	Dependent variable: Reg. Tempt
	EU	Rus+	USA	Chi+	Ind+	Afr+	Lat	Oth
Reg. Temp$_{t-1}$	0.073	0.060	0.046	0.190$^{***}$	0.200$^{***}$	0.123$^{**}$	0.232$^{***}$	0.012
	(0.069)	(0.059)	(0.062)	(0.057)	(0.058)	(0.062)	(0.064)	(0.053)
GMST	1.243$^{***}$	1.766$^{***}$	1.182$^{***}$	1.010$^{***}$	0.709$^{***}$	0.949$^{***}$	0.652$^{***}$	1.604$^{***}$
	(0.118)	(0.124)	(0.090)	(0.077)	(0.055)	(0.071)	(0.060)	(0.094)
Constant	0.006	0.059	−0.029	−0.008	−0.020	−0.020	−0.005	0.018
	(0.039)	(0.037)	(0.029)	(0.019)	(0.015)	(0.015)	(0.015)	(0.026)
Observations	168	168	168	168	168	168	168	168
R2	0.583	0.764	0.691	0.831	0.795	0.856	0.785	0.827
Adjusted R2	0.578	0.761	0.687	0.829	0.792	0.854	0.783	0.825
Residual Std. Error	0.394	0.363	0.288	0.194	0.154	0.154	0.154	0.259

$^{*}$  p < 0.1; $^{**}$  p < 0.05; $^{***}$  p < 0.01.

In this new version, we compute the global temperature as

$$\begin{align*} RT_{gt} \sim \begin{cases} \mathcal{N}\left(PT_{gt}, \sigma_g\right) & \text{if } t = 1\\ \mathcal{N}(\alpha_g + \beta_g PT_{gt} + \gamma_g (PT_{gt} - PT_{g,t-1} + RT_{g,t-1}), \sigma_g) & \textrm{else} \end{cases} \end{align*}$$

where α_g , β_g , and γ_g are the constant, LOESS, and autoregression terms from model (5) in table C2, and σ_g is the corresponding standard error of the residuals. The $PT_{g,t-1}$ variable is the same as PAGE-ANN. The term $PT_{gt} - PT_{g,t-1}$ removes the lagging of temperatures that would otherwise by observed in RT_gt due to the autoregression.

Next, we compute regional temperatures as

$$\begin{align*} RTL_{rt} \sim \begin{cases} \mathcal{N}\left(\textrm{AF}_r RT_{gt}, \sigma_r\right) & \textrm{if}\ t = 1\\ \mathcal{N}\left(\frac{\textrm{AF}_r}{\beta_r + \gamma_r} \left(\alpha_r + \beta_r RTL_{rt} + \gamma_r (\textrm{AF}_r (PT_{gt} - PT_{g,t-1}) + RT_{g,t-1})\right), \sigma_g\right) \end{cases} \end{align*}$$

where $\textrm{AF}_r$ is the PAGE amplification factor. The $\textrm{AF}_r (PT_{gt} - PT_{g,t-1})$ term is again introduced to remove lagging temperatures. The entire regression-based expression is multiplied by $\frac{\textrm{AF}_r}{\beta_r + \gamma_r}$ to scale the regional temperatures to approximate the original PAGE results, $\textrm{AF}_r RT_{gt}$, accounting for any mismatch between the regional temperature dataset and PAGE's regional calibration.

Under this model, we assume that residual temperature variation is uncorrelated across regions, after accounting for GMST and regional auto-regression. Since the regions used in PAGE are chosen for economic coherence, rather than climatological distinctness, we test this assumption by looking at the correlation in model errors across regions from our model.

Zoom In Zoom Out Reset image size

As shown in figure C3, compared to the correlation across regions in the raw temperatures, controlling for either a smoothed GMST predictor or an auto-regressive term reduces average correlation of the residuals across regions by half. The full model, including both a smooth GMST predictor and auto-regressive terms, approximately halves the cross-region correlation again, to a median 0.15 correlation. The final variance of these residuals is also reduced, from 0.22 for raw temperatures to 0.04 for the full model.

As an additional result, we explore the effects of a simpler specification of temperature variability that takes into account observational uncertainties. We specify this alternative model without autogression as follows:

$$\begin{align} &T_{g,t+1} = \Theta_{gt} + \mathcal{N}(0, \mathcal{N}(\mu_{g}^{2} ,\Delta_{g}^{2})) \end{align} \tag{ C1a }$$

$$\begin{align} &T_{g,t+1} = \Theta_{gt} \cdot AF_{r} + \mathcal{N}(0, \mathcal{N}(\mu_{r}^{2} ,\Delta_{r}^{2})) \end{align} \tag{ C1b }$$

where µ represents the standard deviations on the mean annual temperature for a geographical region, and Δ is the standard deviation capturing the uncertainty in the 10 000 member ensemble. Δ of the two is specified exogenously in the Monte Carlo draw for each run, while µ is randomly drawn every year. The impacts of this alternative model on the SCCO₂ are shown in table A4.