Permafrost response and feedback under temperature stabilization and overshoot scenarios with different global warming levels

Permafrost in the northern high latitudes faces significant degradation risks under global warming and this threatens the achievement of global climate goals. This study explores the response and feedback of permafrost and the associated carbon loss under temperature stabilization scenarios, where the global mean temperature stabilizes at various global warming levels, and overshoot scenarios, and where the global mean temperature temporarily exceeds the 1.5 °C warming target. Under the 1.5 and 2 °C stabilization scenarios, permafrost area is projected to decrease by 4.6  (4.5 to 4.7) and 6.6  (6.4 to 6.8) million km², respectively, from a pre-industrial level of 17.0 million km². Corresponding permafrost carbon losses are estimated at 54  (32 to 79) and 72  (42 to 104) PgC, relative to a pre-industrial carbon stock of 484  (383 to 590) PgC. In overshoot scenarios, permafrost area shows effective recovery, with additional losses of only 0.6  (0.3 to 1.1) million km² compared to the 1.5 °C stabilization scenario. In contrast, permafrost carbon loss remains largely irreversible, with additional loss of 24  (4 to 52) PgC compared to the 1.5 °C stabilization scenario. Both stabilization and overshoot scenarios show that additional warming due to permafrost carbon feedback rises with higher global warming levels, and the additional warming in overshoot scenarios becomes most pronounced during the cooling phase. The additional permafrost area loss due to permafrost carbon feedback, which accounts for 5  (2 to 11) % of the total loss, is influenced by both the magnitude of additional warming and the sensitivity of permafrost area to global warming. Moreover, the responses of permafrost area, permafrost carbon, and associated radiative forcing to a broad range of global warming exhibit near-linear relationships under stabilization scenarios. Based on the simulations presented, permafrost carbon feedback is unlikely to initiate a self-perpetuating global tipping process under both stabilization and overshoot scenarios. These findings have significant implications for long-term climate change and mitigation strategies.

1 Introduction

Permafrost soils in the northern high latitudes contain an estimated 1100–1700 Pg of carbon, primarily in the form of frozen organic matter, which is roughly twice the amount of carbon in the atmosphere (Hugelius et al., 2014; Schuur et al., 2015). As the climate warms, both the gradual and abrupt permafrost thaw processes, along with subsequent microbial decomposition, would release carbon dioxide (CO₂) and methane (CH₄) into the atmosphere, thereby amplifying the warming effect (Koven et al., 2011; Feng et al., 2020; Smith et al., 2022). The positive feedback mechanism, combined with the fact that warming rates in the Arctic exceed the global average (Fyfe et al., 2013; Liang et al., 2022; Rantanen et al., 2022), underscores the critical role of permafrost as a key tipping element in the climate system (Armstrong McKay et al., 2022). However, current Earth system models inadequately represent or omit the permafrost carbon processes, which has become one of the largest sources of uncertainty in future climate projections (Schädel et al., 2024). Therefore, researching the permafrost carbon loss and its feedback on the climate system is crucial for accurately assessing climate risks and formulating effective emission reduction strategies.

The Paris Agreement aims to limit global average temperature rise to well below 2 °C above pre-industrial levels, with efforts to keep it below 1.5 °C. Despite these goals, global warming has already exceeded 1 °C and is on track to surpass 3 °C by the end of the 21st century, primarily due to increased anthropogenic CO₂ emissions (Haustein et al., 2017; Hausfather and Peters, 2020). If current emission rates persist, the remaining carbon budgets compatible with the 1.5 °C target will be critically tight and likely exhausted within the next few years (Rogelj et al., 2015; Goodwin et al., 2018; Masson-Delmotte et al., 2018; Forster et al., 2023; Smith et al., 2023). It is unlikely that the 1.5 °C target set by the Paris Agreement will be met (Raftery et al., 2017). However, it might still be achievable after a period of temperature overshoot, by compensating for excessive past and near-term emissions with net-negative emissions at a later time – i.e., through on-site CO₂ capture at emission sources and carbon dioxide removal from the atmosphere (Gasser et al., 2015; Sanderson et al., 2016; Seneviratne et al., 2018; Drouet et al., 2021; Schwinger et al., 2022).

Several existing studies have assessed the climate response to overshoot pathways. Many components of the physical climate system have been identified as reversible, although typically with some hysteresis behavior (Boucher et al., 2012; Wu et al., 2015; Tokarska and Zickfeld, 2015; Li et al., 2020; Cao et al., 2023). In this context, reversibility refers to a partial recovery of climate conditions in an overshoot scenario toward an Earth system state without overshoot. These studies demonstrate that global mean temperature, sea surface temperature, and permafrost area can recover within decades to centuries in response to net negative emissions. Carbon loss from permafrost has been shown to be irreversible on multi-decadal to millennial timescales (MacDougall et al., 2013; Schwinger et al., 2022). The presence or absence of hysteresis effect in the permafrost processes is influenced by multiple factors, including the thermal inertia of permafrost soils, potential shifts in vegetation composition, and the extent to which irreversible permafrost carbon losses are offset by gains in vegetation and non-permafrost (i.e. usual) soil carbon reservoirs (MacDougall, 2013; Schwinger et al., 2022). Furthermore, the soil carbon loss under overshoot scenarios significantly affects the hydrological and thermal properties of soils (Zhu et al., 2019), which in turn modulate the processes involved. The interactions between physical and biophysical processes can potentially stabilize the carbon, water, and energy cycles at distinct post-overshoot equilibria (de Vrese and Brovkin, 2021). Therefore, a temporary warming of the permafrost regions entails important legacy effects and lasting impacts on its physical state and carbon cycle. However, these existing studies have yet to assess permafrost carbon feedback under overshoot scenarios.

Few studies have examined permafrost carbon response and feedback under a long-term climate stabilization scenario, which have received less attention than overshoot scenarios. However, such stabilization scenarios could prove more realistic than overshoot scenarios if carbon dioxide removal ultimately falls short. Most carbon dioxide removal methods have not yet been proven at large scales (Anderson et al., 2023), and sustainability concerns further limit land-based carbon dioxide removal options (Deprez et al., 2024). Multi-model ensemble mean projections suggest that global temperature change following the cessation of anthropogenic greenhouse gas emissions will likely be close to zero in the following decades, although individual models show a range of temperature evolution after emissions cease from continued warming for centuries to substantial cooling (MacDougall et al., 2020; MacDougall, 2021; Jayakrishnan et al., 2024). MacDougall (2021) demonstrated that under a cumulative emission scenario of 1000 PgC, permafrost carbon feedback could raise global temperature by approximately 0.06 °C within 50 years after emissions cease. However, systematic research is still lacking on how permafrost carbon will be released and further amplify global warming under temperature stabilization scenarios, as well as how these processes differ from temperature overshoot scenarios reaching the same peak warming. Additionally, most post-Paris Agreement studies have focused on understanding the impacts of 1.5 and 2 °C global warming levels, with little attention given to the higher global warming levels (Rogelj et al., 2011; Comyn-Platt et al., 2018; King et al., 2024). Assessment of permafrost carbon feedback under various global warming levels would provide a more comprehensive understanding of the critical role of permafrost in the climate system.

This study aims to fill these gaps using an Earth system model of intermediate complexity to systematically assess the permafrost response and feedback under temperature stabilization or overshoot scenarios achieving various global warming levels. The structure of the remainder of this paper is as follows: Sect. 2 introduces the model description, design of temperature stabilization and overshoot scenarios at different global warming levels, as well as the methodology for perturbed parameter ensemble simulations. Section 3 provides the results of these simulations. Finally, Sect. 4 offers conclusions and discussion of our findings.

2 Methods

2.1 Model Description

This study uses the University of Victoria Earth System Climate Model version 2.10 (UVic ESCM v2.10), an intermediate-complexity Earth system model with a uniform horizontal resolution of 3.6° longitude by 1.8° latitude, to simulate permafrost carbon response and feedback under temperature stabilization and overshoot scenarios. The atmospheric component of UVic ESCM is a single layer moisture-energy balance model. The oceanic component incorporates a fully three-dimensional general circulation model, with a vertical resolution of 19 levels, coupled to a thermodynamic and dynamic sea-ice model. A more detailed description of UVic ESCM v2.10 can be found in Weaver et al. (2001) and Mengis et al. (2020).

The terrestrial component of UVic ESCM v2.10 uses the Top-down Representation of Interactive Foliage and Flora Including Dynamics (TRIFFID) vegetation model to describe the states of five plant functional types (PFT): broadleaf tree, needleleaf tree, C₃ grass, C₄ grass, and shrub (Cox, 2001; Meissner et al., 2003). A coupled photosynthesis-stomatal conductance model is used to calculate carbon uptake via photosynthesis, which is subsequently allocated to vegetation growth and respiration. The resulting net primary productivity drives changes in vegetation characteristics, including areal coverage, leaf area index, and canopy height for each PFT. The UVic ESCM v2.10 utilized in this study does not account for nutrient limitations in the terrestrial carbon cycle, leading to an overestimation of global gross primary productivity and an enhanced capacity of land to take up atmospheric carbon (De Sisto et al., 2023). However, the model reasonably represents the dominant PFTs of C₃ grass, shrub, and needleleaf tree at northern high latitudes, although it underestimates vegetation carbon density over this area (Mengis et al., 2020).

The UVic ESCM v2.10 represents the terrestrial subsurface with 14 layers, extending to a total depth of 250.3 m to correctly capture the transient response of permafrost on centennial timescales. The top eight layers (10.0 m) are hydrologically active, while the deeper layers are modeled as impermeable bedrock (Avis et al., 2011). The carbon cycle is active in the top 6 layers (3.35 m), where organic carbon from litterfall, simulated by the TRIFFID vegetation model, is allocated to soil layers with temperatures above 1 °C according to an exponentially decreasing function with depth. If all soil layers are below 1 °C, the organic carbon is added to the top soil layer. The soil respiration is calculated for each layer individually as a function of temperature and moisture, but the respiration ceases when the soil layer temperature falls below 0 °C (Meissner et al., 2003; Mengis et al., 2020). In permafrost regions the model applies a revised diffusion-based cryoturbation scheme to redistribute soil carbon within the soil column. Compared to the original diffusion-based cryoturbation scheme proposed by Koven et al. (2009), the revised cryoturbation scheme calculates carbon diffusion using an effective carbon concentration that incorporates the volumetric porosity of the soil layer, rather than the actual carbon concentration, thereby resolving the disequilibrium problem of the permafrost carbon pool during model spin-up (MacDougall and Knutti, 2016). However, as the UVic ESCM v2.10 only simulates permafrost carbon in the top 3.35 m of soil, the current cryoturbation scheme cannot initiate the formation of Yedoma. As a result, soil carbon stored in deep deposits of Yedoma regions is omitted in our simulations. It is worth noting that the permafrost region for carbon-related variables is defined as the area where permafrost carbon exceeds zero. This differs from the permafrost area, which is defined as the area where soil temperature remains below 0 °C for at least two consecutive years.

In the UVic ESCM v2.10, the usual soil carbon and the permafrost carbon are depicted as two distinct depth-resolved carbon pools within the upper 6 soil layers. Soil carbon that is diffused downward and crosses the permafrost table (i.e. the depth of the active layer) is transformed into permafrost carbon. Conversely, permafrost carbon that is diffused upward and crosses the permafrost table is converted back into usual soil carbon. However, as the active layer deepens, thawed permafrost carbon is not transferred to the usual soil carbon pool; it remains in the permafrost carbon pool and decomposes at a rate distinct from that of usual soil carbon. The coexistence of thawed and frozen permafrost carbon is attributed to the depth-resolved soil carbon scheme. Permafrost carbon can only be decomposed into CO₂, as the UVic ESCM v2.10 does not include a methane production module (MacDougall and Knutti, 2016). The permafrost carbon pool is characterized by four key parameters: (1) a decay rate constant; (2) the available fraction for decomposition, which represents the combined size of both the fast and slow carbon pools subject to decay; (3) a passive pool transformation rate, which governs the rate at which passive permafrost carbon transitions into the available fraction; and (4) a saturation factor used to calibrate the total size of the permafrost carbon pool, which is linked to soil mineral porosity and accounts for the decreasing concentration of soil carbon with depth in permafrost regions (Hugelius et al., 2014). These four key parameters determine the size and vulnerability to decay of the permafrost carbon pool. A more detailed description of the UVic ESCM's permafrost carbon parameterization scheme and its simulated permafrost carbon characteristics can be found in MacDougall and Knutti (2016).

As an Earth system model of intermediate complexity, the UVic ESCM offers relatively low computational costs and serves as an ideal instrument for performing a large number of simulations that are not yet feasible with state-of-the-art Earth system models (Weaver et al., 2001; MacIsaac et al., 2021). In recent years, the UVic ESCM has played a key role in assessments of carbon-climate feedbacks (Matthews and Caldeira, 2008; Matthews et al., 2009; Tokarska and Zickfeld, 2015; Zickfeld et al., 2016) and long-term climate change projection (Ehlert and Zickfeld, 2018; Mengis et al., 2020; MacDougall, 2021). UVic ESCM v2.10 realistically reproduces historical changes in temperature and carbon fluxes (Mengis et al., 2020), compared to the observation-based dataset for temperature change (Haustein et al., 2017) and Global Carbon Project 2018 for carbon fluxes (Le Quéré et al., 2018). Moreover, it has also been validated for simulating permafrost area and permafrost carbon stocks (Mengis et al., 2020).

Figure 1 Temperature stabilization and overshoot scenarios designed through a simple proportional control scheme (Zickfeld et al., 2009) with permafrost carbon module deactivated (corresponding to the NPFC simulations). (a) Cumulative CO₂ emission and (b) temperature change relative to the pre-industrial levels (1850–1900), in the historical and SSP5-8.5 scenario (black), the stabilization (dashed lines) and overshoot (solid lines) scenarios for 1.5 °C (orange), 2.0 °C (brown), 3.0 °C (green), and 4.0 °C (purple) global warming levels.

2.2 Experimental Design

To quantify the impact of permafrost carbon response under different climate conditions, we designed a series of idealized climate scenarios, including four stabilized warming levels (SWL) and three overshoot (OS) trajectories, spanning a long period from 1850 to 2300. In these scenarios, the global mean temperature increases to a maximum of 1.5, 2.0, 3.0, and 4.0 °C above pre-industrial levels (1850–1900). Following this initial warming phase, the scenarios diverge into two groups following distinct temperature trajectories. In the stabilization scenarios, temperature stabilizes at the respective global warming levels and remains steady until 2300, referred to as SWL-1.5, SWL-2, SWL-3, and SWL-4. In contrast, the overshoot scenarios allow temperatures to temporarily reach these warming peaks before entering a cooling phase symmetrical to the warming phase, gradually reducing the global mean temperature to 1.5 °C above pre-industrial levels, where it remains steady until 2300. These overshoot pathways are referred to as OS-2, OS-3, and OS-4, respectively.

A simple proportional control scheme (Zickfeld et al., 2009) was used to derive the CO₂ emissions pathway corresponding to each stabilization and overshoot scenarios, adjusting emissions to align the global mean temperature with the prescribed trajectory. The proportional control equation is given by:

$$\begin{array}{}\text{(1)}& E_{i+\mathrm{1}}=k_{\text{PE}}(T_i-T_{i,\text{goal}}),\end{array}$$

Where E_i+1 represents the carbon emission at year i+1, T_i represents the simulated annual mean global mean temperature at year i, and T_i,goal represents the prescribed global mean temperature trajectory at year i for each scenario. The time-invariant coefficient k_PE is set to 1515 Pg C K⁻¹, ensuring that the control scheme responds neither too quickly nor too slowly to the diagnosed temperature derivation. Using the proportional control scheme, a set of emission trajectories were produced from initial simulations of stabilization and overshoot scenarios using the UVic ESCM v2.10 with the permafrost carbon module deactivated, thereby excluding the influence of permafrost carbon emissions. Specifically, all these initial simulations of stabilization and overshoot scenarios firstly follow a unified emissions trajectory, which is based on historical CO₂ emissions (Friedlingstein et al., 2022) until 2021, and thereafter follow the emissions trajectory of the Shared Socioeconomic Pathway 5-8.5 (SSP5-8.5) (O'Neill et al., 2016) until the highest temperatures of each scenario are approached. As the highest temperatures of each scenario are approached, the proportional control gradually steps in to mediate CO₂ emissions and realize the warming trajectory of each scenario. The CO₂ emission trajectories obtained from these initial simulations were used later to drive formal emissions-driven simulations of stabilization and overshoot scenarios (Fig. 1a). Once these emission trajectories were established, there is no further application of the proportional control in the formal simulations.

Figure 2 Probability distribution functions of the four key permafrost carbon parameters perturbed in the UVic ESCM v2.10 to represent uncertainty in permafrost carbon response. This figure is adapted from MacDougall (2021).

To isolate the contribution of permafrost carbon feedback, two parallel sets of formal simulations were conducted for each scenario, both driven by the same CO₂ emission trajectories obtained from the initial simulations. One set activated the permafrost carbon module and is referred to as PFC simulations, while the other set deactivated the permafrost carbon module and is referred to as NPFC simulations. Since the emission trajectories were obtained from the initial simulations in which the proportional control scheme was applied to achieve the desired temperature pathways, applying them in the formal NPFC simulations can effectively achieve the same temperature trajectories (Fig. 1b). However, applying the obtained emissions in the formal PFC simulations results in any permafrost carbon fluxes being effectively added on top of the applied emissions, thereby causing additional warming. In other words, to achieve the intended climate targets under the same emission pathways, removals equivalent to the permafrost carbon emissions would be required. Therefore, the comparison between the PFC and NPFC simulation sets provides a robust framework to isolate and quantify the additional warming and radiative forcing effects due to permafrost carbon emissions under stabilization and overshoot scenarios. Comparing the PFC simulations with the NPFC simulations, rather than with the initial simulations, eliminates any potential effects of the proportional control scheme. For comparison, two parallel simulations with permafrost carbon module activated or deactivated were also conducted for the high-emissions SSP5-8.5 scenario.

To evaluate the uncertainty in permafrost carbon response under stabilization and overshoot scenarios, we perturbed the four key permafrost carbon parameters following the methodologies of MacDougall and Knutti (2016) and MacDougall (2021). The permafrost carbon decay constant was derived from the mean residence time (MRT) of the slow soil carbon pool at 5 °C in permafrost soils and adjusted to reflect decay at 25 °C using the method proposed by Kirschbaum (2006). The probability distribution function of mean residence time was taken as a normal distribution with a mean of 7.45 years and a standard deviation of 2.67 years from Schädel et al. (2014). The available fraction of permafrost carbon for decomposition was derived from the size of the fast, slow, and passive soil organic carbon pools separately for organic, shallow mineral, and deep mineral soils measured by Schädel et al. (2014). The probability distribution function of available fraction of permafrost carbon for decomposition was described by weighted gamma distributions, with each distribution, respectively describing the probability distribution function of available fraction of permafrost carbon for decomposition in organic, shallow mineral, and deep mineral soils. The passive carbon pool transformation rate was estimated from the ¹⁴C age of the passive carbon pool in midlatitude soils (Trumbore, 2000). Its mean residence time at 5 °C is 300 to 5000 years, yielding a passive carbon pool transformation rate of 0.25 × 10⁻¹⁰ to 4 × 10⁻¹⁰ s⁻¹ after adjustment to 25 °C. The mean residence time of the passive carbon pool was assumed to follow a uniform distribution (MacDougall and Knutti, 2016). For the initial quantity of total permafrost region soil carbon, its probability distribution function was taken as a normal distribution with a mean of 1035 PgC and a standard deviation of 75 PgC, informed by Hugelius et al. (2014). A series of 5000 year sensitivity runs were performed under preindustrial steady conditions with varying saturation factors to determine their relationship with the quantity of the total permafrost region soil carbon. This relationship was then utilized to tune the total permafrost region soil carbon pool, ensuring alignment with observational data (Hugelius et al., 2014). Figure 2 illustrates the probability distribution function for each perturbed parameter. MacDougall and Knutti (2016) and MacDougall (2021) additionally perturbed two physical climate parameters controlling climate sensitivity and Arctic amplification, but they are not perturbed in this study due to their limited influence on global mean temperature in stabilization scenarios.

The Latin hypercube sampling method (McKay et al., 1979) was used to explore the effects of parameter uncertainty on projections of permafrost carbon change. In this study, the probability distribution function of each key permafrost carbon parameter was divided into 25 intervals of equal probability. One value was randomly selected from each interval for a given parameter, and then randomly matched with values of the other three key parameters selected in the same manner to generate parameter sets. This sampling procedure was repeated 10 times, resulting in 250 unique parameter sets (i.e., 250 model variants). For each parameter set, the UVic ESCM v2.10 was first run through a 10 000 year spin-up phase under pre-industrial conditions to achieve a quasi-equilibrium state. For these spin-up runs, the atmospheric CO₂ concentration was fixed at 284.7 ppm and the solar constant was set to 1360.747 W m⁻². Following the spin-up, emission-driven transient simulations were conducted under the stabilization, overshoot, and SSP5-8.5 scenarios. The results are presented as the median across all model variants, with uncertainty quantified as the range between the 5th to the 95th percentiles.

3 Results

The UVic ESCM v2.10 reliably simulates permafrost area, permafrost carbon, and total permafrost region soil carbon stocks. Under pre-industrial conditions, the simulated Northern Hemisphere permafrost area was 17.01  (17.00 to 17.04) million km², the simulated total soil carbon stock in the permafrost regions was 1031  (915 to 1149) PgC, of which 484  (383 to 590) PgC was classified as permafrost carbon and 547  (533 to 559) PgC was classified as usual soil carbon. For the period 1960–1990, the model simulated Northern Hemisphere permafrost area at 16.8  (16.7 to 16.9) million km², which falls within the reconstructed range from 12.0 to 18.2 million km² (Chadburn et al., 2017) and the observation-based extent from 12.21 to 16.98 million km² (Zhang et al., 2000). Additionally, the simulated total soil carbon stock in the permafrost regions for this same period was 1034  (919 to 1151) PgC, with 483  (382 to 587) PgC classified as permafrost carbon, accounting for 47 % (42 % to 51 %) of the total permafrost region soil carbon stock, in agreement with Hugelius et al. (2014). Notably, the area with non-zero permafrost carbon is 17.2 million km² and changes only minimally throughout the simulations.

The UVic ESCM v2.10 also realistically simulates historical temperature changes and the partitioning of anthropogenic carbon emissions among the atmosphere, ocean and land. During the period 2011–2020, the model estimated a global mean temperature increase of 1.14  (1.13 to 1.15) °C relative to preindustrial levels, which is closely aligned with the observed rise of 1.09  (0.91 to 1.23) °C (Gulev et al., 2021). From 2010 to 2019, the model estimated that anthropogenic carbon emissions of 11 Pg C yr⁻¹ were distributed as follows: 5.5  (5.4 to 5.6) Pg C yr⁻¹ to the atmosphere, 3.0  (2.98 to 3.03) Pg C yr⁻¹ to the ocean, and 2.5  (2.4 to 2.6) Pg C yr⁻¹ to terrestrial ecosystems. These estimates are broadly consistent with the global anthropogenic CO₂ budget assessment by the Global Carbon Project (GCP) with figures of 5.1 ± 0.02 Pg C yr⁻¹ for the atmosphere, 2.5 ± 0.6 Pg C yr⁻¹ for the ocean, and 3.4 ± 0.9 Pg C yr⁻¹ for terrestrial ecosystems (Friedlingstein et al., 2020).

Figure 3 Timeseries of annual mean (a) temperature change since 1850–1900, (b) permafrost area, (c) permafrost carbon loss since 1850–1900 and (d) total permafrost region soil carbon loss since 1850–1900 under stabilization (dashed lines) and overshoot (solid lines) scenarios at 1.5 °C (orange), 2.0 °C (brown), 3.0 °C (green), and 4.0 °C (purple) global warming levels, as well as the SSP5-8.5 scenario. Square markers indicate the time points when the temperature overshoot reaches its peak or stabilized warming begins, while circle markers indicate when the overshoot returns to 1.5 °C. All changes are relative to the pre-industrial period (1850–1900). Results represent the ensemble median of 250 simulations based on the PFC simulations. Dots on the adjacent right panels represent values in the year 2300, with uncertainty ranges estimated as the 5th to 95th percentiles. Notably, both permafrost carbon and total permafrost region soil carbon are calculated for the regions where the permafrost carbon is greater than zero.

Figure 4 Similar to Fig. 3, but showing timeseries of changes relative to the SWL-1.5 scenario in (a) permafrost area loss, (b) permafrost carbon loss, and (c) total permafrost region soil carbon loss under stabilization and overshoot scenarios for 2.0 °C (brown), 3.0 °C (green), and 4.0 °C (purple) global warming levels. Square markers indicate the time points when the temperature overshoot reaches its peak or stabilized warming begins, while circle markers indicate when the overshoot returns to 1.5 °C. Results represent the ensemble median of 250 simulations based on the PFC simulations. Dots on the adjacent right panels represent values in the year 2300, with uncertainty ranges estimated as the 5th to 95th percentiles.

Figure 5 Timeseries of changes in (a) permafrost carbon inputs, (b) total permafrost region soil carbon inputs, (c) permafrost carbon decomposition and (d) total permafrost region soil carbon decomposition, under the stabilization (dashed lines) and overshoot (solid lines) scenarios for 1.5 °C (orange), 2.0 °C (brown), 3.0 °C (green), and 4.0 °C (purple) global warming levels, along with the SSP5-8.5 scenario (black). Square markers indicate the time points when the temperature overshoot reaches its peak or stabilized warming begins, while circle markers indicate when the overshoot returns to 1.5 °C. Results represent the ensemble median of 250 simulations based on the PFC simulations. Dots on the adjacent right panels represent values in the year 2300, with uncertainty ranges estimated as the 5th to 95th percentiles.

3.1 Permafrost Response

The Northern Hemisphere high-latitude permafrost area is strongly correlated with changes in global mean temperature (Figs. 3a, b and S1a, b in the Supplement). As global warming increases from 1.5 to 4 °C, the permafrost area declines from 13.9 to 8.3 million km² under the SSP5-8.5 scenario by year 2300. In the stabilization scenarios, when global warming is stabilized at the Paris Agreement targets of 1.5 or 2 °C, permafrost degradation is effectively suppressed compared to the SSP5-8.5 scenario. By 2300, permafrost area decreases by 4.6  (4.5 to 4.7) and 6.6  (6.4 to 6.8) million km² from the pre-industrial level of 17.0 million km² under SWL-1.5 and SWL-2 scenarios, respectively, accounting for 39  (38 to 40) % and 56  (54 to 58) % of the reduction observed under SSP5-8.5 scenario. The additional permafrost degradation under the SWL-3 compared to SWL-2 is significantly larger than that under SWL-4 compared to SWL-3. This is mainly because the remaining permafrost available for degradation becomes progressively limited under higher global warming levels, and the permafrost area under higher stabilization scenarios has not yet reached a steady state in our simulations. Additionally, the permafrost area under the SWL-3 and SWL-4 scenarios exceeds that under the SSP5-8.5 scenario by only 2.0  (1.9 to 2.1) and 0.9  (0.8 to 1.0) million km² by 2300, respectively. This is primarily because the permafrost area in the SSP5-8.5 scenario is smaller and, as a transient simulation, is further from equilibrium compared to SWL-3 and SWL-4. However, during the cooling phase of the overshoot scenarios, as the global mean temperature returns to 1.5 °C above pre-industrial levels, the permafrost area gradually recovers to that under the SWL-1.5 scenario (Figs. 4a and S2a). By 2300, it converges to similar levels of 11.1–12.4 million km² under the OS-2, OS-3 and OS-4 scenarios, with an additional loss of only 0.2–1.2 million km² compared to the SWL-1.5 scenario. This indicates that permafrost area is nearly reversible and largely follows the global mean temperature trajectory, recovering as temperature reduction, consistent with previous studies (MacDougall, 2013; Lee et al., 2021; Schwinger et al., 2022). The incomplete recovery partly results from slightly higher regional temperatures in permafrost regions. By 2300, the northern high-latitude permafrost regions are 0.01–0.13 °C warmer compared to the SWL-1.5 scenario. It is noteworthy that permafrost area loss under overshoot scenarios typically exhibits a hysteresis effect (Boucher et al., 2012; MacDougall, 2013; Eliseev et al., 2014), peaking about 10–30 years after global mean temperature reaches its maximum in our simulations.

Under both temperature stabilization and overshoot scenarios, permafrost carbon declines monotonically over time (Figs. 3c and S1c), driven by the imbalance between weaker permafrost carbon inputs and relatively faster decomposition. The primary source of permafrost carbon input is a very slow physical process of downward diffusion through the permafrost table due to cryoturbation mixing effect (MacDougall and Knutti, 2016). As a result, the rates of permafrost carbon decomposition significantly exceed its inputs (Figs. 5a, c and S3a, c). Under the SWL-1.5, SWL-2, SWL-3, and SWL-4 scenarios, permafrost carbon losses are projected to be 54  (32 to 79), 72  (42 to 104), 106  (64 to 152), and 127  (74 to 180) PgC by 2300, respectively, while under the OS-2, OS-3, and OS-4 scenarios, the losses are 60  (35 to 87), 78  (50 to 111), and 97  (63 to 135) PgC. This indicates that despite global mean temperature returning to 1.5 °C in overshoot scenarios, considerable permafrost carbon losses still occur, accounting for 82  (80 to 85) %, 73  (71 to 78) %, and 76  (72 to 85) % of the losses in the corresponding stabilization scenarios achieving same global warming levels. Additionally, the extra permafrost carbon losses under the SWL-2, SWL-3, and SWL-4 scenarios relative to the SWL-1.5 scenario continue to increase over time; whereas in the overshoot scenarios, these additional losses decrease during the stabilization phase (Figs 4b and S2b) and reduce to 24  (4 to 52) PgC by 2300. This decrease occurs because, after temperature change stabilizes at 1.5 °C, the permafrost carbon decomposition rate under the overshoot scenarios closely aligns with that of the SWL-1.5 scenario, while the permafrost carbon inputs under the overshoot scenarios surpass that of the SWL-1.5 scenario.

Our simulations show that permafrost carbon inputs (Fig. 5a) do not follow the same trajectory as total permafrost region soil carbon inputs (Fig. 5b), especially under overshoot scenarios. Permafrost carbon inputs continue rising for some time after peak warming, whereas total permafrost region soil carbon inputs closely track the temperature trajectory (Fig. 3a). This partly stems from the unrealistic parameterization currently employed in the model. As noted in the model description (Sect. 2.1), litterfall is allocated to soil layers with temperatures above 1 °C according to an exponentially decreasing function of depth. When all soil layers are below 1 °C, organic carbon from the litterfall is added to the top soil layer. Meanwhile, permafrost carbon and usual soil carbon are both represented as depth-resolved carbon pools within the top 6 soil layers. The movement of permafrost carbon due to cryoturbation mixing is parameterized as being proportional to the gradient of total soil carbon with depth. Soil carbon that diffuses downward through the permafrost table is converted to permafrost carbon. During the cooling phase of overshoot scenarios, increased litterfall and a rising permafrost table lead to elevated carbon concentrations in surface soil layers, resulting in enhanced vertical diffusion and a surge in permafrost carbon inputs. Conversely, under the SSP5-8.5 scenario, permafrost carbon inputs exhibit only a minor peak around the 2150s, followed by a sharp decline (Figs. 5a and S3a). This is due to the continuous reduction in permafrost area and the deepening of the permafrost table, both of which reduce carbon concentrations in the upper soil layers and weaken vertical diffusion, despite the increasing litter flux under a strong CO₂ fertilization background. We note that the approach adopted in the model may not accurately describe natural processes of vertical carbon movement, which are influenced by soil porosity heterogeneity, freeze-thaw cycles, and ice expansion upon freezing.

Total permafrost region soil carbon shows a strong tendency to recover after temperature overshoot (Fig. 3d). In contrast, total permafrost region soil carbon continues to decrease under temperature stabilization scenarios, but the rate of decrease is gradually slowing down. The total permafrost region soil carbon loss in the overshoot scenarios is significantly mitigated compared to stabilization scenarios at same global warming levels (Figs. 3d and S1d). By 2300, total permafrost region soil carbon losses under OS-2, OS-3, and OS-4 scenarios are projected to be 41  (15 to 72), 61  (29 to 98), and 81  (44 to 124) PgC, respectively, with reductions of 14  (8 to 21), 41  (25 to 57), and 62  (41 to 82) PgC compared to the SWL-2, SWL-3, and SWL-4 scenarios, respectively. During the stabilization phase in the overshoot scenarios, additional total permafrost region soil carbon losses compared to the SWL-1.5 scenario continue to decrease, but the total permafrost region soil carbon in the overshoot scenarios does not fully recover to the level under the SWL-1.5 scenario by 2300 (Figs. 4c and S2c). Notably, in all stabilization and overshoot scenarios simulated in this study, the permafrost region soil serves as a net cumulative carbon source for atmospheric CO₂ by 2300. However, during the stabilization phase of OS-3 and OS-4, the permafrost region soil turns into a carbon sink, as soil carbon inputs surpass the reduced decomposition activity due to the depletion of soil carbon stocks and reduced warming levels. Total permafrost region soil carbon inputs primarily originate from processes related to vegetation litterfall, with their intensity influenced by warming levels and CO₂ fertilization effects, while total permafrost region soil carbon decomposition is closely tied to global mean temperature. In the overshoot scenarios, the peak for total permafrost region soil carbon inputs occurs slightly later than the global mean temperature peak, whereas the peak for total permafrost region soil carbon decomposition happens marginally earlier than the global mean temperature peak (Figs. 5b, d and S3b, d).

Figure 6 Timeseries of annual mean areal fraction (left column; a, d, and g), net primary productivity (middle column; b, e, and h) and vegetation carbon (right column; c, f, and i) under stabilization (dashed lines) and overshoot (solid lines) scenarios for 1.5 °C (orange), 2.0 °C (brown), 3.0 °C (green), and 4.0 °C (purple) global warming levels, as well as the SSP5-8.5 scenario (black). Each row represents one of the three dominant plant functional type (PFT): (a–c) needleleaf tree, (d–f) C₃ grass and (g–i) shrub. Square markers indicate the time points when the temperature overshoot reaches its peak or stabilized warming begins, while circle markers indicate when the overshoot returns to 1.5 °C. Results represent the ensemble median of 250 simulations based on the PFC simulations, and the shadings denote the 5th and 95th percentile uncertainty ranges. Notably, broadleaf trees and C₄ grasses are not shown because their presence in the permafrost region is negligible.

The total permafrost region soil carbon inputs generally track the trajectory of litter flux across the same area, with an approximate delay of 10–20 years (not shown). To attribute the contribution of total permafrost region soil carbon inputs, we examined how dominant vegetation types (needleleaf tree, C₃ grass and shrub) over the permafrost region adapt to temperature and atmospheric CO₂ concentrations in both overshoot and stabilization scenarios (Fig. 6). Needleleaf trees expand slowly and continuously in the permafrost region in both overshoot and stabilization scenarios, whereas the expansion of shrubs closely follows the trajectory of global mean temperature. Upon 1.5 °C warming (projected for the 2040s relative to 1850–1900), trees and shrubs are projected to cover approximately 62 % of the area where permafrost carbon is non-zero, slightly higher than the 24 %–52 % range projected for 2050 using a statistical approach that links climate conditions to vegetation types under two distinct emission trajectories (Pearson et al., 2013). During the warming and cooling phases of overshoot scenarios, the expansion and reduction of shrubs correspond with the degradation and expansion of C₃ grasses, respectively. Among the three dominant PFTs, only shrubs show a nearly reversible response in areal coverage, net primary productivity (NPP) and vegetation carbon with respect to global mean temperature under overshoot scenarios. In contrast, the continuous reduction of C₃ grasses and the expansion of needleleaf trees suggest a degree of irreversibility in the structure and vegetation carbon density of northern high latitude terrestrial ecosystems under overshoot scenarios. Our results are in line with an earlier study by Tokarska and Zickfeld (2015), but contrast with Schwinger et al. (2022) who reported only minor differences in vegetation carbon after the overshoots compared to the reference simulation with no overshoot by prescribing vegetation distributions. In our study, the shifts in vegetation composition and changes in living biomass, especially those associated with woody vegetation, are key drivers of total permafrost region soil carbon inputs.

The uncertainty in total permafrost region soil carbon loss is nearly the same as that of permafrost carbon loss (Figs. 3c, d and S1c, d). For example, the 5th to 95th percentile range of total permafrost region soil carbon loss under the OS-2 and OS-4 scenarios is 58 and 81 PgC, respectively by 2300, compared to 52 and 72 PgC for permafrost carbon loss. This indicates that the uncertainty in total permafrost region soil carbon loss is largely driven by the uncertainty in permafrost carbon loss. Therefore, we evaluate the relative importance of perturbed permafrost carbon parameters on total permafrost region soil carbon loss under different temperature pathways through calculating their correlations across all ensemble simulations. The influence of model parameters on the uncertainty of permafrost carbon losses by 2300 is relatively consistent across the SSP5-8.5, OS-4, and SWL-4 scenarios, with the strongest correlations observed for the permafrost passive carbon pool transformation rate (R = 0.81–0.85), followed by the initial quantity of total permafrost region soil carbon (R = 0.55–0.61). This finding aligns with Ji et al. (2025), who highlights the critical role of these two parameters in the uncertainty of total permafrost region soil carbon loss under temperature overshoot and 1.5 °C warming stabilization scenarios.

Figure 7 Additional (a) radiative forcing, (b) global warming and (c) permafrost area loss due to permafrost carbon feedback, calculated as the difference between the PFC and NPFC simulations. Results are shown under stabilization (dashed lines) and overshoot (solid lines) scenarios for 1.5 °C (orange), 2.0 °C (brown), 3.0 °C (green), and 4.0 °C (purple) global warming levels, along with the SSP5-8.5 scenario (black). Square markers indicate the time points when the temperature overshoot reaches its peak or stabilized warming begins, while circle markers indicate when the overshoot returns to 1.5 °C. Results represent the ensemble median of 250 simulations. Dots on the adjacent right panels represent values in the year 2300, with uncertainty ranges estimated as the 5th to 95th percentiles. In panel (a), the additional radiative forcing is calculated using the simplified expressions (Etminan et al., 2016) based on simulated CO₂ concentrations. In panel (c), the additional permafrost area loss is smoothed using a 5 year rolling average to eliminate interannual variability.

3.2 Radiative Impacts of Permafrost Carbon Loss

The permafrost carbon loss increases global mean radiative forcing and surface temperature. By comparing two parallel sets of simulations with the permafrost carbon module activated (PFC) or deactivated (NPFC), we were able to quantify the additional radiative forcing and warming caused by permafrost carbon loss. The time evolution of additional global warming closely resembles the additional radiative forcing (Figs. 7a, b and S4a, b) due to approximately linear relationship between radiative forcing and temperature change based on the energy balance of the climate system (Forster et al., 1997; Myhre et al., 2014). Both in the stabilization and overshoot scenarios, the magnitude of additional radiative forcing and warming increases with higher global warming levels. By 2300, the additional warming in the overshoot scenarios steadily rises from 0.10  (0.06 to 0.15) °C in OS-2 to 0.14  (0.09 to 0.20) °C in OS-3 and 0.18  (0.11 to 0.25) °C in OS-4. Similarly, in the stabilization scenarios, the additional warming increases from 0.10  (0.06 to 0.14) °C in SWL-1.5 to 0.13  (0.07 to 0.19) °C in SWL-2, 0.21  (0.12 to 0.31) °C in SWL-3, and 0.24  (0.15 to 0.35) °C in SWL-4. This is because higher global warming levels lead to more significant reductions in permafrost carbon and total permafrost region soil carbon, and further intensifying global warming. Furthermore, in the stabilization scenarios, the additional warming continues to increase over time due to delayed permafrost degradation and positive permafrost carbon feedback. In overshoot scenarios, the additional warming tends to stabilize once the temperature returns to 1.5 °C above pre-industrial levels. By 2300, the additional warming under SWL-2, SWL-3, and SWL-4 exceeds that of OS-2, OS-3, and OS-4 by 0.03  (0.01 to 0.04) °C, 0.07  (0.02 to 0.11) and 0.07  (0.03 to 0.10) °C, respectively, which is 22 % to 56 % higher than in the corresponding overshoot scenarios.

However, the additional warming during the cooling phase is most substantial in overshoot scenarios, and it is also significantly greater than that in stabilization scenarios over the same period. This is primarily due to the sustained reduction in atmospheric CO₂ concentration during the cooling phase, which amplifies the radiative forcing caused by permafrost carbon loss. Specifically, because of the logarithmic relationship between CO₂ concentration and radiative forcing (Etminan et al., 2016), the decline of background CO₂ concentration to low levels causes the additional increases in CO₂ concentration due to permafrost carbon loss to produce more significant changes in radiative forcing. Similarly, the global mean warming under the SSP5-8.5 scenario is 2.43  (2.41 to 2.44) °C higher than under the SWL-1.5 scenario by 2100, but the additional warming difference due to permafrost carbon loss is minimal at 0.01  (0 to 0.03) °C. Despite global mean warming reaching 8.20  (8.12 to 8.28) °C by 2300 in the SSP5-8.5 scenario, the additional warming is limited to 0.14  (0.08 to 0.20) °C due to the profoundly higher background CO₂ concentration, the additional warming is only comparable to that of SWL-2 and OS-3 scenarios.

Figure 8 Timeseries of permafrost feedback factor under stabilization (dashed lines) and overshoot (solid lines) scenarios for 1.5 °C (orange), 2.0 °C (brown), 3.0 °C (green), and 4.0 °C (purple) global warming levels, as well as the SSP5-8.5 scenario (black). Square markers indicate the time points when the temperature overshoot reaches its peak or stabilized warming begins, while circle markers indicate when the overshoot returns to 1.5 °C. Results represent the ensemble median of 250 simulations. Dots on the adjacent right panels represent values in the year 2300, with uncertainty ranges estimated as the 5th to 95th percentiles. The permafrost feedback factor is calculated as the ratio of additional global warming caused by the permafrost carbon feedback (i.e., the difference between the PFC and NPFC simulations) to the global mean temperature change in the NPFC simulations.

The additional warming caused by permafrost carbon loss can be utilized to assess whether the permafrost carbon feedback could be classified as a global tipping point process. This means it is not only positive but also sufficiently strong to sustain itself. To qualify, an initial rise in global mean temperature would need to trigger permafrost carbon emissions that result in a further increase in global mean temperature surpassing the initial warming. As a result, the positive permafrost carbon feedback would induce sufficient additional thawing to initiate a self-sustaining feedback loop (Nitzbon et al., 2024). We employed the permafrost feedback factor, which is defined as the ratio of the additional warming to the initial warming simulated with the permafrost carbon module deactivated, to determine if the permafrost carbon feedback can be considered as a global tipping process. In all perturbed parameter ensemble simulations for the stabilization, overshoot and SSP5-8.5 scenarios, the maximum permafrost feedback factor is 0.21 °C °C⁻¹. By 2300, the permafrost feedback factor for the OS4 and SSP5-8.5 scenarios are estimated at 0.12  (0.08 to 0.17) and 0.02  0.01–0.03 °C °C⁻¹, respectively. The permafrost feedback parameter is the highest under the OS4 scenario, while it is the lowest under the SSP5-8.5 scenario (Figs. 8 and S5). Interestingly, the feedback factors are quite similar across the stabilization scenarios, with values of 0.064  (0.037 to 0.096), 0.064  (0.036 to 0.095), 0.069  (0.040 to 0.103) and 0.061  (0.038 to 0.089) °C °C⁻¹ for the SWL-1.5, SWL-2, SWL-3, and SWL-4 scenarios by 2300, respectively. Although the feedback factor in the overshoot scenarios is substantially larger than the recent estimate of 0.035  (0.004–0.110) °C °C⁻¹ based on the Sixth Assessment Report of the Intergovernmental Panel on Climate Change (Nitzbon et al., 2024), our findings indicate that the positive permafrost carbon feedback is unlikely to result in enough additional thawing and corresponding carbon emissions to initiate a self-perpetuating tipping process. Since this study only models the gradual thawing of permafrost through the deepening of the active layer, we cannot rule out the possibility of tipping points associated with the abrupt thawing of talik development, thermokarst and thermo-erosion processes.

Figure 9 (a) Relationship between SPAW and global warming derived from the SSP5-8.5 scenario based on the PFC simulations. Black curve represents SPAW calculated as the slope of a local regression (Cleveland, 1979; Cleveland and Devlin, 1988; Cleveland and Grosse, 1991) between changes in permafrost area and global mean temperature under the SSP5-8.5 scenario. Colored circle markers represent transient SPAW for 1.5 °C (orange), 2.0 °C (brown), 3.0 °C (green), and 4.0 °C (purple) global warming levels. (b-i) Time series of original (solid line) and reconstructed (dashed line) additional permafrost area loss (ΔPF) due to permafrost carbon feedback for each scenario. The constructed results are derived by multiplying the transient SPAW with the additional warming (Figs. 7b and S4b) under each scenario. Smaller RMSE values indicate higher reconstruction accuracy.

The permafrost carbon feedback also causes additional permafrost area loss under stabilization, overshoot and SSP5-8.5 scenarios. During the cooling phase of OS-2, OS-3, and OS-4, there is a significant additional permafrost area loss of 0.3  (0.1 to 0.4), 0.4  (0.2 to 0.6), and 0.5  (0.3 to 0.7) million km², respectively (Figs. 7c and S4c), contributing 6  (2 to 8) %, 8  (4 to 12) % and 9  (5 to 12) % of total permafrost area loss by 2300. Although global warming has been identified as the primary driver of permafrost degradation (McGuire et al., 2016; Lawrence and Slater, 2005), the temporal evolution of additional permafrost area loss does not align with changes in the additional warming (Figs. 7b and S4b), particularly in the stabilization and SSP5-8.5 scenarios (Figs. 7c and S4c). To better understand this puzzle, we conducted sensitivity analysis based on the SSP5-8.5 scenario. The results show the sensitivity of permafrost area to global warming (SPAW, units of million km² °C⁻¹) is not constant. The maximal sensitivity occurs below a global warming level of 1.5 °C. The magnitude of SPAW decreases as global warming surpasses 1.5 °C (Fig. 9a), suggesting the response of permafrost area to rising temperature weakens and indicating a diminishing feedback effect. Interestingly, by multiplying transient SPAW by the additional warming, the temporal evolution of the additional permafrost area loss can be well reconstructed across various scenarios (Fig. 9b–i). The root mean square error (RMSE) between original and reconstructed additional permafrost area loss, ranging from 0.04 to 0.09 million km², indicates high reconstruction accuracy.

Permafrost carbon feedback usually induces greater additional warming at higher global warming levels, but the lower SPAW diminishes its impact on permafrost area loss. Under the SSP5-8.5 scenario, the additional permafrost area loss peaks around the 2060s, reaching approximately 0.19 million km², coinciding with a strong sensitivity of 3.4 million km² °C⁻¹. After this peak, despite continued rise in additional warming, the permafrost area loss declines notably due to the gradual weakening of SPAW. Similarly, in the stabilization scenario, SPAW decreases markedly at higher global warming levels, even as additional warming increases, preventing a positive correlation between the additional permafrost area loss and global warming levels (Figs. 7c and S4c). For example, the additional warming in the SWL-4 scenario reaches 0.24  (0.15 to 0.35) °C by 2300, nearly twice that of the SWL-2 scenario (0.13  (0.07 to 0.19) °C). However, due to its SPAW being only 0.8 million km² °C⁻¹, approximately one quarter of 3.4 million km² °C⁻¹ under SWL-2, the additional permafrost area loss under SWL-4 scenario is only 0.2  (0.1 to 0.3) million km², significantly lower than 0.5  (0.2 to 0.7) million km² under SWL-2. This suggests the maximal SPAW occurring near 1.5 and 2 °C global warming levels has significant implications for the Paris Agreement's targets of limiting global warming at the same levels.

Figure 10 Relationship between global warming levels and three permafrost metrics: (a) permafrost area loss, (b) permafrost carbon loss, and (c) permafrost radiative forcing in the stabilization (colored dashed lines) and overshoot (colored solid lines) scenarios for 1.5 °C (orange), 2.0 °C (brown), 3.0 °C (green), and 4.0 °C (purple) global warming levels, along with the SSP5-8.5 scenario (black). Square and circle markers indicate values in the year 2300 for the stabilization and overshoot scenarios, respectively. All results are based on the PFC simulations. Grey solid lines show linear fits of permafrost metrics to global warming levels in stabilization scenarios by 2300, while black dashed lines show corresponding fits for the SSP5-8.5 scenario. Note that in panel (a), both the stabilization scenarios and the corresponding SSP5-8.5 points included in the linear fit are limited to global warming levels between 1.5 and 3.0 °C, whereas in panels (b) and (c), the fits include points with global warming levels ranging from 1.5 to 4.0 °C. For stabilization scenarios, only the results from the year 2300 are used for fitting, while for the SSP5-8.5 scenario, all results within the specified global warming level ranges are used for fitting. Shaded regions represent the 5th to 95th percentile ranges across 250 ensemble simulations.

3.3 Linearity of Permafrost Response and Feedback

After exploring the response and feedback of permafrost under temperature stabilization and overshoot scenarios at various global warming levels, it is natural to question the linearity of these responses and feedbacks as functions of global warming levels. Our results show that the responses of permafrost area, permafrost carbon feedback and associated radiative forcing to a broad range of global warming are nearly linear (Fig. 10). The permafrost area change exhibits a strongly nonlinear relationship with global warming below 1.5 °C level, then a quasilinear relation between them in the global warming ranges from 1.5 to 3 °C. Above 3 °C global warming, the sensitivity of permafrost area to global warming decreases nonlinearly, and it is evident in both stabilization and SSP5-8.5 scenarios (Fig. 10a). In contrast, permafrost carbon loss and associated radiative forcing exhibit a nearly linear response to increasing global warming levels, especially above 1 °C, for both stabilization and SSP5-8.5 scenarios (Fig. 10b and c).

Meanwhile, the sensitivities of permafrost area, permafrost carbon loss, and associated radiative forcing to global warming under stabilization scenarios are all stronger than those under the SSP5-8.5 scenario. Specifically, based on the simulated permafrost area in the year 2300 under stabilization scenarios with global warming levels between 1.5 and 3 °C, the sensitivity of permafrost area to global warming is −3.19  (−3.01 to −3.36) million km² °C⁻¹. In comparison, a linear fit of permafrost area change against global warming levels over the same temperature range in the SSP5-8.5 scenario yields a sensitivity of −2.85  (−2.77 to −2.89) million km² °C⁻¹. Similarly, the permafrost carbon feedback per degree of global warming derived from a linear fit based on the total permafrost carbon loss in the year 2300 under stabilization scenarios, is −27.6  (−16.5 to −38.2) Pg C°C⁻¹. In contrast, the corresponding value under the SSP5-8.5 scenario, estimated from a linear fit over the 1.5 to 4.0 °C warming range, is −19.3  (−15.7 to −24.1) Pg C°C⁻¹. Applying the same approach, the associated radiative forcing per degree of global warming is estimated to be 0.08  (0.05 to 0.12) $\mathrm{W}{\mathrm{m}}^{-\mathrm{2}}\mathrm{°}{\mathrm{C}}^{-\mathrm{1}}$ for the stabilization scenarios and 0.04  (0.03 to 0.05) $\mathrm{W}{\mathrm{m}}^{-\mathrm{2}}\mathrm{°}{\mathrm{C}}^{-\mathrm{1}}$ for the SSP5-8.5. These differences between the stabilization and SSP5-8.5 scenarios are mainly attributable to the differing response time scales represented by the two scenarios: SSP5-8.5 reflects a typical transient response, while the stabilization scenarios maintain stabilized temperatures over extended periods and thus approximate a quasi-equilibrium response of the climate-carbon system. Furthermore, the smaller sensitivity of permafrost radiative forcing per degree of global warming under the SSP5-8.5 can be partially attributed to its higher background atmospheric CO₂ concentration compared to the stabilization scenarios. The same amount of CO₂ emissions would produce smaller additional radiative forcing under a higher background atmospheric CO₂ concentration, due to the logarithmic relationship between CO₂ concentration and radiative forcing (Etminan et al., 2016).

To a certain extent, our findings align with those of Nitzbon et al. (2024), who suggested that the accumulated response of Arctic permafrost to climate warming is approximately quasilinear. Nitzbon et al. (2024) reported a quasilinear decrease in the equilibrium permafrost extent to global warming, with a rate of approximately 3.5 million km² °C⁻¹. This quasilinear relation holds for global warming ranges from 0 to 4 °C, derived from the empirical relationship between the local permafrost fraction and the annual mean global temperature. However, when applying a local regression method across the full SSP5-8.5 trajectory (Fig. 9a), we find that the quasilinear relationship only holds between 1.5 and 2 °C. In contrast, when applying a simple linear regression, a quasilinear behavior still holds over the 1.5 to 3 °C warming range (Fig. 10a). Furthermore, the permafrost carbon feedback and the associated radiative forcing per degree of warming, as derived from our simulations of both stabilization and SSP5-8.5 scenarios, are within the ranges of −18  (−3.1 to −41) Pg C°C⁻¹ and 0.09  (0.02 to 0.20) $\mathrm{W}{\mathrm{m}}^{-\mathrm{2}}\mathrm{°}{\mathrm{C}}^{-\mathrm{1}}$, respectively, reported by Canadell et al. (2021). Our estimates also align with the estimated range of equilibrium sensitivity of permafrost carbon decline to global warming, which is −21  (−4 to −48) Pg C°C⁻¹ (Nitzbon et al., 2024). This may represent an upper bound on the permafrost carbon feedback per degree of global warming, as the estimate reflects the maximum amount of carbon that could be lost from frozen soils.

Under overshoot scenarios, permafrost area responds nearly reversibly and presents an almost closed loop (Fig. 10a). In contrast, permafrost carbon loss exhibits an open loop with respect to global warming levels. In other words, permafrost carbon loss does not reverse as temperatures decline, indicating irreversible permafrost carbon radiative forcing. Among the three metrics investigated here, only permafrost area exhibits strong reversibility under the overshoot scenarios. This also explains why, in Fig. 9a, the permafrost area sensitivity derived from the SSP5-8.5 scenario, when multiplied by additional warming, can reasonably reconstruct permafrost area loss in the stabilization and overshoot cases.

4 Conclusions and discussion

This study utilizes the intermediate-complexity Earth system model UVic ESCM v2.10 and perturbed parameter ensemble modelling approach to study the response and feedback of permafrost under temperature stabilization and overshoot scenarios across different global warming levels. The UVic ESCM v2.10 has been validated against observational and reconstructed datasets, demonstrating its ability to reproduce observation-based permafrost area and permafrost carbon stocks for the present day. In addition to presenting the changes in permafrost area, permafrost carbon, and total permafrost region soil carbon under various warming trajectories, this study also quantifies the additional radiative forcing, global warming, and permafrost area loss induced by permafrost carbon loss and suggests that permafrost carbon feedback is unlikely to initiate a self-perpetuating global tipping process under both stabilization and overshoot scenarios. In addition, this study reveals quasilinear responses of permafrost area, permafrost carbon, and associated radiative forcing over a broad range of global temperature change, providing insights into the implications of permafrost carbon feedback for long-term climate change and mitigation strategies.

Reductions in permafrost area and carbon exhibit a strong correlation with global warming for both stabilization and overshoot scenarios. In stabilization scenarios, lower global warming levels effectively mitigate permafrost degradation compared to the SSP5-8.5 scenario, whereas higher global warming levels lead to substantial permafrost degradation due to cumulative warming effects. In overshoot scenarios, permafrost area largely recovers as global warming returns to 1.5 °C levels, though this recovery is delayed by hysteresis effects, with degradation persisting for decades after temperature peaks. Permafrost carbon declines under both stabilization and overshoot scenarios, driven by the dynamic balance between soil carbon inputs and decomposition. The overshoot scenarios partially mitigate permafrost carbon losses compared to the stabilization scenarios at the same global warming levels. Significant carbon losses persist during the cooling and stabilization phases of overshoot scenarios, highlighting the essentially irreversible nature of this process. Total permafrost region soil carbon exhibits a certain degree of recovery under the overshoot scenarios. In fact, soil in these regions even transitions into a carbon sink during the stabilization phase of the overshoot scenarios with high global warming levels, supported by reduced decomposition rates and sustained inputs from vegetation litterfall. However, the higher the overshoot levels, the less recovery there is. These findings underscore the critical role of temporary temperature overshooting levels in affecting long-term permafrost carbon loss and recovery potential.

The responses of permafrost area, permafrost carbon, and associated radiative forcing to a broad range of global warming are nearly linear. The permafrost area and global mean temperature exhibit a quasilinear relation for the global warming ranges from 1.5 to 3 °C in both the stabilization and SSP5-8.5 scenarios. The permafrost carbon loss and associated radiative forcing exhibit a quasilinear relation to global warming ranges from 1 to 4 °C under the stabilization and SSP5-8.5 scenarios. The sensitivity of permafrost area to global warming derived from the stabilization and SSP5-8.5 scenarios is substantially higher than the estimates obtained using the Stefan solution (Peng et al., 2023), but our result is close to that derived from an observation-constrained equilibrium projection, approximately 3.5 million km² °C⁻¹ (Nitzbon et al., 2024) (See Fig. 9a). Nitzbon et al. (2024) noted that permafrost area decreases quasi-linearly with increasing global mean temperature. However, we found that this relationship holds most strongly within 1.5 to 3 °C global warming range when using a linear regression approach (Fig. 10a). According to the SSP5-8.5 simulations, the sensitivity of permafrost area to global warming reaches its peak below 1.5 °C global warming. It then decreases to a relatively stable level between 1.5 and 3 °C warming levels, and continues to decline beyond the 3 °C warming level. This is in line with Comyn-Platt et al. (2018), who found the feedback processes due to permafrost thaw respond more quickly at temperatures below 1.5 °C. The maximal sensitivity occurring just below 1.5 °C global warming level suggests the fastest permafrost degradation is anticipated to take place within Paris Agreement's warming levels. Our findings have significant implications for the development of mitigation and adaptation strategies addressing permafrost-thaw impacts consistent with keeping global warming at the Paris Agreement's levels.

Our study highlights the substantial permafrost carbon feedback during the cooling phase of overshoot scenarios. Permafrost carbon loss evidently increases global radiative forcing and amplifies global warming. The permafrost carbon feedback can be more profound in temperature stabilization and overshoot scenarios than in high-emissions scenarios. In stabilization scenarios, additional radiative forcing and warming from permafrost carbon loss continue to rise after temperature stabilizes, driven by delayed permafrost thaw and the resulting positive permafrost carbon feedback. In overshoot scenarios, additional warming almost stabilizes once global warming drops to 1.5 °C levels. During the cooling phase of overshoot scenarios, lower background CO₂ concentrations amplify the warming effect of permafrost carbon loss. In contrast, under the high-emissions SSP5-8.5 scenario, the additional warming is limited due to higher background CO₂ levels reducing the additional radiative forcing from permafrost carbon loss. The additional permafrost area loss due to permafrost carbon feedback is around 5  (2 to 11) % of the total loss by 2300 in the stabilization and overshoot scenarios. This additional permafrost area loss can be well explained by the sensitivity of permafrost area to global warming and the magnitude of additional warming. The complex interactions between global warming and permafrost degradation emphasize the importance of accounting for these nonlinear effects in climate projections, particularly at 1.5–2 °C global warming levels.

Our results show incomplete recovery of permafrost area under the overshoot scenarios, which is influenced by multiple factors: First, the additional permafrost carbon loss leads to greater additional warming under the overshoot scenarios than the SWL-1.5 scenario, causing additional permafrost degradation. Second, the thermal inertia of deep soil layers limits the rate of permafrost recovery. Even after global mean temperatures return to the 1.5 °C target, residual heat accumulated in deeper soil layers during temperature overshoot period continues to inhibit permafrost refreezing, preventing full restoration to its pre-overshoot state. Third, substantially greater soil carbon loss under overshoot scenarios than in the SWL-1.5 scenario alters the hydrological and thermal properties of soil, affecting the processes that govern carbon cycling (Zhu et al., 2019; Avis, 2012; Lawrence and Slater, 2008), which in turn affects the recovery of permafrost area. Moreover, irreversible shifts in vegetation composition of high-latitude terrestrial ecosystems also contribute to the incomplete recovery of permafrost area under overshoot scenarios. For instance, among the two dominant vegetation types, needleleaf trees continue to expand while C₃ grasses decline, even after global temperatures return to the 1.5 °C warming level. These irreversible changes may stabilize the carbon, water, and energy cycles over the permafrost region at different equilibria after overshoot, through the interactions between physical and biophysical processes (de Vrese and Brovkin, 2021), thereby constraining the ability of permafrost to fully recover under the overshoot scenarios.

Although permafrost carbon loss is essentially irreversible, overshoot scenarios exhibit a certain degree of recovery relative to the SWL-1.5 stabilization scenario (Figs. 4b and S2b). It is therefore scientifically interesting to know whether permafrost carbon under overshoot scenarios will eventually converge with that under SWL-1.5. Our results show that permafrost carbon inputs are consistently higher under overshoot scenarios than under SWL-1.5, while permafrost carbon decomposition differ only slightly between the two (Figs. 5a, c and S3a, c). This tends to suggest that the smaller permafrost carbon stocks under overshoot scenarios by 2300 would eventually catch up to the levels under SWL-1.5. To assess this potential convergence, we extended our simulations of both SWL-1.5 and overshoot scenarios to the year 2400 (data not shown). Then we estimated the convergence time by calculating the ratio between the difference in permafrost carbon stocks and the difference in net permafrost carbon inputs (i.e., annual permafrost carbon inputs minus decomposition) for the overshoot scenarios relative to the SWL-1.5 scenario. Based on simulation results for the year 2300, the median estimated convergence times for OS-2, OS-3 and OS-4 are 1076, 1008 and 1433 years, respectively. When using results from the year 2400, the corresponding estimates increase to 1377, 1199 and 1568 years. This means that convergence would take even longer if estimated from later simulation results, mainly due to gradually weakened permafrost carbon inputs. The relatively larger permafrost carbon inputs under overshoot scenarios result mainly from increased litterfall during the overshoot phase. The extra litterfall during the overshoot phase gradually moves through the active layer and is transported to the permafrost zone. Over time, however, the effect of this extra litterfall gradually diminishes. Consequently, it may take extremely long timescales for the overshoot scenarios to fully converge with SWL-1.5 in terms of permafrost carbon stocks. In addition, due to incomplete recovery of permafrost area and persistent changes in surface climate and soil properties, the overshoot scenarios might ultimately fail to converge to SWL-1.5 scenario in terms of permafrost carbon stocks.

Different permafrost carbon loss and associated additional warming under overshoot scenarios confirm the path-dependent fate of permafrost region soil carbon (Kleinen and Brovkin, 2018) and the path-dependent reductions in CO₂ emission budgets (MacDougall et al., 2015; Gasser et al., 2018). Since the permafrost carbon has accumulated very slowly during the last few tens of thousands of years, its loss would contribute to atmospheric CO₂ burdens, and the resulting additional warming poses a challenge to achieving global climate goals by substantially reducing the remaining carbon budget compatible with the Paris Agreement (MacDougall et al., 2015; Natali et al., 2021). In the overshoot scenarios simulated in this study, permafrost carbon loss by 2300 ranges from 60  (35 to 87) to 97  (63 to 135) PgC. The resulting additional warming from the permafrost carbon loss ranges from 0.10  (0.06 to 0.15) to 0.18  (0.11 to 0.25) °C. This permafrost carbon feedback contributes a substantial addition on top of 1.5 °C warming target under overshoot scenarios, and the magnitude of this additional warming rises with the amplitude of overshoot. To accomplish the 1.5 °C target under the OS-2, OS-3, and OS-4 scenarios, anthropogenic carbon emissions would need to be reduced by amounts equivalent to the permafrost carbon loss. The proportion of carbon removal required to offset permafrost emissions corresponds to 27  (16 to 39) % under OS-2, 10  (6 to 15) % under OS-3, and 7  (4 to 9) % under OS-4, relative to the cumulative carbon removal from the peak warming to the 1.5 °C target sustained until 2300. Our findings are consistent with previous research utilizing the Monte Carlo ensemble method to evaluate the response of permafrost carbon and its influence on CO₂ emission budgets under overshoot scenarios targeting a 1.5 °C warming limit (Gasser et al., 2018). Specifically, for overshoot amplitudes of 0.5 °C (peak warming of 2 and 1 °C (peak warming of 2.5 °C), the reductions in anthropogenic CO₂ emissions due to permafrost are estimated to be 130 (with a range of 30–300) Pg CO₂ and 210 (with a range of 50–430) Pg CO₂, respectively, to meet the long-term 1.5 °C target (Gasser et al., 2018). These results are comparable to our estimates of 60  (35 to 87) PgC under OS-2 and 78  (50 to 111) PgC under OS-3. The differences between the two studies can be partly attributed to different warming trajectories to achieve the same 1.5 °C target. Our study further confirms that deploying negative CO₂ emissions to reverse anthropogenic warming would be undermined by substantial post-peak warming arising from permafrost carbon loss, thereby diminishing their overall effectiveness (MacDougall, 2013; Tokarska and Zickfeld, 2015).

This study does not simulate the changes of deep Yedoma carbon under the temperature stabilization and overshoot scenarios. Yedoma deposits represent a significant deep carbon reservoir and are widespread across Siberia, Alaska, and the Yukon region of Canada, having primarily formed during the late Pleistocene, especially in the late glacial period. These deep, perennially frozen sediments are particularly ice-rich, and the freeze-locked organic matter in such deposits can be re-mobilized on short time-scales, representing one of the most vulnerable permafrost carbon pools under future warming scenarios (Schuur et al., 2015; Strauss et al., 2017). According to Zimov et al. (2006), these perennially frozen Yedoma sediments cover more than 1 million km², with an average depth of approximately 25 m. Recent estimates place the organic carbon stock in Yedoma deposits at 213 ± 24 PgC, constituting a significant portion of the total permafrost carbon stocks (Strauss et al., 2017). However, the UVic ESCM v2.10 utilized in this study simulates permafrost carbon only within the top 3.35 m of soil, limiting our ability to directly assess the impacts of temperature overshoot on deep Yedoma carbon. Considering their ice-rich nature and potential susceptibility to rapid-thaw processes, we analyzed the average and maximum active layer thickness (ALT) in Yedoma regions (Strauss et al., 2021, 2022) under the simulated scenarios to approximate potential impacts. We find that the average ALT in Yedoma regions remains below 1 m in all stabilization and overshoot scenarios, while the maximum ALT rarely exceeds 3.35 m in overshoot scenarios but does exceed this depth in some stabilization scenarios. However, in all scenarios, the maximum ALT does not exceed 6 m, which is relatively shallow compared to the average depth (∼ 25 m) of Yedoma deposits (Fig. S6). Consequently, the impact on Yedoma is expected to be minimal in all scenarios, and the effect of overshoot scenarios on the deep Yedoma carbon is relatively minor compared to stabilization scenarios as well.

This study, like previous ones, uncovers considerable uncertainty in projections of permafrost carbon under global warming. The uncertainty represented by perturbed model parameters for each scenario can be interpreted as model uncertainty. We note that model uncertainty in permafrost carbon loss gradually increases with the peak warming level and the duration of overshoot for each scenario (Figs. 3c and S1c). However, the uncertainty ranges in permafrost carbon loss substantially overlap for overshoot and stabilization scenarios with adjacent warming levels, such as OS-2, SWL-1.5 and SWL-2. This is especially evident in low-level warming scenarios, where the uncertainty in projected permafrost carbon loss is mainly driven by model uncertainty due to parameter perturbations, rather than scenario-related uncertainty. Given the significant roles of the permafrost passive carbon pool transformation rate and the initial quantity of total permafrost region soil carbon in determining the uncertainty of total permafrost region soil carbon loss, it is expected that these two parameters contribute significantly to the overlapping uncertainty ranges of permafrost carbon and total permafrost region soil carbon losses across different warming levels. Due to the interaction with soil carbon inputs, the overlapping uncertainty in total permafrost region soil carbon loss tends to differ from that of permafrost carbon loss. For example, the uncertainty ranges in permafrost carbon loss under OS-4 and SSP5-8.5 scenarios show considerable overlap, but the same does not apply to total permafrost region soil carbon loss, which results from significant differences in soil carbon inputs under distinct CO₂ fertilization backgrounds. The large overlapping uncertainty in projecting permafrost carbon loss under low-level warming scenarios, as shown in this study and in previous research (MacDougall et al., 2015; MacDougall and Knutti, 2016; Gasser et al., 2018), constitutes a significant challenge in accurately estimating the remaining carbon budgets consistent with temperature goals of the Paris Agreement.

This study aims to serve as a meaningful supplement to the limited existing research on permafrost response and feedback under temperature stabilization and overshoot scenarios. For example, CMIP6 includes only one overshoot scenario (SSP5-3.4-OS), which is insufficient for comprehensively analyzing the potential contribution of permafrost feedback under different global warming levels (Melnikova et al., 2021). The upcoming CMIP7 is planned to encompass a broader range of overshoot scenarios, from very low to medium and high overshoot scenarios. This expansion aims to delve into the potential for climate restoration, the feasibility of achieving Paris Agreement targets, and the risks of irreversibility and hysteresis in the slower components of the Earth system from beyond 2125 to 2300 (van Vuuren et al., 2025; Sanderson et al., 2024). Our study investigated the permafrost response and feedback under overshoot scenarios with varying levels of warming, which can be used for comparative analysis of the upcoming CMIP7 multi-model results. Additionally, research on permafrost changes under different climate stabilization scenarios is also limited. Some of existing studies are based on equilibrium permafrost models (e.g. Liu et al., 2021), which do not consider transient effects and cannot predict the timing of permafrost loss based on imposed warming scenarios (Smith et al., 2022). Our study utilizes a process-based model and provides insights into the relationship between permafrost degradation and global warming.

This study demonstrates a method to quantify the permafrost carbon feedback under specified global warming levels and warming trajectories with a fully coupled Earth system model driven by specified CO₂ emission pathways. As permafrost carbon feedback is highly sensitive to permafrost states and processes, controlling global warming level or warming trajectory is beneficial to isolate individual contributions from physical climate and carbon cycle processes. Due to tight link between terrestrial carbon feedback and physical climate feedback, Goodwin (2019) found around half the uncertainty in derived terrestrial carbon feedback originates from uncertainty in the physical climate feedback. Therefore, controlling global warming level or warming trajectory is helpful to accurately assess a specific aspect of terrestrial carbon feedback by isolating the uncertainty from the physical climate feedback. We encourage the exploration of potential application for this method in quantifying other aspects of carbon-climate feedbacks. The findings of this study are based on a single Earth system model, and it motivates further studies aimed at comprehensively understanding of permafrost response and feedback.