A total of 2,212 modern pollen assemblages published previously (Figure 1) for northern continental Asia (east of 50°E and north of 45°N) are included in our analyses: most of them are extracted from topsoil, moss samples, or lake core-top samples. The Chinese and Mongolian modern pollen data primarily come from Cao et al. (2014), while Siberian data primarily come from Bordon et al. (2009) and Natalia et al. (2020). Other data include records for the northern and central Yakutia (Müller et al., 2010), the Russian Far-East (Tarasov et al., 2011), and the Khatanga River region (Niemeyer et al., 2015; Klemm et al., 2016). The fossil pollen records (Figure 1) were obtained from lacustrine sediments and peat from the same area of northern Asia, and comprise 6,873 fossil pollen assemblages from 199 records (Cao et al., 2020). For these modern and fossil pollen data, pollen names were taxonomically harmonized at the family or genus level generally and pollen percentages were re-calculated based on the total number of terrestrial pollen grains, following the method described by Cao et al. (2013, 2020). Based on calibrated radiocarbon dates, an age-depth model was established for each record using a Bayesian approach (further details are described in Cao et al., 2013). In this study, we selected records between 12.5 and 0 cal ka BP, covering the Holocene.

Modern permafrost data was downloaded from the open database Circum-Arctic Map of Permafrost and Ground-Ice Conditions (version 2) (Brown et al., 2002), which maps the distribution of permafrost in the Northern Hemisphere from 20°N to 90°N. Based on the estimated percentage of area, the permafrost extent is divided into continuous permafrost (90–100%), discontinuous permafrost (50–90%), sporadic permafrost (10–50%), isolated patches (<10%), and no permafrost. For easier classification, we merged the three transitional states between continuous permafrost and no permafrost into one category of non-continuous permafrost. The resulting three states of continuous permafrost, non-continuous permafrost, and permafrost-free are used throughout our research. Using this distribution map of permafrost, conditions prevailing at the modern pollen sites are assigned.

We fitted a logistic regression model to estimate the probability of the presence of permafrost in northern Asia by relating modern pollen percentages from 74 pollen taxa shared by modern and fossil pollen records to the permafrost data. To obtain more reliable regression results, we only selected samples from regions with permafrost-free and continuous permafrost. The association between the pollen percentage (predictor variable) and the presence of permafrost (response variables) was assessed by the Wald test. A p-value less than 0.05 is considered statistically significant. The logistic regression was implemented using the glm function in the stats package in R version 4.0.3 (R Core Team, 2020).

Random Forest is a machine-learning algorithm that can be used to solve regression and classification problems (Breiman, 2001). RF operates by constructing a multitude of decision trees. Although factors such as the number of trees in RFs may bias the regression results to some extent (Strobl et al., 2007; Arlot and Genuer, 2014), RFs have been successfully used in the field of earth science to predict future species distributions, and to reconstruct local and even global past tree species distributions (Prasad et al., 2006; Benito Garzón et al., 2007; Lindgren et al., 2021; Qin, 2021).

Our RF model was trained on 74 pollen taxa from the 2,212 modern samples that were matched to permafrost conditions. For the 2,212 samples, 70% of samples served as the training set for RF, while the remaining 30% of samples served as a test-set. We conducted three separate training runs using RF to assess the stability of the model and selected the training model with the highest overall statistical accuracy and Kappa value. The trained RF was then applied to down-core palynological records from 199 boreholes, sampled in northern Asia to reconstruct past permafrost conditions. These steps were implemented collectively using the R version 4.0.3 built-in package randomForest (version 4.6-14; Liaw, 2018).

Pollen data from topsoil together with current permafrost extent offer a unique opportunity to understand the relationship between the presence of permafrost and the occurrence of specific pollen taxa. Among the 74 pollen taxa analyzed, 20 taxa significantly correlate with the presence of permafrost (P < 0.05). Alnus (shrub), Betula (shrub), Larix, and another 13 taxa are positively correlated with the presence of permafrost, while Abies, Pinus (Diploxylon), Picea, and another 7 taxa are negatively correlated with the presence of permafrost (Supplementary Figure 1). Using Larix and Abies as examples, the probability of permafrost being present increases as the percentage of Larix pollen increases but decreases as the percentage of Abies pollen decreases (Figure 2). This suggests that the variability of pollen taxa can reflect shifts between the presence or absence of permafrost and thus be used to reconstruct permafrost conditions.

Based on assessments of the test-sets of three separate training runs, the RF model has great stability with a mean accuracy of 0.84 and a mean Kappa value of 0.76 (Table 1) and almost no variance between the runs. Mean Decrease Accuracy and Mean Decrease Gini of a taxon indicate the importance of that taxon to the accuracy of permafrost state classification. Of the 74 taxa used in our classification model, 30 had a Mean Decrease Accuracy value of more than 10 and 29 had a Mean Decrease Gini value of more than 10 (Figure 3). Obviously, these taxa play a key role in classifying the state of permafrost.

The precision and recall rate are two important statistical metrics to evaluate the classification quality, measuring the fraction of relevant instances among the retrieved instances and the fraction of relevant instances that were retrieved, respectively. The precision and recall rate of our RF classification for different states of permafrost differ slightly (Table 2) being best for continuous permafrost and worst for non-continuous permafrost.

To gain an understanding of the overall condition in different time periods, we take the most frequently reconstructed state of each record in each 1000-year window. For different time slices of the Holocene, the distribution of permafrost changed significantly, especially in the early-to-middle Holocene (Figure 4 and Supplementary Figure 2). From 12 to 5 cal ka BP, in addition to Asia north of 60°N, which is currently the main distribution area of continuous permafrost, the vast northern Asia region south of 60°N, such as the West Siberian Plain, southeastern Siberia, Lake Baikal region, and Kamchatka peninsula, sporadically experienced continuous permafrost state. Focusing on the West Siberian Plain, despite the widespread permafrost conditions before 9 cal ka BP, no continuous permafrost was reconstructed for this region at 7 cal ka BP. The main distribution of non-continuous permafrost between 12 and 5 cal ka BP is found in northern Asia south of 60°N, especially in the Tianshan-Altai region and the Mongolian Plateau region. However, there is no obvious spatial pattern in the distribution of sites with non-continuous permafrost. In contrast to the area north of 60°N, a permafrost-free state is frequently reconstructed from sites in the area south of 60°N. Additionally, on the West Siberian Plain between 9 and 7 cal ka BP, permafrost-free states increased significantly.

The percentages of the three states at different time slices portray the evolution of permafrost in northern Asia during the Holocene (Figure 5). Continuous permafrost is always the highest among the three states throughout the Holocene, ranging from 41.6 to 67.1%. Continuous permafrost clearly decreases in the early Holocene (12–8 cal ka BP), with the rate of decline slowing down after reaching 48.8% at 8 cal ka BP. A slight increase in this state is observed in the most recent 1 ka. The proportion of non-continuous permafrost ranges from 19.0 to 24.7%, with an average value of 21.5%. Although the proportion of non-continuous permafrost has a stable trend during the Holocene, there was a slight increase at 1 cal ka BP followed by a decrease. In corollary, the percentage of permafrost-free state shows a sharply increasing trend from 12 to 8 cal ka BP and a gradually increasing trend after 8 cal ka BP until the late Holocene when it decreased. The decrease during 2–1 cal ka BP is caused by the formation of more non-continuous permafrost, while the decrease during 1–0 cal ka BP is caused by an increase in continuous permafrost.

Reconstructing past permafrost conditions using pollen data is credible as confirmed by logistic regression and statistical metrics from test-set. The relationship between pollen taxa and permafrost has been identified by logistic regression, and the results presented here imply that we can reconstruct the history of permafrost using pollen data. Additionally, the high average accuracy (0.84) based on the test-set means that the state of permafrost down-core can be reconstructed with a reasonably low error, while the Kappa statistic (up to 0.76) also suggests a substantial classification quality of the RF model (Landis and Koch, 1977). The precision and recall rate of continuous permafrost reach 0.90 and 0.88, respectively, indicating a good ability of the model. Although some of the samples in the test-set were mismatched by our trained RF model, it is worth noting that these misallocated samples are mainly distributed in marginal regions where the permafrost is in a state of flux (Supplementary Figure 3) and the pollen taxa show transitional characteristics. Overall, the classification results of samples within the different permafrost states are very reliable.

The synchrony of changes in permafrost conditions with the Northern Hemisphere climatic signal and consistency with other permafrost proxies in northern Asia further support the reliability of pollen-based permafrost reconstructions (Figure 6). The freezing and thawing of permafrost have a close relationship to climate, and our reconstruction demonstrates the variation in permafrost state in northern Asia could depend on suborbital variations in summer solar insolation (Laskar et al., 2004). During the early Holocene, permafrost thawed rapidly (Figure 6C, from 12 to 8 cal ka BP), coinciding with a sharp increase in Northern Hemisphere temperature and North Atlantic air temperature after the Younger Dryas (Figures 6A,B; Andersen et al., 2004; Shakun et al., 2012; Marcott et al., 2013). After the northern hemisphere temperature reached a maximum in the middle Holocene, the rate of permafrost thawing in northern Asia slowed as the temperature and summer solar insolation decreased. The decreasing trend in proportion of permafrost-free during the late Holocene (Figure 6D) can also be explained by variation in the Northern Hemisphere temperature anomalies (Figure 6B). Continuous permafrost expanded during 2–0 cal ka BP when the stacked temperature anomalies show a return to under 0°C. The temporal and spatial variation of permafrost in western Siberia is very significant, which is consistent with the traditional research results of permafrost. Research on the palaeo-permafrost and palaeo-periglacial evidence shows that the severe freezing conditions that developed during Last Glacial time persisted of western Siberia in early Holocene, but until the Holocene optimum (7.5–4 cal ka BP), the southern boundary of surficial permafrost receded northward to near the Arctic Circle (Supplementary Figure 4; Velichko et al., 1984). Additionally, thermokarst lake formation history (Figure 6H), which can indicate permafrost degradation on a large scale, has recently been published for northern mid- to high-latitudes (Brosius et al., 2021) and the result of northern Asia suggests a rapid permafrost thawing period during the early Holocene (12–8 cal ka BP) and a relatively slow decreasing rate of thaw in the mid-late Holocene (7–2 cal ka BP): this result is consistent with our reconstruction based on pollen data.

A reanalysis of vegetation turnover (Figure 7), including evergreen conifer tree, summer-green conifer tree, summer-green broad-leaved tree, and non-tree vegetation from northern Asia based on pollen records (Cao et al., 2019) is instructive for our understanding of the relationship between permafrost and vegetation. Although an early-to-middle Holocene (12–7 cal ka BP) decreasing trend is observed in the non-tree vegetation component, the other vegetation components do not show a corresponding increasing trend. It was not until 8 ka that evergreen conifer species, dominated by pine and spruce, increased significantly, which is consistent with the increase of permafrost-free conditions. A re-analysis of the spatiotemporal distribution of the four vegetation types (Figure 8) shows that summer-green conifers (larch) were the main trees in northern Asia in the early Holocene, and the transition from summer-green conifer trees to evergreen conifer trees occurred around 8 cal ka BP, especially on the Western Siberian Plain (G2, G3, G8, G9), which coincided with the establishment of permafrost-free conditions.

It is reasonable to assume that the presence of permafrost in the early Holocene inhibited the spread of dark coniferous taiga. Despite increases in temperature, permafrost could not thaw deeply enough across Siberia to support dark coniferous taiga (dominated by pine) during the early Holocene. Traditional proxy-based climate reconstructions in the Northern Hemisphere suggest abrupt warming occurred at the onset of the Holocene, followed by a long-term cooling trend throughout the middle-to-late Holocene (Shakun et al., 2012; Marcott et al., 2013), which should be consistent with afforestation. However, the turnover of vegetation is not simply driven by changes in climate as traditionally thought, such as seed dispersal strategy of specific species (Travis et al., 2013), competition to the existing dominant vegetation (Meier et al., 2012), the existence of permafrost (Tian et al., 2018; Cao et al., 2019), and extensive fires (Schulze et al., 2012), which may affect the progress of plant migration. Research in northeast Asia has challenged the general view that the vegetation-climate lag is last no more than a few centuries, suggesting that during the Plio-Pleistocene, interglacial vegetation was influenced by the persistence of permafrost, mainly reflects conditions of the previous glacial period, lagging the climate by several millennia (Herzschuh et al., 2016). Therefore, compared with other non-permafrost regions, the extensive and deep permafrost in our study area, northern Asia, may have a profound impact on vegetation. Previous studies have demonstrated that larch can survive on permafrost with an active layer depth of less than 40 cm (Osawa et al., 2010), while pine requires at least 1.5 m of active layer to grow (Tzedakis and Bennett, 1995). In Alaska, permafrost with a shallow active layer is assumed to limit the northward extent of Picea (Lloyd, 2005). A pollen-based biochemical study also finds that during the early Holocene, only larch forests could survive on these shallow active-layer permafrost regions in northern Asia (Tian et al., 2018). A large amount of heat is needed to thaw permafrost and it takes time for deep permafrost to thaw (Galushkin, 1997; Vaks et al., 2013). However, our studies have not taken all relevant geoecological causal chains and their interactions into account. Our knowledge of permafrost and vegetation delay is thus still limited.

The freezing and thawing of permafrost can change the hydrological processes, and the turnover in vegetation controlled by permafrost could affect the land-surface albedo. The thawing of frozen soil will lead to changes in the ice/water ratio, altering the thermal and hydraulic qualities of the soil, and also affecting ground surface temperature (Xin et al., 2012). Numerical simulations with a supercooled soil-water model show that the thawing of frozen soil has a great impact on winter and spring soil moisture as well as temperature in Eurasia (Li et al., 2011; Xin et al., 2012). Landscapes with trees generally have a lower albedo than landscapes without trees and this difference may be magnified when there is snow cover. The albedo values of tundra and grassland when covered by fresh snow, can reach 0.76 and 0.65, respectively, while the average albedo value of three other vegetation types with trees is only 0.27 (Schaeffer et al., 2006). In addition, larch, being a deciduous coniferous species, lacks a canopy in winter and has a poor ability to cover snow compared to evergreens (Betts and Ball, 1997), with a peak winter albedo value of up to twice that of evergreen coniferous forest (dark taiga forest; Shuman et al., 2011).

A series of global atmospheric circulation models has suggested that the Eurasian snow cover has a significant effect on the Asian summer monsoon (Yang and Xu, 1994; Liu and Yanai, 2002; Wu et al., 2009; Xu et al., 2021), with a recent study showing that when there was excessive spring snowmelt in Siberia between 1981 and 2014, the EASM was weaker, resulting in less summer precipitation in north China (Xu et al., 2021). The inter-decadal variation of snow cover and summer monsoon can essentially be attributed to altered soil hydrology and albedo, which affect the surface and atmospheric temperature and soil moisture, leading to changes in thermal contrast and atmospheric circulation and thus regulating the intensity of the Asian summer monsoon rainfall (Cohen, 1994; Liu and Yanai, 2002; Jacob et al., 2005). Given the rapid thawing of permafrost and the high albedo caused by sparse evergreen coniferous trees in northern Asia during the early Holocene (Figures 6C–E), similar hydrological and atmospheric responses may have occurred, reducing the intensity of the early Holocene EASM (Figure 6G; Chen et al., 2015; Herzschuh et al., 2019; Kang et al., 2020). This mechanism, however, requires a mechanism model for further validation.