Topographic Controls on Vegetation Changes in Alpine Tundra of the Changbai Mountains

The vegetation of alpine tundra is undergoing significant changes and topography has played a significant role in mediating such changes. The roles of topography varied at different scales. In this study, we intended to identify topographic controls on tundra vegetation changes within the Changbai Mountains of Northeast China and reveal the scale effects. We delineated the vegetation changes of the last three decades using the normalized difference vegetation index (NDVI) time series. We conducted a trend analysis for each pixel to reveal the spatial change and used binary logistic regression models to analyze the relationship between topographic controls at different scales and vegetation changes. Results showed that about 30% of tundra vegetation experienced a significant (p < 0.05) change in the NDVI, with 21.3% attributable to the encroachment of low-altitude plants resulting in a decrease in the NDVI, and 8.7% attributable to the expansion of tundra endemic plants resulting in an increase in the NDVI. Plant encroachment occurred more severely in low altitude than in high altitude, whereas plant expansion mostly occurred near volcanic ash fields at high altitude. We found that plant encroachment tended to occur in complex terrains and the broad-scale mountain aspect had a greater effect on plant encroachment than the fine-scale local aspect. Our results suggest that it is important to include the mountain aspect in mountain vegetation change studies, as most such studies only use the local aspect.


Introduction
Alpine tundra is located in relatively isolated regions and is very sensitive to global environmental change [1,2].Alpine tundra in many places around the world has been undergoing dramatic changes in vegetation communities over the past few decades [3,4].Some changes are considered to be the results of climate warming, while others are attributed to nitrogen deposition.Hallinger et al. (2010) found a sharp expansion of shrubs in the alpine tundra of European Alps due to increasing temperature [5].Climate warming in the Swiss Alps reduced snow cover and advanced snowmelt, which altered aboveground biomass allocation [4].In subalpine tundra in northern Sweden, increased nitrogen caused the replacement of shrubs by fast-growing grasses in this nutrient-poor system [6].Nitrogen deposition also led to changes of species composition in the alpine tundra of the southern Rocky Mountains [7].Changes in the structure of the plant community may be strongly mediated by topography in alpine tundra, a biome that typically has complex topography.Altitude, and its tight correspondence with temperature, is the obvious topographic variable that can affect vegetation growth [8,9].Slope and aspect, features that regulate solar radiation and soil water content, are also important topographic controls [10,11].Topography often causes spatial heterogeneity in soil temperature and moisture by affecting the spatial redistribution of atmospheric temperature, moisture, and solar radiation [8,12].Studies have shown that topography affects the rate and magnitude of vegetation changes in alpine areas [3,13,14].Carmel and Kadmon (1999) observed that the rate of vegetation changes was affected by slope and aspect in the Northern Galilee Mountains of Israel [13].Similarly, Gottfried et al. (1999) found that the distribution of plant species in the European Alps shifted to higher elevations under specific terrain conditions in response to climate warming [15].
Topographic variation exhibits strong effects on alpine plant communities.High mountains present a significant barrier to atmospheric circulation, which may alter airflow, resulting in different precipitation and temperature regimes depending on the aspect of the mountain [16,17].In addition, different aspects of the mountain receive different amounts of solar radiation.Therefore, at the mountain-wide scale, topography factors, especially the mountain aspect, correspond to different environments for plant growth.Meanwhile, at local scales, the aspect of the immediate surroundings also affects habitat conditions by affecting the microclimate [18][19][20].Studies have suggested that highly complex topographies could create favorable conditions for plants to cope with changing environmental conditions, such as climate warming [21][22][23].Vegetation changes are closely related to changes in environmental conditions [24][25][26] and, thus, patterns of vegetation changes may be affected by both the local aspect and mountain aspect.However, most studies about the relationship between topography and changes of mountainous vegetation were conducted at local, micro-topographical scales, rarely accounting for mountain-scale effects [14,15].Particularly, the mountain aspect as a mountain-scale topographical variable is often ignored in related researches, whose effect should be worthy of attention.
The tundra communities of the Changbai Mountains in China are typical of temperate alpine tundra.When comparing historical records with recent survey results, there have been obvious plants community changes since the 1980s.[27,28].The most significant changes include rapid upward shifts of herbaceous plants represented by Narrowleaf small reed (Deyeuxia angustifolia) on the west side of the mountain [27] and colonization of birch (Betula ermanii)-A treeline tree species-On the north side of the mountain [28].The encroachment of low-altitude plants has dramatically changed the pattern of vegetation in the tundra, and studies indicate that the distribution of the encroaching plants was affected by topography.For example, the abundance of encroaching herbs decreased with increasing elevation [29], and the colonization by B. ermanii was influenced by the curvature of the terrain and aspect [28].These results indicate that there are differences in how vegetation changes, not only on different aspects of a single mountain, but also within different topographic niches on the same mountain aspect.Thus, the alpine tundra of the Changbai Mountains is an ideal platform to explore the relationship between vegetation changes and multiscale topographic factors.
In this study, we analyzed the temporal and spatial distributions of vegetation within the alpine tundra of the Changbai Mountains using a vegetation index time series derived from satellite images.We particularly intended to reveal the topographic controls on vegetation changes and explore the importance of the mountain aspect in mountainous vegetation researches.The conclusions of our study are expected to offer useful information regarding alpine vegetation dynamics.

Study Area
Our study area covers the entire alpine tundra of the Changbai Mountains, which is located along the southern edge of the Eurasian alpine tundra zone, and developing on the upper slopes of a Forests 2018, 9, 756 3 of 13 volcanic cone (Figure 1).The climate of the alpine tundra is characterized by low temperatures, intense precipitation, and strong winds, with an average annual temperature of about −7.3 • C and annual precipitation of over 1400 mm [30].About 70% of the annual precipitation is concentrated between June and September [30].The tundra experiences strong westerly winds year-round, with the prevailing winds coming from the west-southwest.Tundra plants are dominated by dwarf shrubs, including Rhododendron (Rhododendron aureum), Blueberry (Vaccinium uliginosum), and Mountain avens (Dryas octopetala var.asiatica).In addition, there are some moist herbs in the tundra, e.g., Sanguisorba sitchensis (Sanguisorba stipulata) and Sanguisorba parviflora (Sanguisorba tenuifolia), mainly occurring within low-lying topographies, such as gullies and depressions.Recently, it has been found that the upward shift of herbaceous plants from low altitudes is seriously encroaching upon the habitat of native shrubs [27,29].Topographic features are highly heterogeneous, and the topographical characteristics of the mountain aspect differ greatly.The terrain on the east side is relatively flat, while the other three sides are highly rugged and characterized by crisscross gullies.Our study area covers the entire alpine tundra of the Changbai Mountains, which is located along the southern edge of the Eurasian alpine tundra zone, and developing on the upper slopes of a volcanic cone (Figure 1).The climate of the alpine tundra is characterized by low temperatures, intense precipitation, and strong winds, with an average annual temperature of about −7.3 °C and annual precipitation of over 1400 mm [30].About 70% of the annual precipitation is concentrated between June and September [30].The tundra experiences strong westerly winds year-round, with the prevailing winds coming from the west-southwest.Tundra plants are dominated by dwarf shrubs, including Rhododendron (Rhododendron aureum), Blueberry (Vaccinium uliginosum), and Mountain avens (Dryas octopetala var.asiatica).In addition, there are some moist herbs in the tundra, e.g., Sanguisorba sitchensis (Sanguisorba stipulata) and Sanguisorba parviflora (Sanguisorba tenuifolia), mainly occurring within low-lying topographies, such as gullies and depressions.Recently, it has been found that the upward shift of herbaceous plants from low altitudes is seriously encroaching upon the habitat of native shrubs [27,29].Topographic features are highly heterogeneous, and the topographical characteristics of the mountain aspect differ greatly.The terrain on the east side is relatively flat, while the other three sides are highly rugged and characterized by crisscross gullies.In our study, we regarded the area above 1800 m asl in the Changbai Mountains as the study area (Figure 1).However, due to frequent landslides near the crater and numerous manmade buildings in the canyon on the north side, we excluded these two areas from the study to insure that results could be attributed to the same causes.

General Approach
We analyzed vegetation changes in the study area using normalized difference vegetation index (NDVI) values derived from Landsat Thematic Mapper (TM)/Enhanced Thematic Mapper Plus (ETM+) time series during 1988-2017 (Figure 2).We acquired six topographic variables from a digital elevation model (DEM).Finally, we used binary logistic models to quantify relationships between vegetation changes (encroachment and expansion) and topographic variables.In our study, we regarded the area above 1800 m asl in the Changbai Mountains as the study area (Figure 1).However, due to frequent landslides near the crater and numerous manmade buildings in the canyon on the north side, we excluded these two areas from the study to insure that results could be attributed to the same causes.

General Approach
We analyzed vegetation changes in the study area using normalized difference vegetation index (NDVI) values derived from Landsat Thematic Mapper (TM)/Enhanced Thematic Mapper Plus (ETM+) time series during 1988-2017 (Figure 2).We acquired six topographic variables from a digital elevation model (DEM).Finally, we used binary logistic models to quantify relationships between vegetation changes (encroachment and expansion) and topographic variables.

Data Source and Processing
Landsat TM/ETM+ datasets were used in the study, which were obtained from the USGS (United States Geological Survey) [31].The spatial resolution was 30 m. Due to similar spectral response profiles on each spectral band, the radiometry of images acquired by TM and ETM sensors was comparable.Eleven Landsat TM/ETM+ scenes were selected from 1988 to 2017 at the senescence periods (late September to early October) for our study (Table 1), and all images of the study area were clear.Firstly, we conducted radiometric calibration and atmospheric correction to acquire surface reflectance values using ENVI 5.0 (Harris Geospatial Solutions, Broomfield, CO, USA).Subsequently, we performed projection transformation and mask extraction using ArcGIS 10.2 (ESRI, Redlands, CA, USA).We used a GaoFen-2 (GF-2) image with high resolution (0.8 m) to analyze vegetation pattern characteristics in the study area.The cloud-free image of our study area was acquired on 23 September 2017.

Data Source and Processing
Landsat TM/ETM+ datasets were used in the study, which were obtained from the USGS (United States Geological Survey) [31].The spatial resolution was 30 m. Due to similar spectral response profiles on each spectral band, the radiometry of images acquired by TM and ETM sensors was comparable.Eleven Landsat TM/ETM+ scenes were selected from 1988 to 2017 at the senescence periods (late September to early October) for our study (Table 1), and all images of the study area were clear.Firstly, we conducted radiometric calibration and atmospheric correction to acquire surface reflectance values using ENVI 5.0 (Harris Geospatial Solutions, Broomfield, CO, USA).Subsequently, we performed projection transformation and mask extraction using ArcGIS 10.2 (ESRI, Redlands, CA, USA).We used a GaoFen-2 (GF-2) image with high resolution (0.8 m) to analyze vegetation pattern characteristics in the study area.The cloud-free image of our study area was acquired on 23 September 2017.The DEM data were generated from the PRISM (panchromatic remote-sensing instrument for stereo mapping) sensor attached to the ALOS (advanced land observing) satellite, and are considered to be the most precise global-scale elevation data.The spatial resolution is 5 m, and the horizontal and vertical accuracy is also up to 5 m.

Landsat Vegetation Index
NDVI values are widely used in the research of vegetation ecology.They are calculated using red and near-infrared band reflectance: where ρ NIR and ρ R are the spectral reflectance values in the near-infrared band and red band, respectively.

Topographic Variables
We derived six topographic variables from digital elevation model (DEM) data with ArcGIS 10.2 (Table 2).Altitude, slope, and aspect are the three topographic variables most frequently used in previous studies of the relationships between vegetation and topography.In addition, we calculated the coefficient of variation in elevation (CVE) to characterize topographic complexity, and distance to the edge of the nearest volcanic ash deposit (DVA) to quantify the extent of interference by volcanic ash.In this study, we used two indices relating to aspect.The first described the aspect of the mountain as determined by the pixel orientation in relation to the summit of the mountain.The second index of aspect described the local aspect as indicated by the orientation of each pixel.Particularly, the aspect was transformed according to Beers et al. (1996) [32].DEM represents digital elevation model, CVE represents the coefficient of variation in elevation, DVA represents the distance to the nearest volcanic ash deposit, α represents the azimuth of the mountain aspect, which is the pixel orientation in relation to the summit of the mountains, and β represents the azimuth of the pixel aspect.

Trend Analysis
To analyze the NDVI changes observed between 1988 and 2017, we applied ordinary least-squares regression (OLS), which has been widely used in other studies of vegetation dynamics.The slope of the regression shows the extent of changes in the NDVI.
where θ slope is the slope of OLS, n is the number of years in the study period (from 1988 to 2017), i is the serial number of the year in the study period, and p i is the NDVI value in year i.For each pixel, Forests 2018, 9, 756 6 of 13 the slope and statistical significance of the regression were calculated using a Student's t-test at the 95% confidence level.

Binary Logistic Regression
Using the six topographic variables, we built a binary logistic regression model to predict changes in NDVI values over time.
where p is the probability of a significant change in NDVI, x i is the i-th explanatory variable, a is a constant, and b i is the regression coefficient of explanatory variable x i .All explanatory variables were normalized before entering the logistic regression model, so as to assess the relative contribution of each variable.
The regression was calculated in a forward stepwise manner using the SPSS 21.0 (IBM Corporation, Armonk, NY, USA) 'binary logistic regression' procedure, which is based on a maximum-likelihood ratio test.In our study, we divided samples into 70% training samples to build the model and 30% verification samples to test the model.The reliability of the regression model was determined with a classification table, in which the number of correctly predicted samples for each class indicates overall predictive accuracy of the model.The classification table was constructed based on validation samples to provide relative trustworthy accuracy metrics.The area under the receiver operating characteristic (ROC) curve was assessed to quantify the performance of the regression model: When the area under the curve (AUC) equals or exceeds 0.7, the model is considered to be accurate.We reported odds ratios (OR) to convey the relative impact of the explanatory variables.
To eliminate multicollinearity caused by the correlation among explanatory variables, we performed correlation analysis on six topographic variables, and the weak relationship between the variables allowed us to keep them.In the regression model, spatial autocorrelation may inflate statistical significance [33].For the high-resolution analysis of our study, however, we expected that the similarity between neighboring pixels was caused by similar topographic conditions [34].To avoid the correlations resulting from spatial proximity, we built the model based on 5000 random points using a spatially stratified sampling scheme by controlling the minimum distance between points within the study area.

Extraction of Volcanic Ashes and R. aureum Distribution Areas
The spatial distribution of volcanic ash deposits and R. aureum were determined using true color composite images derived from the GF-2 image.In this image, volcanic ash showed as white, while R. aureum, as the only evergreen plant in the tundra, showed a distinctive green.Therefore, we extracted the distribution areas of volcanic ash and R. aureum through unsupervised classification and visual interpretation of the GF-2 image (Figure 3).We consider this extraction to be reliable because of the high spatial resolution of the image.Meanwhile, the areas of volcanic ash we extracted were similar to the areas extracted by other published results [35].because of the high spatial resolution of the image.Meanwhile, the areas of volcanic ash we extracted were similar to the areas extracted by other published results [35].

Vegetation Changes
There are obvious spectral differences between typical tundra plants at the senescence period, due to their distinct phenological characteristics.During this period, R. aureum has the highest Landsat NDVI value as the only evergreen shrub in the tundra.The remaining plants (e.g., herbaceous plants) had a significantly lower NDVI than R. aureum due to yellowing of the leaves at the senescence period.This was confirmed by the distinct reflectance curves of R. aureum and Deyeuxia angustifolia at the senescence period (Figure 4).In addition, volcanic ashes or areas without vegetation cover have the lowest Landsat NDVI, which clearly differed from the areas covered by vegetation.Accordingly, we determined the vegetation pattern in the tundra using remote sensing data.Through overlaying the distribution of R. aureum and volcanic ashes (Figure 3b) on the Landsat image acquired on 5 October 2017 nearest in time to the GF-2 image, we found that Landsat NDVI values indicating areas of pure R. aureum (>80% of pixels covered by R. aureum) were greater than 0.2 (Figure 5a), and Landsat NDVI values indicating the area of extensive volcanic ash (>80% of pixels

Vegetation Changes
There are obvious spectral differences between typical tundra plants at the senescence period, due to their distinct phenological characteristics.During this period, R. aureum has the highest Landsat NDVI value as the only evergreen shrub in the tundra.The remaining plants (e.g., herbaceous plants) had a significantly lower NDVI than R. aureum due to yellowing of the leaves at the senescence period.This was confirmed by the distinct reflectance curves of R. aureum and Deyeuxia angustifolia at the senescence period (Figure 4).In addition, volcanic ashes or areas without vegetation cover have the lowest Landsat NDVI, which clearly differed from the areas covered by vegetation.Accordingly, we determined the vegetation pattern in the tundra using remote sensing data.
because of the high spatial resolution of the image.Meanwhile, the areas of volcanic ash we extracted were similar to the areas extracted by other published results [35].

Vegetation Changes
There are obvious spectral differences between typical tundra plants at the senescence period, due to their distinct phenological characteristics.During this period, R. aureum has the highest Landsat NDVI value as the only evergreen shrub in the tundra.The remaining plants (e.g., herbaceous plants) had a significantly lower NDVI than R. aureum due to yellowing of the leaves at the senescence period.This was confirmed by the distinct reflectance curves of R. aureum and Deyeuxia angustifolia at the senescence period (Figure 4).In addition, volcanic ashes or areas without vegetation cover have the lowest Landsat NDVI, which clearly differed from the areas covered by vegetation.Accordingly, we determined the vegetation pattern in the tundra using remote sensing data.Through overlaying the distribution of R. aureum and volcanic ashes (Figure 3b) on the Landsat image acquired on 5 October 2017 nearest in time to the GF-2 image, we found that Landsat NDVI values indicating areas of pure R. aureum (>80% of pixels covered by R. aureum) were greater than 0.2 (Figure 5a), and Landsat NDVI values indicating the area of extensive volcanic ash (>80% of pixels Through overlaying the distribution of R. aureum and volcanic ashes (Figure 3b) on the Landsat image acquired on 5 October 2017 nearest in time to the GF-2 image, we found that Landsat NDVI values indicating areas of pure R. aureum (>80% of pixels covered by R. aureum) were greater than 0.2 (Figure 5a), and Landsat NDVI values indicating the area of extensive volcanic ash (>80% of pixels covered by volcanic ash) were less than 0.0 (Figure 5b).Therefore, we divided NDVI values into three classes: <0.0 (volcanic ash or no vegetation cover), >0.2 (evergreen vegetation cover of R. aureum), and 0.0-0.2(largely non-evergreen, herbaceous plants cover).Based on the NDVI pattern maps during the 1990s-2010s derived from the averages of Landsat NDVI from the relevant decade, we found the volcanic ash was mostly distributed in the southeast side of the mountain, decreasing in total area from 10.7% to 9.7% between 1988 and 2017 (Figure 6).The areas covered by R. aureum had a marked decrease from 70.6% in the 1980s to 31.7% in the 2010s, especially in low-altitude tundra (Figure 6).Non-evergreen plants had a significant expansion, increasing from 18.7% to 58.6% of the total area (Figure 6).covered by volcanic ash) were less than 0.0 (Figure 5b).Therefore, we divided NDVI values into three classes: <0.0 (volcanic ash or no vegetation cover), >0.2 (evergreen vegetation cover of R. aureum), and 0.0-0.2(largely non-evergreen, herbaceous plants cover).Based on the NDVI pattern maps during the 1990s-2010s derived from the averages of Landsat NDVI from the relevant decade, we found the volcanic ash was mostly distributed in the southeast side of the mountain, decreasing in total area from 10.7% to 9.7% between 1988 and 2017 (Figure 6).The areas covered by R. aureum had a marked decrease from 70.6% in the 1980s to 31.7% in the 2010s, especially in low-altitude tundra (Figure 6).Non-evergreen plants had a significant expansion, increasing from 18.7% to 58.6% of the total area (Figure 6).Based on the rate of change between 1988 and 2017, and levels of statistical significance, we classified NDVI change into three types: Increased, decreased, and unchanged.About 30% of the area showed statistically significant changes in NDVI values (p < 0.05), with NDVI values decreasing (encroachment of low-altitude plants) in 21.3% of areas and increasing (expansion of endemic plants) in 8.7% of areas (Figure 7).The remaining area (70%) had no changes.The encroachment of low-altitude plants occurred mainly in the low-altitude areas on the south and west sides of the study area.The expansion of endemic plants largely occurred near the volcanic ash areas.The areas of unchanged vegetation were concentrated on the east side of the mountain (Figure 7).covered by volcanic ash) were less than 0.0 (Figure 5b).Therefore, we divided NDVI values into three classes: <0.0 (volcanic ash or no vegetation cover), >0.2 (evergreen vegetation cover of R. aureum), and 0.0-0.2(largely non-evergreen, herbaceous plants cover).Based on the NDVI pattern maps during the 1990s-2010s derived from the averages of Landsat NDVI from the relevant decade, we found the volcanic ash was mostly distributed in the southeast side of the mountain, decreasing in total area from 10.7% to 9.7% between 1988 and 2017 (Figure 6).The areas covered by R. aureum had a marked decrease from 70.6% in the 1980s to 31.7% in the 2010s, especially in low-altitude tundra (Figure 6).Non-evergreen plants had a significant expansion, increasing from 18.7% to 58.6% of the total area (Figure 6).Based on the rate of change between 1988 and 2017, and levels of statistical significance, we classified NDVI change into three types: Increased, decreased, and unchanged.About 30% of the area showed statistically significant changes in NDVI values (p < 0.05), with NDVI values decreasing (encroachment of low-altitude plants) in 21.3% of areas and increasing (expansion of endemic plants) in 8.7% of areas (Figure 7).The remaining area (70%) had no changes.The encroachment of low-altitude plants occurred mainly in the low-altitude areas on the south and west sides of the study area.The expansion of endemic plants largely occurred near the volcanic ash areas.The areas of unchanged vegetation were concentrated on the east side of the mountain (Figure 7).Based on the rate of change between 1988 and 2017, and levels of statistical significance, we classified NDVI change into three types: Increased, decreased, and unchanged.About 30% of the area showed statistically significant changes in NDVI values (p < 0.05), with NDVI values decreasing (encroachment of low-altitude plants) in 21.3% of areas and increasing (expansion of endemic plants) in 8.7% of areas (Figure 7).The remaining area (70%) had no changes.The encroachment of low-altitude plants occurred mainly in the low-altitude areas on the south and west sides of the study area.The expansion of endemic plants largely occurred near the volcanic ash areas.The areas of unchanged vegetation were concentrated on the east side of the mountain (Figure 7).

Topographic Controls on Vegetation Changes
The logistic regression model that predicted the encroachment of low-altitude plants had high accuracy (80.2%), with an AUC of 0.78 (Table 3).Four of the six topographic variables (altitude, mountain aspect, CVE, and local aspect) contributed significantly to the model.Altitude was the most dominant factor; the negative coefficient indicates that encroachment occurs at low altitudes.The positive effect of CVE indicates that encroaching plants are mostly distributed in complex terrain.The intensity and direction of effects due to the mountain aspect and local aspect differed.The mountain aspect was the second most important topographic factor; the negative correlation reflected the observation that encroaching plants were mostly distributed on the south and west sides of the mountain.However, the local aspect had a positive effect, showing that plant encroachment occurred commonly in the pixels facing north and east.Although the variable were statistically significant in the regression model, odds ratios (OR) near 1 indicated small effect sizes of local aspect on the encroachment of low-altitude plants.The model that predicted endemic plant expansion also had high predictive accuracy (91.4%), with an AUC of 0.77 (Table 3).The regression model incorporated two topographic variables: DVA and CVE showed negative effects on endemic plants expansion.DVA was the main contributor to

Topographic Controls on Vegetation Changes
The logistic regression model that predicted the encroachment of low-altitude plants had high accuracy (80.2%), with an AUC of 0.78 (Table 3).Four of the six topographic variables (altitude, mountain aspect, CVE, and local aspect) contributed significantly to the model.Altitude was the most dominant factor; the negative coefficient indicates that encroachment occurs at low altitudes.The positive effect of CVE indicates that encroaching plants are mostly distributed in complex terrain.The intensity and direction of effects due to the mountain aspect and local aspect differed.The mountain aspect was the second most important topographic factor; the negative correlation reflected the observation that encroaching plants were mostly distributed on the south and west sides of the mountain.However, the local aspect had a positive effect, showing that plant encroachment occurred commonly in the pixels facing north and east.Although the variable were statistically significant in the regression model, odds ratios (OR) near 1 indicated small effect sizes of local aspect on the encroachment of low-altitude plants.The model that predicted endemic plant expansion also had high predictive accuracy (91.4%), with an AUC of 0.77 (Table 3).The regression model incorporated two topographic variables: DVA and CVE showed negative effects on endemic plants expansion.DVA was the main contributor to endemic plant expansion, which mostly occurred near the volcanic ash areas.In addition, the expansion commonly occurred in relatively flat areas, e.g., the east side of the mountain.

Discussion
Our results showed that altitude was a dominant environmental control in the encroachment of low-altitude plants, with the degree of encroachment gradually weakening with altitude.This was confirmed in a field investigation of D. angustifolia [29].Our finding that low altitudes were generally more susceptible to plant encroachment than higher altitudes was in general agreement with findings reported in other mountainous areas [36,37].Becker et al. (2005) observed that the numbers of alien species declined steeply with altitude at 232 sites in the Swiss Alps [38].Ropars and Boudreau (2012) discovered that shrub expansion was more intense on low altitude than on exposed hilltops at the forest-tundra ecotone of Quebec, Canada [39].High alpine zones with harsh environmental conditions (e.g., cold temperatures, strong winds, poor soil nutrients) may be resistant to the encroachment by plant species with low environmental tolerances [38,40,41].
In previous studies, the local aspect was considered to be an important factor affecting vegetation changes [13,19].Gottfried et al. (1999) found that the upward shift of native plants was closely related to their abilities to find suitable micro-topographical situations [15].However, we found that the mountain aspect had an overwhelming effect on vegetation changes.In alpine ecosystems, the mountain aspect strongly affects the intensity of several abiotic factors (e.g., wind, solar radiation, precipitation), resulting in significant differences in soil characteristics (soil temperature, water, and nitrogen availability) on different mountain aspects [16,17,42].Therefore, the composition and dynamics of mountainous vegetation (e.g., treeline) are directly linked to the mountain aspect [42][43][44][45].Our results showed that plant encroachment on the south-and west-facing aspects of the mountain was more severe than on other mountain aspects.This was likely because the southand west-facing mountain aspects receive greater precipitation from the dominant southwestern winds, and consequently have high soil moisture and nitrogen deposition [46].Zong et al. (2014) also reported heavy encroachment by herbaceous plants on the western mountain aspects [29].Although micro-topography also provoked changes in vegetation, our study showed that the mountain-scale effects of topographic controls should be worthy of more attention in future research on mountain ecology.
Previous studies have reported the effect of topographic complexity on mountainous vegetation changes [19,47].They showed that strong environmental heterogeneity induced by complex terrain could provide diverse habitats for plants to respond to climate change.Scherrer and Koerner (2010) found that spatial variation in microhabitat temperatures, mediated through complex terrain, provides refuge for plants under warming climate conditions [23].Refuges were the sites that could offer favorable environmental features, decoupled from external climate conditions [21,48].Dobrowski (2011) found that refuges, topographic depressions (e.g., valleys, sinks) with relatively cold temperatures and high moisture content, were widely distributed in areas with complex terrain [21].In addition to providing refuge for native plants, encroaching plants also utilize sites with better habitats [49,50].Thus, areas with better environmental conditions are vulnerable to encroachment by highly competitive nonnative plants [51].Our study confirmed the above findings by showing that severe encroachment occurred in topographically complex areas by competitive species, such as D. angustifolia.The encroaching plants were mostly distributed in low-lying areas, such as valleys [52].
Our study area experienced a destructive volcanic eruption about 1000 years ago, and still retains a large amount of volcanic ash [53].In general, areas with severe volcanic erosion require a longer pioneer phase than areas with weaker erosion [54,55].For the mountainous areas affected by the volcano, volcanic ash was the primary abiotic factor constraining recovery of vegetation [56].Del Moral and Wood (1988) found that rates of vegetation recovery on tephra plots were relatively slower than at other sites on Mount St. Helens [57].We found an expansion of endemic plants near the volcanic ash fields in the alpine tundra.This may be because these areas have experienced a longer period of vegetative recovery than the areas located far from volcanic ash fields, and endemic plants did not appear until recent decades.A study based on image interpretation also reported the expansion of the tundra plants in the uppermost area of the mountain (where volcanic ash and rocks were mainly distributed) [35].However, most of the alpine tundra is currently rarely affected by volcanic erosion, and the vegetation changes are mostly attributed to environmental changes [4][5][6].Therefore, volcanic erosion was not a common factor affecting alpine vegetation changes.

Conclusions
We found significant changes in vegetation in the alpine tundra of the Changbai Mountains using Landsat time series data at the senescence periods.The changes did not occur in a random fashion, suggesting that environmental heterogeneity may have affected the change of vegetation.Our results showed that different topographic variables affected encroachment by low-altitude plants and expansion of endemic plants.In our study, the mountain aspect was a significant factor affecting vegetation changes.We anticipate that future environmental conditions will exacerbate the effects of the mountain aspect on changes in alpine plant communities.In alpine systems, topographic factors at the mountain scale may play important roles in vegetation changes, and their effects deserve attention in future research of changes in alpine vegetation.Finally, we found that topographic complexity also strongly affects vegetation changes.The suitable habitats offered within a complex terrain facilitate inhabitation and competition among species, especially under climate warming.Therefore, vegetation changes in mountainous regions are predicted to occur in the areas with high topographic complexity under future climate conditions.

Figure 1 .
Figure 1.Location of the study area, CMNR representing the Changbai Mountains Nature Reserve.

Figure 1 .
Figure 1.Location of the study area, CMNR representing the Changbai Mountains Nature Reserve.

Figure 2 .
Figure 2. Flow diagram illustrating quantification of topographic controls on vegetation changes.CVE represents the coefficient of variation in elevation, DVA represents the distance to the nearest volcanic ash deposit.DEM = digital elevation model.NDVI = normalized difference vegetation index.TM/ETM+ = Thematic Mapper/Enhanced Thematic Mapper Plus.

Figure 2 .
Figure 2. Flow diagram illustrating quantification of topographic controls on vegetation changes.CVE represents the coefficient of variation in elevation, DVA represents the distance to the nearest volcanic ash deposit.DEM = digital elevation model.NDVI = normalized difference vegetation index.TM/ETM+ = Thematic Mapper/Enhanced Thematic Mapper Plus.

Figure 3 .
Figure 3. (a) GF-2 image of the study area during the senescence period (the green areas represent Rhododendron aureum in Subfigure ①, the white areas represent volcanic ashes in Subfigure ②), (b) the distributions of R. aureum and volcanic ashes.

Figure 4 .
Figure 4. Reflectance curves of R. aureum and Deyeuxia angustifolia at the senescence periods.

Figure 3 .
Figure 3. (a) GF-2 image of the study area during the senescence period (the green areas represent Rhododendron aureum in Subfigure 1 , the white areas represent volcanic ashes in Subfigure 2 ), (b) the distributions of R. aureum and volcanic ashes.

Figure 3 .
Figure 3. (a) GF-2 image of the study area during the senescence period (the green areas represent Rhododendron aureum in Subfigure ①, the white areas represent volcanic ashes in Subfigure ②), (b) the distributions of R. aureum and volcanic ashes.

Figure 4 .
Figure 4. Reflectance curves of R. aureum and Deyeuxia angustifolia at the senescence periods.

Figure 4 .
Figure 4. Reflectance curves of R. aureum and Deyeuxia angustifolia at the senescence periods.

Figure 6 .
Figure 6.Vegetation changes as reflected by NDVI values during the 1980s, 1990s, 2000s, and 2010s, respectively.The NDVI value represents the average of the Landsat NDVI for the relevant decade.NDVI < 0.0 represents volcanic ashes, 0.0-0.2represents non-evergreen plants, and >0.2 represents R. aureum.

Figure 6 .
Figure 6.Vegetation changes as reflected by NDVI values during the 1980s, 1990s, 2000s, and 2010s, respectively.The NDVI value represents the average of the Landsat NDVI for the relevant decade.NDVI < 0.0 represents volcanic ashes, 0.0-0.2represents non-evergreen plants, and >0.2 represents R. aureum.

Figure 6 .
Figure 6.Vegetation changes as reflected by NDVI values during the 1980s, 1990s, 2000s, and 2010s, respectively.The NDVI value represents the average of the Landsat NDVI for the relevant decade.NDVI < 0.0 represents volcanic ashes, 0.0-0.2represents non-evergreen plants, and >0.2 represents R. aureum.There were two directional changes in tundra vegetation.In one, the plant community shifted from one dominated by R. aureum to increasing dominance of non-evergreen plants (mostly due to encroachment of low-altitude plants), reflected by a gradual decline in NDVI values.The other change was the expansion of tundra endemic plants, such as R. aureum, reflected by an increase in NDVI values.Based on the rate of change between 1988 and 2017, and levels of statistical significance, we classified NDVI change into three types: Increased, decreased, and unchanged.About 30% of the area showed statistically significant changes in NDVI values (p < 0.05), with NDVI values decreasing (encroachment of low-altitude plants) in 21.3% of areas and increasing (expansion of endemic plants) in 8.7% of areas (Figure7).The remaining area (70%) had no changes.The encroachment of low-altitude plants occurred mainly in the low-altitude areas on the south and west sides of the study area.The expansion of endemic plants largely occurred near the volcanic ash areas.The areas of unchanged vegetation were concentrated on the east side of the mountain (Figure7).

Figure 7 .
Figure 7. Spatial distribution of NDVI change from 1988 to 2017.

Table 1 .
List of satellite images.

Table 1 .
List of satellite images.

Table 2 .
Topographic variables and ecological meanings.

Table 3 .
The results of binary logistic regression.

Table 3 .
The results of binary logistic regression.