Volcanic eruptions are among the most impactful natural disturbances on Earth. They perturb ecosystems1, landscapes2,3, and climate4,5,6. Their deposits range from subtle but very widespread dust fallout to proximal tephra fall and flow sequences that can mantle topography and infill valleys to depths of tens of meters or more. Globally, more than 1300 volcanoes are known or thought to have erupted in the Holocene, and in the modern period, around 80 eruptions a year are recorded7. While there have been many studies of the impacts of volcanism on vegetation dynamics, their focus has been on the aftermath of a small number of eruptions, spanning relatively short timescales of up to decades, including Krakatau 18838, Mount Usu 19109, Surtsey 196310, and Mount St Helens 198011. Factors that influence ecological recovery include depth of burial of the seed bank in pre-eruption soils, nature of the deposit, local climate conditions, herbivore activity, and proximity to refugia12. For example, the type, amount, and distribution of biological legacies sparse on the landscape after a volcanic disturbance are key drivers of vegetation succession13,14.

Despite this progress in understanding the ecological impacts of volcanic eruptions, little remains known about longer-term post-eruption ecological succession and the influences of local and regional factors (e.g., hydroclimate, pedogenesis) as well as eruption characteristics (e.g., depth of burial, deposit characteristics)15,16,17,18. Under harsh environmental conditions (e.g., cold and/or dry habitats), the ecological repercussions of volcanic eruptions may play out over centuries and millennia19. The factors that influence a post-eruption landscape, vegetation succession, and forest recovery over such extended time periods can be challenging to disentangle, as they can reflect the history of land-use and changing climate as well as volcanological and geomorphic factors that influence recruitment and recovery20,21.

One site where the long-term consequences of a major volcanic eruption appear remarkably evident is Changbaishan (Mt. Paektu), situated on the border between PR China and DPR Korea (Fig. 1). Precisely dated to late 946 CE22 and involving rhyolitic and trachytic phases23, its so-called Millennium Eruption (ME) ranks among the largest explosive eruptions of the Common Era. It disgorged an estimated 7–36 km3 (dense-rock equivalent) of magma23 as 96 ± 19 km3 24 (40–98 km3 25) of tephra fallout and pyroclastic current deposits, the latter infilling valleys dissecting the mountain to depths in excess of 100 m. The volcanic plume dispersed eastwards more than 1000 km, resulting in widespread deposition of up to several cm thickness of ash in parts of Japan26. Today, the thickness of deposits, their hydraulic properties, the short growing season on the mountain, and low annual temperatures of the region sustain an alpine tundra environment around the summit of the volcano. These same deposits, however, contain the stems of trees, some that lived for centuries until they were killed by the eruption, suggesting profound and enduring changes in vegetation composition resulting from the ME27.

a Location of Changbaishan (orange box) on the border between PR China and DPR Korea. Image source: Google Earth. White squares indicate one tundra site (T) located at ~2050 m a.s.l.; three forest sites (F1, F2, and F3) located at ~1000 m a.s.l., where carbonized trees were sampled; the star represents the Chichi peatland, located at ~1800 m a.s.l. in a B. ermanii forest. b Topographic profile of the transect shown as red line in (a). The vegetation zones are conifer-broadleaf mixed forests (1000–1200 m a.s.l.), conifer forest (1200–1500 m a.s.l.), conifer-birch (B. ermanii) forest ecotone zone (1500–1700 m a.s.l.), birch (B. ermanii) forest (1700–1900 m a.s.l.), birch-tundra ecotone zone (1900–2060 m a.s.l.), and alpine tundra (>2060 m a.s.l.) up to Tianchi Lake (shown in blue). F1–F3 and T represent the sample sites of carbonized trees as shown in (a, c, d). The Chichi peatland is not located directly on this transect; instead, it has been plotted in (b) at its corresponding elevation for context. c Examples of carbonized trees found from the tundra area (T) and the current alpine tundra view near the site T. d Examples of carbonized trees found in the field at sites F1, F2, and F3, and current conifer-broadleaf mixed forest located at ~1000 m a.s.l. near site F3.

Full size image

The immediate impacts of the emplacement of tephra fall deposits23 were surely extreme on terrestrial ecosystems, within a radius of around 50 km from the summit of the volcano. While, despite its high magnitude, the ME appears to have muted impacts on northern hemisphere climate28, syn-eruptive volatile emissions (sulfur, halogens, and trace metals29,30) likely affected aquatic and terrestrial environments. Previous modeling studies suggested the importance of remnant trees and seed sources as well as environmental factors (e.g., climate, terrain, and soils) in influencing vegetation dynamics on the northeastern sector of Changbaishan over the past three centuries3,31, but hitherto there has been no focused study of the pre-ME and present vegetation species composition in order to evaluate the millennium-scale picture of recovery.

Here we analyzed of preserved carbonized and partially carbonized tree stems embedded in ME deposits32,33,34,35 and phytolith records in peatland cores36 to compare pre- and post-ME vegetation species composition. Based on 102 pre-eruption carbonized samples we collected at two elevations on Changbaishan (Fig. 1), we reconstructed the pre-ME forest composition on Changbaishan. The elevation of the tundra vegetation where carbonized stems were found lies well above the current alpine treeline. We hypothesize that the vegetation species composition and distribution are observably different between the two periods due to the effects of ME. We also reconstructed paleovegetation for different periods during the past millennium using phytolith records. Additionally, complementing our earlier carbonized tree-based study37, this work provides additional insights into millennial-scale climate dynamics through comprehensive reanalysis. This study builds upon our long-term scientific investigation at the study site, synthesizing previous researches and expanding the knowledge of the ecology of this area after ME. By incorporating both newly discovered materials and historical data through complementary analytical approaches, we aim to characterize and disentangle post-eruption vegetation change and global warming-induced vegetation dynamics on the volcano, and to present novel perspectives, to our knowledge, on post-ME ecological succession.

Comparing pre-ME and present forest structure

We recognize eleven tree species (91 trees) among the carbonized samples at lower elevations (~1000 m a.s.l.), including Pinus koraiensis, Larix olgensis, Picea koraiensis, Picea jezoensis, Abies holophylla, Abies nephrolepis, Betula platyphylla, Tilia amurensis, Populus ussuriensis, Fraxinus mandschurica, and Quercus mongolica (Fig. 2a and Supplementary Fig. 1). At this elevation, distinct differences in tree species composition between the pre-ME and current forest are evident. Before the ME, the forests were dominated by conifers, accounting for ~87% of individuals, particularly by P. jezoensis (22%), P. koraiensis (21%), A. nephrolepis (17%), P. koraiensis (14%), and A. holophylla (12%) (Fig. 2a). The five broadleaf species accounted for 13% of the sample set. Similar results were found for the three individual forest sites (Supplementary Fig. 1). In contrast, conifers account for ~54% of trees (Fig. 2b), and Picea koraiensis and P. jezoensis are also less abundant in the contemporary forest. On the other hand, Acer species (A. mono and A. tegmentosum) account for 17% of contemporary samples, but are not represented in the pre-ME collection. This characterizes the pre-ME low-elevation forest as a climax community of the late-successional period38. The proportion of basal areas shows very similar results with the number of trees (Fig. 2c, d vs Fig. 2a, b).

a, b Tree species composition, c, d proportion of basal areas, and e, f diameter at breast height (DBH) distributions for pre-ME (left column) and present-day forest (right column, DBH in 2013). Numbers in parentheses represent sample size.

Full size image

Moreover, all carbonized P. jezoensis, A. holophylla, and A. nephrolepis trees (accounting for more than half of the samples) have modest diameters at breast height (DBH < 25 cm) (Fig. 2e). Carbonized Pinus koraiensis is represented in all DBH classes and is dominant among larger trees (DBH > 40 cm). In contrast, the modern fir trees (accounting for 43% of samples) are represented in all DBH classes and are dominant among smaller trees (DBH < 20 cm). Moreover, current A. mono, A. tegmentosum, and B. platyphylla are primarily small individuals, whereas T. amurensis and Q. mongolica are mainly large individuals (Fig. 2f).

The most obvious differences in forest structure between pre-ME and present forests are lesser representation of small-DBH broad-leaved trees in the medieval forests. Moreover, Picea accounts for 36% of carbonized tree samples, whereas it is not represented in the contemporary population (Fig. 2). Picea only accounts for 8.7% in the range of 900–1200 m based on a large-scale forest survey of the Changbaishan Nature Reserve range (Supplementary Fig. 2), i.e., far reduced compared with the pre-ME forest. Further, about half of the trees with large DBH (>40 cm) are broad-leaved species today (Fig. 2d), whereas the vast majority of large trees (DBH > 40 cm) are conifers in the pre-ME population (Fig. 2c). The higher abundance of broad-leaved species in the modern period may result from different successional stages and/or recent climate warming.

During forest succession, changes in tree species composition depend on disturbance regime (frequency, severity), species’ traits (e.g., growth rate, shade tolerance, reproduction type), species interactions (competition, facilitation) soil features (texture, organic matter, fertility), hydrological properties of the sites, topography, aspect, and other variables. For instance, at low elevations (~1000 m a.s.l.) Changbaishan forest communities in the early stage of succession are mainly affected by soil nitrogen availability, whereas in the middle and late stage of succession soil phosphorus content plays a stronger role38. Moreover, climate change can also affect species reproduction or competition. Previous studies have shown how climate change alters tree species abundance and community composition39,40,41. For example, the relative abundance of deciduous broad-leaved trees increases near their colder range boundaries under global warming, consistent with decreases in the cover of conifers42.

Most Acer and B. platyphylla in our modern population are less than 40 years old (Supplementary Fig. 3), spanning the recent warming period, and suggesting that recruitment of these small broad-leaved trees has accelerated as regional temperatures have increased43. For example, recent climate warming resulted in increased recruitment rates of Betula ermanii in Changbaishan44. However, other large broad-leaved trees (e.g., Tilia amurensis and Quercus mongolica) exceed 150 years in age (Supplementary Fig. 3), indicating they were growing in the early industrial period and that their recruitment is related to disturbance regimes or ecological strategies, rather than recent climate warming. Moreover, mean air temperature differs little (p = 0.39) between the pre-ME and modern periods (Supplementary Fig. 4). Thus, temperature changes are not the primary reason for the pronounced differences in forest structure and composition between the two periods.

Moreover, the phytolith records from the Chichi peatland attest to mixed conifer-birch forests at 450 cal yr CE, coniferous forest before the ME at 700 cal yr CE a transition to herbaceous plants soon after the ME at 950 cal yr CE, and then a return to mixed conifer-birch forests at 1450 cal yr CE, followed by B. ermanii forest at 1700 cal yr CE (Supplementary Fig. 5). Our results align closely with the paleovegetation reconstruction from the Chichi peatland pollen record45, which identified: Mixed conifer-birch forests dominating from 250 to 1350 (or 1450) CE, and a transition to B. ermanii forest dominating from 1350 to 1950 CE. The last period occurred during the Little Ice Age. Centuries of prolonged low temperatures (Supplementary Fig. 4) led to a downslope biome shift, and the vegetation in this location changed from mixed conifer-birch forests to pure B. ermanii forest. Therefore, these differences in forest structure and composition were due to the combined effects of the ME and climate change.

Pre-ME trees in the present-day tundra zone

We found eleven carbonized trees located at and above the treeline-tundra ecotone (2030–2100 m a.s.l.) on the western side of Changbaishan in 2017 and 2018. These samples are all conifers (mostly spruce), consisting of four Picea jezoensis, four P. koraiensis, and three Pinus pumila (Fig. 3). The calibrated ages of these samples based on radiocarbon dates range from 436–656 CE to 704–1158 CE (Fig. 3 and Supplementary Table 1). That is, some of these samples were centuries old when the ME occurred (946 CE). These findings suggest that a spruce forest existed before the ME in what is now the tundra zone, which is dominated by Rhododendron aureum and other low-growing alpine shrubs and herbs (Fig. 1c). We speculate that the present tundra vegetation (Supplementary Fig.  6b, c) will be succeeded by forests under climate warming, despite the effects of slope instability and tephra erosion46.

Calibrated radiocarbon ages (±2σ, 95.4% based on IntCal20 NH) of the AMS 14C ages of the eleven carbonized trees were dated at the National Taiwan University AMS Laboratory. These carbonized trees without entire stem were found at different elevations within the alpine tundra vegetation zone on the western slope of Changbaishan. Note that the current alpine treeline is located at ~2000 m a.s.l. The error bar represents ±2σ (95.4%) for the calibration age. Different colors indicate different tree species, including Picea jezoensis (green), Picea koraiensis (blue), and Pinus pumila (yellow).

Full size image

The altitudinal vegetation zones of Changbaishan show that the elevation difference between the current treeline and the conifer forest upper limit is ~550 m (Fig. 1b). The elevation of the current treeline is ~2060 m43. Since present-day temperatures are comparable to the period prior to the ME (Supplementary Fig. 4), we assume that the treeline pre-ME would have been situated at ~2600 m. This implies that there was likely almost no alpine tundra present on Changbaishan before the ME unless the summit had been markedly higher. Changbaishan’s volcanic history includes two widely recognized explosive eruptions (Tianwenfeng and ME)47,48,49. The present-day caldera, now filled by Tianchi Lake, may be the product of both. The 1815 eruption of Tambora in is thought to have reduced the summit height by over 1400 m50. Thus, it is possible that, prior to the ME, Changbaishan was markedly higher in elevation.

Observed changes in present-day tundra vegetation

The Changbaishan alpine tundra flora consists of 22 genera and 70 species of lichen, 67 genera and 135 species of moss, and 87 genera and 131 species of vascular plants51. About 89.2% of the genera and 60.1% of the species of the Changbaishan alpine tundra flora are shared with the Arctic tundra51. The native and dominant plant species are represented by dwarf shrubs and herbs, e.g., Dryas octopetala, Vaccinium uliginosum, Rhododendron aureum, Rhododendron confertissimum, Diphasiatrum alpinum, Polygonum viviparum, Salix rotundifolia, Phyllodoce caerulea, Gentiana algida, Empetrum nigrum52,53. The northern slope of Changbaishan has the highest proportion of shrubby tundra vegetation, while the eastern slope still exposes much tephra. This likely reflects the eastward trajectory of the ME plume evident from the thick tephra fall deposits on this flank of the volcano, along with other factors such as wind erosion54. However, the shrub-dominated alpine tundra vegetation (e.g., R. aureum) has synchronously changed on all slopes of Changbaishan over the past four decades (Fig. 4), although the vegetation development after volcanic disturbance differs according to aspect.

The changes are the areal fraction of R. aureum, herbs, and ash (unvegetated areas) in the tundra over the last four decades as a function of elevation and aspect: a–d represent eastern, southern, western, and northern slopes, respectively. Similar to the changes in tundra vegetation, changes in Betula ermanii are consistent across different aspects, thus we consider trends in recruitment rate for B. ermanii trees in the large plot within the western tundra zone during 1985–2015 as representative (purple in (c)) (See “Methods” for details). Two (one) asterisks indicate changes significant at p < 0.05 (p < 0.1).

Full size image

In general, the area dominated by herbs (e.g., Deyeuxia angustifolia) increased55,56 at the expense of R. aureum on all slopes and between 1800 and 2400 m (Fig. 4). Notably, D. angustifolia was not found in the tundra of Changbaisan before the 1980s52, being recorded since the late 1980s51. The R. aureum area slightly increased at elevations above 2400 m on the eastern and northern slopes (Fig. 4a, d). Moreover, the changes in area were greater at lower than higher elevations. For example, the decrease of R. aureum was 28, 9, 4, and 2% per decade, and the increase of herbs was 27.5, 8, 2.5, and 0.8% per decade, from low to high elevations on western slopes (Fig. 4c). The changes in tephra-covered and unvegetated areas were small, but relatively greater at higher than lower elevations.

To further investigate the changes of vegetation at high elevations, we established a large plot of 450 m in slope length (2074–2194 m) × 50 m in width in the tundra on the western slope in 2015 and investigated tree recruitment (Supplementary Fig. 6). Our earlier study found that, 598 B. ermanii seedlings/saplings (See “Methods”) were being recruited in the tundra, nearly all of which (99%) had established since 198543 and coinciding with rapid climate warming. This study found a significant (p < 0.01) increasing recruitment rate (~32% per decade) for all 598 trees in the plot during 1985–2015 (Fig. 4c). We also found that the recruitment of B. ermanii was significantly (p < 0.05) positively related to summer mean temperature during this period (Supplementary Fig. 7). Previous studies have showed that recent warming has accelerated the upward shift of the treeline on Changbaishan43 and that the trend is likely to continue57. Alpine shrubs (e.g., R. aureum) are being replaced by herbs and trees (B. ermanii) (Supplementary Fig. 6b, c). The likely consequence is a decline or even disappearance of tundra vegetation, which constitutes the southern fringe of the Eurasian tundra in eastern Asia58.

Our study documents the dramatic changes in vegetation species composition that remain evident when comparing vegetation cover before and more than 1000 years after the 946 CE Millennium Eruption on Changbaishan. The present-day alpine tundra on Changbaishan, the most southernmost in eastern Asia, is shrinking, and we have documented an ongoing ecological succession that has persisted for a millennium following the ME. Therefore, past, present, and future changes in forests including the alpine treeline elevation and species composition are not only driven by their response to climate change58,59,60 or by recent disturbances61, but can also be partly attributed to the long-term impacts of the Millennium Eruption on Changbaishan. This highlights the need to account for the potentially enduring impacts of major past ecological disturbances when evaluating the effects of recent climate change on vegetation dynamics.

Carbonized trees sampling and laboratory treatment

The ME of the Changbaishan volcano produced a large amount of ash that pushed and buried forest trees, carbonizing them. Some of the carbonized trees were exposed due to the construction of roads or water erosion. We sampled 102 carbonized trees at two elevations, 91samples at ~1000 m a.s.l. and 11 ones at ~2015 m a.s.l. The lower elevation included three forest sites (site F1, 42.14993°N, 127.86569°E, 1025 m a.s.l., with 33 samples collected in 2012; site F2, 41.71495°N, 128.16063°E, 1185 m a.s.l., with 36 samples collected in 2012; and site F3, 42.09444°N, 127.70702°′E, 892 m a.s.l., with 22 samples collected in 2012), and another one was a tundra site (site T, 41.98815°N, 128.00491°E, ~2050 m a.s.l., with 11 samples but without an entire stem collected in 2016) (Fig. 1). The carbonized trees were exposed and then found due to road construction for sites F2 and T, and due to water erosion for sites F1 and F3. Many of the sampled trees were not totally carbonized and were whole stems with bark (Supplementary Fig. 8a–c), indicating limited transport from their growth site. Moreover, the tree species of all of these carbonized trees can be exactly identified (see below). The tree age of all of the carbonized trees for the three forest sites can also be exactly identified. Therefore, analyzing the forest composition and structure of pre-ME based on the carbonized samples avoids common taphonomic biases of plant fossil records62. This study represents the first application of these carbonized samples for analyses of forest composition and vegetation succession dynamics—providing novel insights, to the best of our knowledge, into long-term ecological processes.

We identified the tree species of carbonized trees by analyzing microscopic anatomical features of wood on three planes (cross-sectional, radial, and tangential). For example, the microscopic anatomical features of one carbonized Pinus koraiensis wood are shown in Supplementary Fig. 8d–g. We also dated the samples through radiocarbon analyses of the outermost rings of three carbonized trees from sites F1–F3 at the Accelerator Mass Spectrometry Laboratory in Peking University (Supplementary Table 1) and all eleven samples from site T at the National Taiwan University AMS Laboratory (Fig. 3). The samples subjected to radiocarbon dating were collected during two distinct sampling years: 2012 for forest samples at sites F1–F3 and 2016 for tundra samples at site T. The use of different laboratories primarily reflects the temporal gap between these sampling years (2012 vs. 2016), rather than any methodological inconsistencies. The measured AMS 14 C ages are consistent with previous studies of dating the ME22,63, indicating the carbonized trees we collected were killed during the ME (late 946 CE).

Peat phytolith records

We reconstructed vegetation type over the past 1000 years using peat phytolith records of the Chichi (a pond) peatland collected by our co-author Prof. Dongmei Jie36. The Chichi peatland is located in the current Betula ermanii forest on the north slope of Changbaishan (42.05445°N, 128.05611°E, ~1800 m a.s.l.). The Chichi peatland is relatively small, with a circumference of approximately 220 m and a diameter of approximately 80 m. The pond within it is relatively shallow, and a peatbog with a thickness of about 0.40–0.85 m has developed, making it the highest-elevation peat bog in the Northeast China. The pond area is surrounded by Carex lasiocarpa, along with Sphagnum, Drosera species and other shrubs. The dominant community in the surrounding environment is a Betula ermanii-Rhododendron aureum community. We sampled peat profile for 30 cm in the Chichi peatland. We extracted peat samples at five depths, 30, 25, 20, 10, and 5 cm from bottom to top of the profile. The five depths represent the date of 450 CE, 700 CE, 950 CE, 1450 CE, and 1700 CE, respectively64. We identified phytolith types for the peat samples and analyzed the characteristics of phytolith assemblages in different samples (different periods) (Supplementary Fig. 5). Herbs, conifers, and broad-leaved plants are identified by different phytolith morphotypes (see Supplementary Fig. 5). We then reconstructed palaeovegetation of the Chichi peatland for the four periods from 450 CE onward by analyzing the phytolith assemblages in the peat samples to the parent plants they indicated.

Current forest structure

Near to the three forest sites used for sampling carbonized trees, we established three sample square plots (30 m × 30 m) in old-growth stands to investigate the present forest structure. These forests have not been disturbed by human activities or fires. Because of the strong intensity of the ME, there may be some downslope displacement for carbonized trees. Therefore, the present sampling plots were located within 5 km from the carbonized tree sampling plots in the upslope direction. Plots F1 and F2 were established in 2013, and plot F3 was established in 2023 (Fig. 1). We recorded the species, height, diameter at breast height (DBH, measured at 1.3 m) for all trees with height >2 m. Twelve tree species were found in the three plots, including Pinus koraiensis, L. olgensis, A. holophylla, A. nephrolepis, B. platyphylla, T. amurensis, Albizia kalkora, Ulmus pumila, A. mono, A. tegmentosum, Fraxinus mandshurica, and Q. mongolica. We collected at least two cores from each tree of all species using Pressler increment borers for determining their ages by counting the number of rings. The DBH and tree age of trees in plot F3 were extrapolated to 2013 by subtracting radial growth and age over the last 10 years. All individuals in these three plots were used to compare present and ME forest structures.

Additionally, we used a forest inventory dataset of Changbaishan in 2017. This dataset was derived from 333 forest survey sample plots (50 m × 50 m) and contained a total of 68,137 trees from nine conservation stations of the Changbaishan National Nature Reserve, where the disturbance of forests by human activities is limited. This large-scale survey represents regional forest composition. For comparability with the present study, we extracted data from similar low elevations (900–1200 m) with the sampling sites of the carbonized and present trees in this study, including 28,577 trees from 144 sample plots. We classified these trees into broad-leaved trees and conifers and analyzed the distribution characteristics of their DBH (Supplementary Fig. 2).

Changes in the area ratio of different vegetation types in tundra

We used remote-sensing images (Landsat satellites 5 and 7 and GaoFen-2 (GF-2) images taken at 30 and 0.8 m resolution, respectively; See Supplementary Table 2) to extract the distribution of different vegetation types near the tundra from the 1980s to 2010s (1988 to 2017) in areas with elevation >1800 m, where the forest-tundra ecotone and tundra zone are located in Changbaishan (Fig. 1a). The GF-2 image was purchased from the Beijing Lanyu Fangyuan Information Technology Co., Ltd. The Landsat satellites 5 and 7 images were downloaded at https://earthexplorer.usgs.gov/. All used images were very clear, most cloud free, while others showed very low cloud cover (1–3%, and only one for 9%) (Supplementary Table 2).

Rhododendron aureum is the dominant tundra evergreen shrub. We identified R. aureum from herbaceous plants using spectral differences at specific phenological periods. For example, the reflectance curves of R. aureum and Deyeuxia angustifolia (the dominant herbaceous plant) differ markedly during the yellowing and withering phenological phases in autumn (late September to early October). The NDVI of herbaceous plants was markedly lower than that of R. aureum during the yellowing and withering phenological phases. Volcanic ash and unvegetated areas have the lowest NDVI. Therefore, this study focused on the R. aureum and herbaceous plant abundances.

We first mapped the outline of R. aureum and volcanic ash (containing non-vegetation covers) using true color composite images derived from the GF-2 image. In this image, green is recognized as R. aureum, whereas volcanic ash appears as white patches. The remaining area was considered as the distribution of herbaceous plants. We discriminated vegetation types using unsupervised classification and visual interpretation using the GF-2 image acquired on 23 September 2017. Then, we overlaid the vegetation classes on the near contemporaneous Landsat 7 image (5 October 2017). We found that areas with >80% of pixels covering R. aureum had a NDVI (measured with the Landsat image) value greater than 0.2, whereas for ash it was lower than 0.0. Therefore, the pixels with NDVI values > 0.2, 0.0–0.2, and <0.0 over the historical Landsat images were considered as R. aureum, herbaceous plants, and volcanic ash, respectively. Based on this, we calculated the changes in area ratio (the fraction of pixels for a given class) of R. aureum, herbaceous plants, and volcanic ash spanning the 1980s–2010s for different elevations and slopes of Changbaishan (Fig. 4).

Recruitment rate of Betula ermanii in tundra

A large plot (450 m × 50 m, Supplementary Fig. 6a) was established by our co-authors in 2015 in the alpine tundra (2074–2194 m a.s.l.) to analyze whether and how fast B. ermanii recruited beyond the treeline. 598 B. ermanii trees were found in the large plot, most of which were smaller than 1.3 m in height, and therefore, we measured diameter at base (DB, approximately 5 cm above the root crown) for all detected seedlings and saplings. These trees are too small to determine ages by extracting cores. Therefore, the DB and sampled the trunks of B. ermanii trees (totally 35 ones) within 2 m of a long side of the large plot (blue points in Supplementary Fig. 6a) were measured, and their ages were determined by counting their rings. However, there may exist potential differences in the radial growth-age relationship among different elevations. Therefore, we also measured the DBH and determined tree age by drilling tree cores for trees grown at elevations (1990–2050 m a.s.l.) below the alpine treeline (Supplementary Fig. 9). The tree age–DBH (or DB) regression models did not differ across elevations (Supplementary Fig. 9). Therefore, the age-DB regression model by the 35 B. ermanii trees (age = 5.44 × DB + 0.49, R2 = 0.645, p < 0.01) is a robust model, enabling estimation of the age of all individuals found in the plot. The recruitment time was determined by tree age. We summed the number of new recruits each year and divided these annual cumulative values by the total number of trees (598), then multiplied by 100% to obtain the annual cumulative recruitment percentage. Approximately 99% of the trees were recruited after 1985. We therefore calculated the change in this cumulative percentage to derive the trend in recruitment rate for B. ermanii in the large plot during 1985–2015 (Fig. 4c). We also analyzed the relationship between birch recruitment and summer mean temperature during this period (Supplementary Fig. 7).