Spatiotemporal Analysis of Vegetation Changes Along the Belt and Road Initiative Region From 1982 to 2015

As a key ecological zone of terrestrial ecosystem, the Belt and Road Initiative (BRI) region has experienced a significant change of vegetation coverage in recent years. Using the GIMMS NDVI3g, this study investigated the patterns of spatiotemporal variation of vegetation coverage in the BRI region during the period 1982-2015. The Theil-Sen Median trend analysis and Mann-Kendall method were used to analyze the data, followed by the calculation of Hurst index in order to analyze future trends of vegetation coverage. In addition, possible environmental factors affecting this variation were identified by using the partial correlation analysis and residual analysis methods. The results of the study showed that (1) the normalized difference vegetation index (NDVI) of the study area revealed a slow decrease during 1982-2015, with the linear tendency being −0.1%/10a. During this period, a stable increase was detected before 1997 (linear tendency 1.4%/10a), followed by a sharp decline after 1997 (linear tendency −1.8%/10a). (2) In spatial, the areas with increased vegetation NDVI are mainly distributed in Europe, India and China, whereas the regions with decreased vegetation NDVI are mainly distributed in northern Russia, Central Asia, Southeast Asia, the Malay Islands, and northeast China, of which, the magnitude of decrease in the north of Russia is particularly remarkable. This phenomenon indicates that vegetation activities in high latitude regions declined, such as coniferous forest of subfrigid zone and tundra vegetation. (3) The same characteristic of vegetation coverage change were stronger than the reverse characteristics. A total of 89.2% of the study area’s vegetation will change in the same way as in the past, with a continuously increasing area accounting for 34.0% and a continuously decreasing area accounting for 21.6%. (4) Although climate change may play a key role in vegetation growth trends on a long-term scale, human activities are also an important factor driving vegetation change in the BRI area, especially for areas with increased vegetation coverage such as China, India, and Europe.


I. INTRODUCTION
The Belt and Road Initiative (BRI) was proposed and promoted by China's National Development and Reform Commission, Ministry of Foreign Affairs, Ministry of The associate editor coordinating the review of this manuscript and approving it for publication was Xujie Li .
Commerce of China (2015). Under the dual effects of climate change and human activities, terrestrial ecosystems have been unprecedentedly disturbed and threatened. As a sensitive area of climate change and a fragile area of the ecological environment, BRI area is bound to change or destroy its ecosystem structure and function under the background of global change. Therefore, it is necessary to conduct a comprehensive VOLUME 8, 2020 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ survey and scientific evaluation of the various geographical elements in the region to raise the level of awareness of the impact of climate change on typical terrestrial ecosystems.
As the main body of terrestrial ecosystem, vegetation is not only the recipient of climate change, but also has a feedback effect on climate change. It is a priority index for ecological environment monitoring, ecological risk and vulnerability assessment. Therefore, it is very important and urgent to comprehensively evaluate the characteristics of spatial and temporal changes of vegetation coverage and their attribution in the BRI area.
As the topic of global climate change is gaining more and more attention from academic circles, the study of the relationship between global climate change and terrestrial ecosystems is a scientific issue that is currently of great concern to the international community. Among them, the response of vegetation to climate change has undoubtedly become the focus and core issue of scholars around the world. A lot of research work on the Tibetan Plateau [3], northeastern Siberia [4], northern Scandinavia [5], the Yangtze River [6], [7], Loess Plateau [8], [9], Mongolian [10], the African Sahel [11], Central Asia [12], North America [13], temperate and boreal Eurasia [14], semi-arid areas [15], the Northern Hemisphere [16], high latitudes area [17] and global scale vegetation changes [18]- [20]. The above research results have revealed the changes of vegetation cover and its relationship with climate factors from different scales, and have drawn meaningful conclusions. However, these studies mostly focus on the trend and fluctuation characteristics of vegetation and climate change, and the research on the impact of extreme climate events and human activities on vegetation growth is slightly weaker, especially in the sensitive and fragile areas of climate change such as the BRI region. In addition, studies have shown that from the 1980s to the early twentieth century, global warming caused the start of the growing season (SOS) of temperate vegetation in the northern hemisphere to advance and the end of the growing season (EOS) to be delayed [16], and are in line with the increase in net primary production suggested by modeling [21]. This response of vegetation to climate change has aroused widespread concern in the academic community. Beck and Goetz [22] and Jiang et al. [12] analyzed the trend of vegetation change in North America and central Asia, respectively, and found that global warming increased evaporation and water consumption of vegetation, resulting in decreasing boreal forest productivity in North America and significantly decreased for shrubs and sparse vegetation in central Asia. Another study found that the correlation between vegetation and temperature in the northern hemisphere is weakening, and it is believed that the cause may be due to drought [23]. This shows that extreme climate events not only affect the vegetation ecosystem itself, but also change the relationship between vegetation and climate factors. To sum up, it can be seen that the above research results have deepened our understanding of the changes in vegetation coverage. It is worth noting that there are significant differences in the patterns and processes of geographical elements in different regions. Strengthening research on local areas is conducive to a comprehensive understanding of the evolution characteristics of terrestrial ecosystems. In the context of climate change, vegetation is extremely sensitive to extreme climates, and it is necessary to re-understand and evaluate the changes in vegetation cover and its feedback relationship with climate.
With strong surface heterogeneity and complex ecological structure, the BRI area is extremely sensitive to the response to global climate change, which has already caused a significant impact on the vegetation coverage in this area. At the same time, with the proposal of the ''Green Silk Road'', it is essential to fully understand the basic characteristics of vegetation coverage and its response to climate change in BRI region. Based on the above understanding, under the guidance of the landscape pattern & ecological processes theory [24], using vegetation index data, supplemented by Theil-Sen Median analysis, Mann-Kendall method, Hurst exponent, and residual analysis method, this paper analyzes the variations characteristics and future trend of vegetation coverage changes in BRI area, moreover the driving factors of vegetation changes were studied. The purpose is to deeply understand the characteristics of regional ecological environment changes and promote the coordinated and sustainable development of regional ecological environment and social economy.

A. STUDY AREA
The BRI area is located between 12 • E-180 • E to 11 • S-82 • N, mainly involving the central and eastern Eurasia and northern Africa. As a platform for international cooperation, the GDP of BRI area is about 23.56 trillion US dollars, accounting for 31.79% of the global total, and the population reaches 4.58 billion, accounting for 62.61% (World Bank, 2015). During 1982 to 2015, the annual average maximum temperature was 31.7 • C, which appeared in the western side of the Arabian peninsula ( Figure 1D), and the minimum temperature appeared in the northeast Russia, which reached −22.4 • C ( Figure 1B). The maximum total annual precipitation is about 10534.7 mm, which occurs in the Bengal ( Figure 1C), and the smallest precipitation occurs in the Egypt, with annual precipitation about 3.8 mm ( Figure 1A). The vegetation in this area is diverse, including almost all of vegetation types from evergreen broad-leaved forest to tundra ice and snow.  The NDVI data used for vegetation monitoring in this paper is the latest version of the GIMMS-NDVI3g.v1 dataset (http://glcf.umd.edu/data/gimms/) provided by the National Aeronautics and Space Administration (NASA) global monitoring and modeling group, which is the longest time span NDVI dataset that come from same sensors. Therefore, it can avoid the difference caused by the spectral response function between different sensors. The spatial resolution of GIMMS data is 8 × 8 km, the temporal resolution is half a month, and the time range is from 1982 to 2015. The maximum value combination (MVC) method is used to eliminate the effects of clouds and aerosols, and the DNVI of the current month is obtained by synthesizing two data within the same month. In addition, the data are preprocessed by subset extraction, format conversion, reprojection, and quality control to form a NDVI dataset along the BRI area from 1982 to 2015.
The gridded monthly-mean air temperature (T, • C) and monthly total precipitation (P, mm) data come from the University of Delaware Air Temperature & Precipitation datasets (https://psl.noaa.gov/data/gridded/ data.UDel_AirT_ Precip.html), with the spatial resolution of 0.5 degrees (about 50 × 50 km) and the temporal resolution is months. The gridded fields were estimated from monthly weather-station averages using a combination of spatial interpolation methods: digital elevation model (DEM) assisted interpolation, traditional interpolation, and climatologically aided interpolation (CAI). To indicate (roughly) the spatial interpolation errors, station-by-station cross validation was VOLUME 8, 2020 employed [25]. In order to unify the spatial resolution of temperature, precipitation and NDVI data, the temperature and precipitation data were resampled to 8 × 8 km resolution by bilinear interpolation method.

1) TENDENCY ANALYSIS METHOD
The Theil-Sen Median, as a non-parametric trend analysis method, was used to show the trend of long-term NDVI in the BRI area. The Theil-Sen Median is a robust nonparametric statistical trend identification method with the following formula: where i and j are time series numbers, NDVI i and NDVI j denote the NDVI value at time i and j, respectively; β is the median of the slopes in the n(n−1)/2 combinations. In which, β > 0 indicates that the data shows an increasing trend during this time period, and vice versa. Furthermore, the Mann-Kendall method is used to verify the significance of the trend derived from Theil-Sen Median. The Mann-Kendall method was originally proposed by Herrmann et al. [11], and is only used to detect a change trend of the sequence, and Hamed and Rao [26] have further improved this method so that it can more effectively determine the starting position of various changing trends, and it has a wide detection range, a high degree of quantification and full of vitality [27]. The process of the Mann-Kendall trend test is as follows: where n is the number of data in the series, m the number of tied values, t i the number of ties for the i th value. Z c is the normalized statistic of the Mann-Kendall test, and |Z c | > Z (1−α)/2 indicates that the time series data changes significantly at the α level. The confidence level set in this paper is 0.05, that is, |Z| ≥1.96 indicates that the vegetation data changes significantly, and |Z| <1.96 indicates that the change is not significant.

2) FUTURE TREND ANALYSIS METHOD
Now that we understand the direction and significance of vegetation changes in the past 34 years, this urges us to study the sustainability of this vegetation changes further. The rescaled range analysis (R/S), also known as Hurst exponent, is an effective method to quantitatively describe the long-term dependence of time series information. The Hurst exponent quantifies a sustainable and effective method for representing long-term sequence data, being able to pry into future trends based on past data to a certain extent, and has been used in the fields of hydrology [27] and meteorology [28]. The basic principle is: (1) Considering the NDVI time-series (NDVI (t) ) observed at time t = 1, 2, 3, · · · , n, for any positive integer τ ≥ 1, define the average sequence of the time series as: (2) calculate the cumulative dispersion: (3) calculate the range: (4) calculate the standard deviation: (5) then the Hurst exponent is calculated as: where H is the Hurst exponent and its value range is 0-1 (0 < H < 1). In which, 0 < H < 0.5 indicates that the future trend of the sequence is the opposite of the past. The trend of increasing in the past indicates that the future will decrease, and vice versa, and the closer H is to 0, the more obvious the anti-sustainability; H = 0.5 indicates that the time series information is a random sequence that is independent of each other; 0.5 < H < 1 indicates that the future trend is consistent with the past, and the closer H is to 1, the more obvious this persistence.

3) DRIVING FACTORS ANALYSIS METHOD
In this study, the residual analysis method is used to distinguish the influence of natural factors and human factors in the long-term sequence change of NDVI during 1982-2015.
By establishing the regression model of NDVI, temperature and precipitation, the contribution of temperature and precipitation to NDVI can be predicted grid by grid. Without considering the influence of other non-deterministic factors, the residual between the measured and the simulated value of NDVI is the part that contributed by human activities. This method was first proposed by Evans and Geerken [1] and Geerken and Ilaiwi [2], and is currently widely used in research [29]- [31]. The calculation method is as follows: where NDVI obs is the pixel DNVI value from GIMMS dataset observed by the satellite. NDVI pre is the predicted value of NDVI fitted according to temperature and precipitation, which represents the state of vegetation without human influence. ε is the part contributed by human activities in NDVI. The area where ε > 0 indicates that human activities have a positive effect on promoting vegetation. Conversely, ε < 0 indicates that human activities have an inhibitory effect on vegetation growth, and ε ≈ 0 means that human activities have a weak impact. What needs illustration is that the Interactive Data Language (IDL) is used to process the GIMMS NDVI, the air temperature, and the precipitation data, besides, the ArcGIS software is used to complete the mapping.

A. THE CHARACTERISTICS OF THE TEMPORAL VARIATION OF VEGETATION
From 1982 to 2015, the NDVI of vegetation in the BRI region showed a slow decreasing trend, with a decrease rate of 0.1%/10a, but this trend is not significant. (Figure 2). Through piecewise linear analysis, we found that the vegetation changes in the BRI region in the past 34 years can be divided into two stages. (1) From 1982 to 1997, there was a significant increase trend with the growth rate of 1.4%/10a.
(2) From 1998 to 2015, it turned into a downward trend with a deceleration rate of 1.8%/10a. Research by Park and Sohn [32] also showed that the increasing trend of vegetation cover in Eurasia gradually slowed down, stagnated, and even declined in some areas during the 1980s and 1990s.
(3) In addition, we found that the NDVI of vegetation declined rapidly from 2011 to 2012 in the BRI region. By comparing NDVI of two years, it was concluded that 60.7% of pixels showed a declining trend, which may be caused by the change of atmospheric environment caused by the strong la Nina event in this period, which further affected vegetation cover [33], [34].

B. THE CHARACTERISTICS OF THE SPATIO-TEMPORLA TREND OF VEGETATION
The spatial mean value can represent the overall variation trend of NDVI of vegetation. However, due to the fact that the variation trend of different regions is opposite and then cancels each other, the variation characteristics of different regions cannot be well described. In this study, Theil-Sen Median analysis and the Mann-Kendall test were combined to examine the trend of vegetation cover change in the B&R region from 1982 to 2015 (Figure 3). The results of the Theil-Sen Median analysis showed that, from 1982 to 2015, (1) the regions with increased NDVI were mainly distributed in central and eastern Europe, India and southeast China, accounting for 24.9% of the study area.

C. THE SUSTAINABILITY OF THE VEGETATION VARIATIONS
Although the Theil-Sen Median analysis and the Mann-Kendall test can reveal the variation characteristics of vegetation in the BRI region, the sustainability of this trend is still unknown. Therefore, the Hurst exponent is used in this study to represent the sustainability of the vegetation changes. As shown in Figure 4A, the number of pixels with a Hurst exponent value greater than 0.5 accounted for 89.2%, with the average of 0.72, indicating that the vegetation sustainability in the study area was slightly stronger than the anti-sustainability. The regions with a Hurst exponent below 0.5 are mainly distributed in the Kazakhstan, western part of Malaysia, Arabian Peninsula, southeastern part of the Mongolia Plateau and the Arctic coastal region. The Hurst index is greater than 0.5 in most of the remaining regions, especially in Europe, West Asia, South Asia, Eastern China, and Russia.
As can be seen from Table 1 and Figure 4B, (1) the areas where vegetation coverage sustainability significant and slight improvement cover a larger area, accounting for 34.0% of the total area, mainly distributed in East Asia, South Asia, Asia Minor Peninsula, and Central Europe.
(2) The spatial distribution of the sustainability severe and slight degradation regions is scattered, mainly distributed in the northern Russia, Central Asia (the coast of the Black Sea, Caspian Sea and Aral Sea), Southeast Asia, the Malay Archipelago, and northeast China accounting for 21.6% of the total area. (3) The areas where there is no obvious law of vegetation changes are mainly distributed in Egypt, Arabian Peninsula, Iranian Plateau and western China, accounting for 44.4% of the whole study area.

D. THE RELATIONSHIP BETWEEN VEGETATION CHANGE AND TEMPERATURE AND PRECIPITATION
This paper uses partial correlation analysis to analyze the impact of precipitation and temperature factors on NDVI in the BRI region ( Figure 5). In general, the spatial difference between the relationship between precipitation and temperature and NDVI is obvious. The partial correlation coefficient between the annual average NDVI and temperature in the BRI region is between -0.83 and 0.86, and the partial correlation coefficient with the total precipitation is between -0.79 and 0.86. Vegetation and temperature are closely related in Europe, Asia Minor Peninsula, Pakistan and eastern China. However, the areas with strong correlation  between vegetation and precipitation are mainly distributed in South and Central Asia.

E. THE IMPACT OF HUMAN ACTIVITIES ON VEGETATION COVER
In order to further define the influence of human activities affecting the NDVI in the BRI region, this paper calculated the slop of NDVI residual form 1982 to 2015 ( Figure 6). If the residual value tends to zero, indicating that the vegetation growth in this area has a significant correlation with precipitation and temperature; the larger the residual, the greater the influence of human activities on the vegetation growth in this area. A positive and larger residual value indicates that vegetation growth is disturbed by human activities and tends to improve. On the contrary, a negative and smaller residual value indicates that vegetation growth is disturbed by human interference and presents a deteriorating trend. The regions with significantly increased NDVI residuals are mainly concentrated in Europe, South Asia, and East Asia.
The continuous increase in NDVI residuals indicates that the growth of vegetation in these areas cannot be explained solely by climate change, which largely reflects the impact of human activities, and human activities mainly contribute to the increase of vegetation in these areas.

IV. DISCUSSION
The analysis of the dynamic changes of vegetation and its driving mechanism has always been a core issue of concern in the scientific research of global change and the construction of ecological civilization. The key is to select the appropriate scale and combine multiple data for multi-dimensional analysis. This study explores the change trend of vegetation in the BRI area from three aspects: variation trend, significance, and sustainability, and uses the method of residual analysis to explore the impact of human activities on vegetation change. We found that NDVI in southeastern China, India, and Europe has grown rapidly over the past 34 years, which is consistent with previous studies. For example, by analyzing MODIS and Landsat data, Chen et al. [35] found that more than 65% of the vegetation coverage areas in China and India showed a greening trend, becoming the fastest growing region in the world. On the other hand, studies have pointed out that the increasing trend of vegetation coverage in Eurasia gradually slowed down and stagnated during the 1980s and 1990s, and even declined in some areas [32], [36]. For the entire BRI region, our research found that 1997 was a turning point in the region's vegetation from an upward trend to a downward trend (Figure 2). The main reason for this kind of browning is that the vegetation in the Northern Hemisphere in spring (32.1%) and summer (33.5%) has decreased compared with previous years [37].
Regarding the reduction of NDVI, De Jong et al. [18] found in their research that most coniferous forests in northern Russia showed browning. Our research also found that NDVI in the high latitudes of northern Eurasia has generally decreased over the past 34 years. This may be due to temperature-induced drought, because these areas experienced significant warming without increasing rainfall simultaneously [38]. In addition, we detected that the trend of vegetation degradation in Central Asia is also relatively obvious. Increased temperature and reduced precipitation are the main factor affecting vegetation degradation in the Kyzylkum Desert and the northern Ustyurt Plateau [39]. And the degradation of some shrubs in the southern part of the Karakum Desert, the southern Ustyurt Plateau and the wetland delta of the Large Aral Sea were mainly triggered by human activities: the excessive exploitation of water resources and oil and natural gas extraction [12]. Browning in Indonesia and other parts of Southeast Asia may be related to the expansion of rubber and palm oil plantations at the expense of tropical forests [40]. It has been reported that this transformation has occurred on a large enough scale that browning can be seen [18]. But in general, through residual analysis, we found that the impact of human activities on vegetation in the BRI region is mostly positive ( Figure 6). This positive impact may come from development activities such as land development, intensive agriculture (irrigation, fertilization), plantation and nitrogen deposition [41]. For example, evidence from China shows that the greening of farmland and forest land is mainly caused by agricultural practices and afforestation and reforestation project, such as the Grain for Green project and the Three-North Shelter Forest Program [42]. For India, fertilization and irrigation may have led to greening from 1982 to 2003 [43].
Remote sensing, as the most powerful means of monitoring vegetation activities on a global scale, is widely used in the study of terrestrial ecological changes [44], but there are also certain uncertainties. For example, the saturation effect of NDVI on high-coverage vegetation affects the sensitivity of monitoring vegetation changes [19]. In addition, the spatial resolution of the GIMMS data used in this study is 8 km × 8 km. In this case, it is difficult to ensure the uniformity of the features in the pixels. The problem of mixed pixels has also had a certain impact on this study. Finally, this paper uses precipitation and temperature to estimate the vegetation status under natural conditions, and uses residuals to represent man-made influences. This method is applicable in most cases [23], [45], but the vegetation response to temperature and precipitation has hysteresis and regional differences [46], [47]. Human activities also affect precipitation and temperature, therefore, how to separate the impact of climate change and human activities on the change of vegetation coverage is a problem that needs to be solved in future research.

V. CONCLUSIONS
In this paper, the Theil-Sen Median analysis, Mann-Kendal test, Hurst exponent, and residual analysis method were combined with the latest version of GIMMS NDVI3g 1982-2015 data to analyze the spatial and temporal patterns, trends, sustainability, and the principal driving factors of vegetation along the BRI area. It aims to explain the ecological environment status of the BRI area from the perspective of vegetation, and provide decision support for the ecological environmental protection and construction of the area. This article mainly draws the following conclusions: From 1982 to 2015, the NDVI of vegetation in the BRI region showed a slow decreasing trend, with a decreasing rate of 0.1%/10a, which failed the significance test. Vegetation changes have obvious phase characteristics over the past 34 years. Among them, the vegetation NDVI showed an obvious increasing trend from 1982 to 1997, with an increase rate of 1.4%/10a, and showed a decreasing trend from 1998 to 2015, with a decrease rate of 1.8%/10a. Spatially, vegetation NDVI in Europe, India and China showed a significant increase trend, while vegetation NDVI in key ecological regions and transition zones in Northeast China, Central Asia and northern Russia, Southeast Asia, and the Arabian Peninsula showed a significant downward trend.
The results of Hurst exponent analysis show that the same direction features of vegetation variations in BRI area are greater than the reverse features. A total of 89.2% of the study area's vegetation change is sustainable: NDVI in 34.0% of the study area will continue to increase, and 21.6% of the study area will continue to decrease.
In the BRI region, the change of vegetation NDVI is the result of both climate change and human activities. In general, the impact of human activities on vegetation is positive, and it will promote the increase of vegetation, especially in China, India and Europe.