The Spatio-Temporal Variation of Vegetation and Its Driving Factors During the Recent 20 Years in Beijing

1. Introduction

Terrestrial vegetation is an important part of the earth’s ecosystem. It releases oxygen through photosynthesis, which is an important source of oxygen on earth. Moreover, terrestrial vegetation also has a positive impact and feedback on the carbon cycle, nitrogen cycle and water cycle (DAVID 1995; Hughes et al. 2001; Tong et al. 2016). There are many factors affecting vegetation variations, among which climate changes and human activities are the two most important categories (Jiang et al. 2018; Sun et al. 2015; Wang et al. 2015; Zhao et al. 2018). For one thing, it is generally believed that climate is an important factor in vegetation change. Rainfall and temperature are the main climate variables as claimed in many previous studies (PEI et al. 2012; Wang et al. 2015; Yu et al. 2020). For another thing, with the growth of the economy and population, more and more studies demonstrate that human activities are the key factors causing vegetation changes (Feng et al. 2015; Lü et al. 2015; Wang et al. 2016). One of the most representative human activities is the rapid expansion of urbanization. The population in urban areas has increased rapidly from 14% in 1900 to 50% in 2000 of the total population, in particular, it has reached 70% in developed countries (Kinsella 2001; McIntyre 2000). From regional to global scales, the urbanization process and related suburban development not only have a significant impact on our living environment but also inevitably bring changes to the vegetation ecological environment.

As an effective method to replace time-consuming and laborious field experiments, remote sensing technology has the advantages of fast data acquisition and a wide detection range, which facilitates the use of a wide range of large-area vegetation monitoring (Hunt et al. 2004). Many researchers use the Normalized Difference Vegetation Index (NDVI) calculated by the near-infrared and red bands to monitor the growth and coverage of vegetation (Lamchin et al. 2018), which is helpful to analyze the problems in the ecological environment. Sun et al. (Sun et al. 2013) and Huang et al. (Huang et al. 2016) indicated that climate factors have an important impact on vegetation change with high spatial-temporal heterogeneity. Precipitation has a positive effect on vegetation, while temperature has a negative effect based on the NDVI in the Qinghai Tibet Plateau. Many scholars have used NDVI to conduct in-depth research on how climate (Pettorelli et al. 2012; Xu n.d.), topography (Corenblit and Steiger 2009), land use and human activities (Chen et al. 2019; Xu et al. 2018) affect vegetation growth.

Beijing is located in the transition zone between the mountain and the plain. The plain area is surrounded by mountains on three sides, forming a natural arc barrier, resulting in two distinct climate in front of and behind the mountains. In the past decades, the population and economy have developed rapidly. And the process of urbanization has been accelerated. With the rapid urban expansion in Beijing plain area, the mountainous area of Beijing is ecological conservation, an important biodiversity center, a water source protection area, as well as the main natural vegetation distribution area. According to the ecological situation and the current needs of urban development, the government has formulated several development plans of "Five-Year" to ensure that the eco-environment is fully protected and improved while Beijing’s economy is developing rapidly. Therefore, it becomes more and more important to explore the relative role of climate variables and human activities in vegetation changes and the correlation between urbanization expansion and vegetation growth trends in Beijing.

Under the influence of global climate change, Beijing has also experienced rapid economic development and is facing complex ecological problems. In recent years, many works are studying the patterns and driving factors of vegetation variations in Beijing. Zhao Y et al. (Zhao et al. 2019) found that the vegetation in Beijing had an overall increasing trend and the NDVI fluctuations in several particular years were greatly related to temperature or precipitation anomalies based on geographically weighted regression and ordinary least squares during 2000-2015. Jiang M et al. (Jiang et al. 2017b) calculated the fractional vegetation coverage by the method of dimidiate pixel model based on NDVI, which suggests that human activities are very significant factors to influence and explain the changes in Beijing and they are highly spatially heterogeneous from 2000 to 2015. Chang Y et al. (Chang et al. 2020b) considered that the response of vegetation to urbanization showed obvious differences and geographical heterogeneity in the urbanization gradient based on the nighttime light data. Nevertheless, these works do not continue into the latest year and quantify the role of human activities and climate factors in the process of vegetation restoration and degradation in Beijing. Moreover, there are few studies on integrating impervious surface data and NDVI data to explore how vegetation changes in rapid urbanization.

For a fast-developing city, quantitative and continuous research about the influence of climate variability and rapid urbanization on vegetation change is of great significance. Based on the NDVI dataset, meteorological data, and impervious surface data in recent 20 years, we aim to: (1) investigate the spatiotemporal patterns of vegetation change trends by using Theil Sen and Mann Kendall methods excluding the influence of water; (2) evaluate the dominant factors of affecting vegetation change in the mountains by using residual trend analysis approach; (3) explore the trend of vegetation growth in different stages of urban development. The results help to understand the current vegetation driving mechanism and provide a scientific basis for formulating reasonable vegetation construction, land use, and ecological environment protection strategies in Beijing.

2. Materials and Methods

2.1. Study Area

Surrounded by Taihang Mountains and Yanshan Mountains, Beijing (115.7°-117.4°E, 39.4°-41.6°N) lies on the north of North China Plain (Figure 1). The total area of Beijing is 16410.54 $k{m}^2$, of which the plain area is 6200 $k{m}^2$, accounting for 38%, and the mountainous area is 10200 $k{m}^2$, accounting for 62%. The average altitude of the whole city is 43.5 meters with a plain elevation of 20-60m and a mountain elevation of 1000-1500m generally. The city has a typical North Temperate semi-humid continental monsoon climate which is hot and rainy in summer, cold and dry in winter, and short in spring and autumn. The annual average precipitation is about 450mm with 80% of the annual precipitation occurring in June, July, and August. The temperature in July is the highest, and that in January is the lowest. Although it is cold and dry in winter, there is more sunshine, with the average sunshine being more than 6 hours per day. Beijing is divided into 16 districts(counties), including high-intensive building area, low-intensive building area and mountain area. The area within Fourth-Ring Road is referred to as the high-intensive building area. The area with slope greater than 8 degrees is the mountain area. And the remaining area is the low-intensive building area.

2.2. Data

2.2.1. MODIS-NDVI Dataset and Preprocessing

We used Moderate-resolution Imaging Spectroradiometer (MODIS) Terra 16-day vegetation index product–250m NDVI (MOD13Q1) in the manuscript from Land Processes Distributed Active Archive Center. The product has minimized the effect of cloud, cloud shadows, and noise through 16 days of synthetic data. Three MODIS tiles (h26v04, h26v05, and h27v05), 240 images in total, are used in this study from 2000 to 2019.

The three tiles are spliced together, and the NDVI data of Beijing is cut out using the Beijing boundary map as a mask. According to the previous research (Xueying 2010), the vegetation in Beijing reached the peak growth stage in July and August. The Maximum Value Composites (MVC) (Wang and Zhu 2019) was utilized to process the data. The maximum value of the data was taken in July and August respectively. The mean value of July and August is calculated to represent the vegetation growth in that year and further used in the subsequent time series analysis.

2.2.2. Landsat Imagery and Preprocessing

Water body has a great influence on NDVI time series, even sub-pixel water body has a great influence on NDVI. In order to avoid the influence of water body as much as possible, the water body is extracted by the Landsat images in Beijing. Beijing is covered by two (123/032 and 123/033) WRS2 paths/rows of Landsat images. A total of 717 Landsat surface reflectance images from 2000 to 2019 are collected in this study (http://earthexplorer.usgs.gov/). All the data details are listed in Table 1. CFMask algorithm was used to detect and remove snow, cloud, and cloud shadow for Landsat Image (Foga et al. 2017; Zhu et al. 2015; Zhu and Woodcock 2012). NDVI (Tucker 1979) and Enhanced Vegetation Index (EVI) (Huete et al. 2002, 1997) respond rapidly to canopy structure changes and are sensitive to high biomass areas (Huete et al. 2002). Modified Normalized Differences Water Index (MNDWI) is widely used to evaluate surface water information (Xu 2006).

According to (Zou et al. 2017, 2018), if the conditions satisfied ((MNDWI > NDVI or MNDWI > EVI) and (EVI < 0.1)), the observation results could be determined as open surface water. The interannual frequency of open surface water was calculated for each pixel and classified the pixels as open surface water if the interannual frequency was greater than or equal to 0.75 (Zou et al. 2017, 2018). Therefore, we generate 20 water body maps from 2000 to 2019 and select the map with the largest open surface water area as the water mask. The nearest neighbor method is implemented to resample to 250m in order to match MODIS vegetation index data.

2.2.3. Climate Data

The climate data from 2000 to 2019 come from the National Earth System Science Data Center, National Science & Technology Infrastructure of China (http://www.geodata.cn). Based on the global 0.5 ° climate data released by Climatic Research Unit (CRU) and the global high-resolution climate data released by WorldClim, the dataset was downscaled in China through Delta spatial downscaling scheme (Peng et al. 2019). The data of 496 independent meteorological observation points were used for verification, and the verification results are credible. To match MODIS data, these datasets were bilinear interpolated from 1km to 250m resolution.

2.2.4. Impervious Surface Data

The long-term impervious surface data is used to quantify city expansion and investigate the response of vegetation change trends to urbanization in this study. Impervious surface is an important part of the urban foundation surface, revealing the changes in the conurbation area better (Shahtahmassebi et al. 2012). Based on Landsat images, Gong et al. (Gong et al. 2020) developed impervious surface data of Chinese provinces with 30m resolution spatially, covering the period from 1985 to 2018 (Impervious surface data download website: http://data.ess.tsinghua.edu.cn). The impervious surface data of Beijing from 2000 through 2018 were selected to keep consistent time with the MODIS-NDVI product in this study. We aggregated the impervious surface data into 250 m resolution to match MODIS data spatially and calculated the impervious surface ratio (ISR) in each pixel. The impervious surface ratio (ISR) of conurbation demonstrates the complicacy of the city eco-environment and the extent of urban development (Arnold Jr and Gibbons 1996; Ridd 1995).

2.2.5. SRTM DEM Data

InSAR technology was applied in the SRTM mission to obtain elevation information of surface features for its efficiency and economical characteristics. The SRTM data used the WGS84 coordinate system as the horizontal reference and the EGM96 geoid as the vertical reference. By far it’s the best global digital elevation data comprehensively evaluated in terms of coverage, elevation accuracy, and publicity. Since 2003, SRTM3 with a ground resolution of 90 meters has been released. We also resample it to 250m, calculate the slope, and select the slope threshold value as 8 to distinguish mountainous areas (Niu et al. 2009).

2.3. Methods

As shown in Figure 2, the study first calculates the slope and Student’s t test value of the NDVI variable for each pixel during 2000-2019 to obtain change trends of vegetation by Theil-Sen and Mann-Kendall method. In the mountain area, the meteorological data is utilized to analyze the influence of climate factors on vegetation variations by correlation analysis method. And residual analysis is used to distinguish the impact of climate factors and human activities on vegetation. In high and low-intensive building areas, the impervious surface data is leveraged to study the response of vegetation to urbanization by correlation analysis.

2.3.1. Variation Trend Judgement

In this study, the Theil-Sen slope is acquired to get the trend of greenness and meteorological variables (temperature and precipitation). Compared to the least square linear fitting method, TS trend analysis is not affected by noise and outliers in the data (Fensholt et al. 2012; Sen 1968; Theil 1950), and it’s a widespread method in time series trend analysis of climate variables (Martinez et al. 2012; Tabari et al. 2011). The Theil-Sen slope is estimated based on the median value of the observations $X_j$ and $X_i$ at all pairwise time steps j and i. The formula is as follows:

$$slope=\mathrm{Median}\left(\frac{x_j-x_i}{j-i}\right)$$

(1)

where $x_i$ and $x_j$ represent the value of NDVI, annual accumulated precipitation or annual average temperature in the i-th and j-th year respectively. For NDVI variation, $slope$ >0 means greenness restoration and $slope$ <0 means greenness degradation. For meteorological variables, $slope$ >0 represents precipitation or temperature increase , and $slope$ <0 represents precipitation or temperature decrease.

As a common method, we calculate the value of the Mann-Kendall method to test significance because the method has no requirement on the distribution (Kendall 1938). The inspection process is as follows:

For the sequence, the test statistic S is:

$$S=\sum _{i=1}^{n-1}\sum _{j=i+1}^{n}sgn(x_j-x_i)$$

(2)

$$sgn(x_j-x_i)=\left\{\begin{array}{c}+1,x_j-x_i>0\hfill \\ 0,x_j-x_i=0\hfill \\ -1,x_j-x_i<0\hfill \end{array}\right.$$

(3)

where n is the length of the sequence, $x_i$ and $x_j$ are the observations at i and j time respectively. When the data value has independent identically distribution ,S is the normally distributed approximately, and the variance is given by the following formula (Kendall 1948; Mann 1945):

$$VAR\left(S\right)=\left(n(n-1)(2n+5)\right)={\sigma }^2$$

(4)

where $\sigma$ is the standard deviation.The significance of the test method is given by the statistical value Z:

$$Z=\left\{\begin{array}{c}\frac{S-1}{\sqrt{VAR\left(S\right)}},S>0\hfill \\ 0,S=0\hfill \\ \frac{S+1}{\sqrt{VAR\left(S\right)}},S<0\hfill \end{array}\right.$$

(5)

where $\left|Z\right|\ge 1.96$ (equivalent to $p\le 0.05$) is considered as of significance.

According to the vegetation greenness changes and the significance analysis, vegetation change results are classified into four groups: significant restoration, restoration, degradation and significant degradation. Mann-Kendall tau ($\tau$) coefficient is adopted to investigate the correlations between two random variables as known as Kendall rank correlation coefficient.

2.3.2. Residual Analysis

Residual analysis is utilized to separate human-induced vegetation trends from the vegetation changes caused by climate through "NDVI-Climate" model. Considering the correlations between NDVI and meteorological variables, the "NDVI-Climate" model is established by multiple linear regression in every pixel (Qu et al. 2018).

$$\mathrm{NDVI}\left(n\right)=a\times T\left(n\right)+b\times P\left(n\right)+C$$

(6)

where n is year; C is a constant; a and b are correlation coefficients of annual average temperature (T) and accumulated precipitation (P) of that year respectively. a, b and C are calculated by the least square method.

The NDVI of climate-driven vegetation change is predicted by the "NDVI-Climate" model. And the residual represents the NDVI of human-driven vegetation trend (Evans and Geerken 2004; Ibrahim et al. 2015; Jiang et al. 2017a).

$${\mathrm{NDVI}}_{\mathrm{climate}}(m,n)=a\times T(m,n)+b\times P(m,n)$$

(7)

$${\mathrm{NDVI}}_{\mathrm{human}}(m,n)=\mathrm{NDVI}(m,n)-{\mathrm{NDVI}}_{\mathrm{climate}}(m,n)$$

(8)

where ${\mathrm{NDVI}}_{\mathrm{climate}}$ is the predicted value of the model in a specific year with temperature and precipitation impacts, and ${\mathrm{NDVI}}_{\mathrm{human}}$ is the value of anthropogenic NDVI variation during 2000-2019.

According to Theil-Sen method, we calculate the slopes of ${\mathrm{NDVI}}_{\mathrm{human}}$ and ${\mathrm{NDVI}}_{\mathrm{climate}}$ as $S_{\mathrm{human}}$ and $S_{\mathrm{climate}}$, respectively. We also calculate the slope of actual NDVI as $S_{\mathrm{NDVI}}$. Therefore, the six types of impact factors on vegetation dynamics were determined (Table 2). Specifically, $S_{\mathrm{NDVI}}$ >0 indicates vegetation increases and $S_{\mathrm{NDVI}}$ <0 indicates vegetation decreases. $S_{\mathrm{climate}}$ >0 represents that climate factors have a positive impact on vegetation, while $S_{\mathrm{climate}}$ <0 represents that climate factors have a negative impact on vegetation; $S_{\mathrm{human}}$ >0 represents that human activities have a positive impact on vegetation, and $S_{\mathrm{human}}$ <0 represents that human activities have a negative impact on vegetation (Wang et al. 2019a).

2.3.3. Correlation Analysis

Considering the spatial distribution of vegetation change trend results, Mann-Kendall coefficient correlation analysis is adopted to study the correlation between NDVI and ISR in high and low-intensive building areas and the correlation between NDVI and climate variables in the mountain area. The $\tau$ coefficient of the Mann-Kendall method is also called the Kendall rank correlation coefficient, which is used to test the statistical correlation of the observation values of two random variables, especially for small sample size.

3. Results

3.1. Spatial-temporal Characteristics of Vegetation Change Trend

The change trends of vegetation in Beijing during four periods are shown in Figure 3 including 2000-2010,2000-2013,2000-2016,2000-2019. The corresponding statistical results of vegetation trends are shown in Figure 4. It can be seen that the vegetation in Beijing is totally improved. The vegetation increase is mainly concentrated in mountainous areas and high-intensive building areas. Vegetation degradation, especially the significantly degraded vegetation, is mainly concentrated in low-intensive urban areas. The significant restoration area changes from 27.8% in 2010 to 52.8% in 2019, and the restoration area changes from 55.5% in 2010 to 27.4% in 2019, which indicate that a large area of vegetation restoration has changed into significant restoration in recent years. The increase of vegetation in some urban forest parks (Figure 3 1 2 3) shows that the construction of urban parks plays a positive role in vegetation growth. At the same time,the significant degradation area changes from 2.4% in 2010 to 5.6% in 2019,and the degraded area changes from 13% in 2010 to 12.9% in 2019, which is because that the vegetation growth situation in low-intensive urban areas of Beijing is not optimistic.

We counted the vegetation changes in high-intensive building areas, mountainous areas and low-intensive building areas respectively from 2000 to 2010, 2000 to 2013, 2000 to 2016 and 2000 to 2019. In the high-intensive building area, significant restoration is increasing continuously with 57.4% in 2010, 67.3% in 2013, 84.5% in 2016, and 90.2% in 2019 (Figure 5). The mountainous area is similar to the high-intensive building area, where significant restoration is increasing continuously with 34.2% in 2010, 57.9% in 2013, 64.8% in 2016, and 79.5% in 2019 (Figure 6). Unlike mountainous areas and high-intensive built-up areas, there has been a gradual increase in the degradation of vegetation in low-intensive built-up areas of Beijing (Figure 7), especially significantly degraded vegetation.

3.2. Driving Factors of Vegetation Change in Mountainous Area

The precipitation and temperature data are used to obtain the driving factors of vegetation change in mountain area. The trend of precipitation and temperature change in Beijing is shown in Figure 8. Temperature trends show a more pronounced east-west variation. Precipitation shows an upward trend, but with clear north-south variability. The regions with significant consistency between NDVI and temperature (p<0.05) accounted for 2.8%, while the regions with significant consistency between NDVI and precipitation (p<0.05) accounted for 23.6%, indicating that vegetation in Beijing responded more strongly to precipitation than to temperature. The correlation between NDVI and climate change in mountain area is shown in Figure 9. Precipitation shows a positive correlation with vegetation in most of mountain area(97.7%), while temperature shows a negative correlation in most of mountain area, reaching 78.8%. Temperature and precipitation showed significant differences in response to vegetation, especially in the northern mountainous areas.

The trend analysis of NDVI residuals is shown in Figure 10(a). The NDVI residuals in the northern mountainous areas of Beijing show an obvious upward trend, while those in the southern mountainous areas and areas closer to human activities show a downward trend, which indicates that the trend of NDVI residuals has obvious regional characteristics. The TS slope results of NDVI, ${\mathrm{NDVI}}_{\mathrm{climate}}$ and NDVI residuals are superimposed to obtain the NDVI dominant factor results, as shown in Figure 10(b), and the NDVI driving factor statistical results are shown in Table 3. The results show that vegetation restoration accounts for the majority of the areas in the non-impervious surface area, reaching 92.9%, of which 41.5% is caused by human factors, 45.3% is caused by a combination of climatic and human activity factors, and rest is caused by climatic factors. The areas with significant increases of NDVI residuals and vegetation are mainly located in the northern mountainous areas, where human activities dominate the vegetation growth. The increasing trend of residuals in the southern mountainous area is not significant, but the increasing trend of vegetation is significant. Human activities and climate factors jointly dominate the vegetation change in this area.

3.3. Driving Factors of Vegetation Change in High and Low-Intensive Building Area

Impervious surface data from 2000 to 2018 is utilized to illustrate the spatial-temporal characteristics of urbanization expansion in Beijing over this period (Figure 11(a)). The ISA in Beijing has increased by 138.2% from 1501.6km² in 2000 to 3577.4km² in 2018 (Figure 11(b)). In light of interannual dynamic changes of the newly increased ISA shown in Figure 11(a) and Figure 11(b), the ISA increases rapidly during 2000-2016. In the early period of 2000-2016, the expansion of impervious surface was mainly distributed in peri-urban areas centered on high-intensive urban areas. And in the later period, it was mainly distributed in the distant suburban areas of Beijing, covering more rural areas. After 2016, the expansion of impervious surfaces tended to be almost stable.

In the impervious surface area, the correlation results between vegetation and impervious surface are shown in Figure 12. The positive correlation indicates the growth of vegetation restoration in the process of urban expansion, while the negative correlation indicates the growth of vegetation degradation. The proportion of negative correlation increased each year from 37.5% during 2000-2010 to 53.7% during 2000-2016 (Table 4) due to the continuous urban expansion of Beijing, while the positive correlation continued to decrease, indicating a gradual increase in vegetation degradation area caused by urban expansion during this period. After 2016, the urban expansion tends to be stable, the negative correlation decreases, and the vegetation restores. The spatial distribution of correlation results shows an obvious spatial heterogeneity in vegetation growth in different regions. The negative correlation is mainly concentrated in low-intensive building areas and tends to gradually expand outwards. The positive correlation is concentrated in high-intensive building areas of the city, where the region remains stable, but the positive correlation is gradually enhanced.

4. Discussion

4.1. The Relative Role of Driving Forces in the Process of Vegetation Restoration

Above data analysis indicates a trend of vegetation increasing in Beijing. And the restoration area accounts for 80.2%. The driving factors show spatial variability. In mountainous areas, rainfall is the main climatic variable affecting vegetation growth and change. The northernmost mountainous area of Beijing is the Beijing-Tianjin Sandstorm Source Control Project (BTSSCP), where vegetation restoration is predominantly caused by human activities. This shows that the effect of the BTSSCP project is remarkable. Under adverse climate conditions, human activities play a positive role and significantly improve the vegetation growth situation here, which is also consistent with previous research (Xu et al. 2011; Zhang et al. 2016). The NDVI is positively correlated with precipitation in the south of the mountain area, and the joint determinations of precipitation and human activities dominate the vegetation increases in the south of Beijing. The policy of closing hillsides for afforestation and returning farmland to forest has been implemented in the southern mountainous area. Meanwhile, some protected areas have been established. The construction of these protected areas and the implementation of the policy, together with the positive effect of precipitation, are conducive to vegetation restoration and ecological construction.

In addition to the mountainous areas, vegetation restoration is also concentrated in the high-intensive impervious surface area in Beijing. In the advanced stage of urbanization, the speed of urban expansion slows down. Due to the demand for greening, the vegetation growth trend presents a phenomenon of recovery (Liu et al. 2014; Luck et al. 2009). This phenomenon can also be observed in developed countries such as Japan, the United Kingdom, and America (Liu et al. 2015). The greening work in the area is getting better and better in the regions, which benefited from the urban park construction policy that Beijing has long persisted in (Chang et al. 2020b; Sun et al. 2012). In a word, active human activities play an important role in vegetation restoration (Huang et al. 2016; Li et al. 2015).

4.2. The Relative Role of Driving Forces in the Process of Vegetation Degradation

Vegetation degradation in Beijing is mainly distributed in the low-intensive buildings of Beijing. The area is in the early stage of urbanization, where urban expansion has a negative role in vegetation growth. The ISA has expanded rapidly in the low-intensive building area, where the rapid transformation of cities from non-urban areas to urban areas leads to vegetation degradation (Chang et al. 2020a). In developing countries such as Bangladesh (DewanY 2009), and South Africa (Burgoyne et al. 2016), this feature of urbanization to vegetation dynamics is common. Under the large-scale urban expansion, the ecological environment of vegetation is very fragile (Shan et al. 2015; Wang et al. 2019b) with significant degradation. Therefore, vegetation degradation in Beijing is mainly caused by human activities (Anzhou et al. 2017; Tian et al. 2015).

The negative correlation between vegetation and the impervious surface has been spreading outwards, which is consistent with the trend of urbanization expansion, indicating that human factors are more important in explaining the trend of vegetation change (Li et al. 2016; Lü et al. 2015; Tian et al. 2014). The development of economic and social production may be conducive to ecological restoration and vegetation growth (Lü et al. 2015; Madu 2009; Salvati and Zitti 2008). However, these areas lack ecological protection, restoration projects and policies related to vegetation protection. Under the effects of adversely climatic conditions and unreasonable human activities, vegetation degeneration has occurred in the region, which is opposite to the situation in the core urban areas. Beijing’s urban expansion tends to be stable after 2016 and gradually attaches importance to the protection of suburban vegetation, such as the implementation of a large-scale plain afforestation project (2012-2015). Moreover, given the fragile ecological environment of vegetation in these areas, the government should comprehensively consider the socio-economic development and ecological protection, and formulate more scientific and reasonable policies.

5. Conclusions

In this manuscript, we utilize MODIS NDVI data to analyze the pattern of vegetation growth trend temporally and spatially in Beijing. By leveraging climate data, Landsat images, DEM data, and impervious surface data, we evaluate the relative effects of human activities and climate changes on vegetation restoration and degradation in Beijing. Since 2000, vegetation increase has accounted for 80.2% of the total area in Beijing, mainly located in the high-intensive building and mountainous areas of Beijing. Human activities play a positive role in vegetation growth. Moreover, some large-scale ecological restoration projects have achieved good results. Beijing’s core area is located in the advanced stage of urbanization, and vegetation growth is promoted by urban greening projects in this area. The vegetation degradation area in Beijing accounts for 18.5%, mainly distributed in the low-intensive building areas of Beijing. The combined effects of these two factors and climate change dominate the degradation of vegetation. The current vegetation ecological environment is very fragile. Beijing suburbs are still in the early stage of urbanization, and a large number of land converted into impervious surface which led to the reduction of vegetation. Therefore, we suggest that policymakers should formulate different ecological restoration policies for different ecologically vulnerable areas and areas with different urbanization gradients, and establish sustainable development strategies to ensure coordinated economic and ecological development.

Our research has a few limitations. MODIS13Q1 vegetation index product is mainly used in this study. Although the product is a 16-day synthetic data, it is also found that the observed value is high or low due to cloud interference. Meanwhile, we focus on the impact of urbanization on vegetation growth, which is a representative factor of human activities. In future research, we can combine the vegetation index data obtained by other types of sensors and field observation methods to further optimize NDVI. Human activity factors can examine more data, such as human footprint data, night light data, and so on, to evaluate the extent of man-made vegetation change more scientifically and comprehensively.