Vegetation disturbances characterization in the Tibetan Plateau from 1986 to 2018 using Landsat time series and field observations

Disturbances in vegetated land could dramatically affect the process of vegetation growth and reshape the land cover state. The overall greenup of vegetation on the Tibetan Plateau (TP) has almost served as a consensus to date. However, we still lack consistent acquisitions on the timing, the spatial patterns, and the temporal frequency of vegetation disturbance over the TP, limiting the capacity for planning land management strategies. Therefore, we explored the spatiotemporal pattern and variation of vegetation disturbances across the TP during the past decades and analyzed the disturbance agents. We utilized 37-year Landsat time series images and field observations coupled with a temporal segmentation approach to characterize the spatiotemporal pattern of vegetation disturbances across the TP for the period 1986–2018. The results from this study revealed that 75.71 M ha (accounting for 29.34% of TP’s area) vegetation area underwent at least one disturbance, of which 8.44 M ha area ever experienced large-scale disturbances (disturbance area greater than 0.9 ha and disturbance magnitude (the difference between the spectral value of pre-disturbance and that of post-disturbance) over 0.2). Further, the spatial distributions of these large-scale disturbances varied over time: before 2002, the disturbed sites were evenly distributed over the southeast part of the TP probably induced by overgrazing and unscientific livestock management, while after 2002, most disturbances were concentrated in the south of the Yarlung Tsangpo, mainly caused by anthropogenic activities, such as urban area, roadways, railway, and water control projects. This study presents an effort to characterize vegetation disturbances and their variations over the past decades on the TP, which provides crucial insights toward a complete understanding of vegetation dynamics and its causal relationship with human activities.


Introduction
As an integral part of the terrestrial ecosystem, vegetation regulates the exchange of carbon, water, impetus, and energy between the earth's surface and the atmosphere (Bonan et al 1992, Haberl et al 2007. Change in vegetation structure and function could infer the effects of climate change and human interventions because vegetation is highly sensitive to various natural and anthropogenic factors (Parmesan and Yohe 2003). To tackle urgent issues related to global climate change, carbon emissions and losses of biodiversity, detailed information about long-term vegetation dynamics over a large spatial extent and a long temporal scale is urgently needed (Dixon et al 1994, DeFries et al 1999, Brun et al 2019, Chen et al 2021, Smith et al 2022. In the past several decades, considerable attention has been dedicated to the changes in vegetation at regional to global scales using forest inventory and satellite observation data (Hansen et al 2013, Ceccherini et al 2020, Piao et al 2020. However, most of these only focused on the long-lasting trend of vegetation growth, lacking more detailed depiction of the inter-annual and intra-annual changes in vegetation, especially for the disturbance change within them (Verbesselt et al 2010, Watts andLaffan 2014). Referring to other related studies (Verbesselt et al 2010, Zhu et al 2020, we define the vegetation disturbance as the abrupt change where there is a significant decrease in vegetation volume within vegetated land caused by anthropogenic or natural processes. Vegetation disturbances caused by logging, urbanization, and fire could potentially alter the process of vegetation growth and affect local ecosystems. Accurate detection of such changes and attribution to the causes are therefore crucial for ecosystem management, biodiversity protection, and conservation policies (Krawchuk et al 2020, Ye et al 2021b. The Tibetan Plateau (TP), referred to as the Earth's 'Third Pole' , is one of the most vulnerable areas under global climate change. Over the past decades, a series of studies have investigated and revealed the trend of vegetation greening in the TP (e.g. Shen et al 2015a, Zhu et al 2016. Climate change, particularly the warming (Qiu 2008) and wetting (Chen et al 2013) trends, are thought to driving the long-term greening tendency of the vegetation on the TP (Piao et al 2014, Zhu et al 2016. Nevertheless, detailed information on inter-annual vegetation fluctuation was often ignored from the general greening trend in these previous studies (Piao et al 2011). The disturbance change of vegetation is a significant source of the inter-annual fluctuations especially for the ecological fragile zones like the TP. In recent years, the TP has suffered from more and more vegetation disturbance such as fire, overstocking (Harris 2010), active rodents (Yang et al 2017), increasing attention has therefore been paid on this regard. However, most of the previous studies were developed at a local scale or based on field observations (Qin et al 2014). Few studies have attempted to comprehensively characterize the disturbance-driven changes in vegetation over the TP with regards to their timing, duration, frequencies, and intensity.
Over the past decades, remote sensing technology has been highly sought after and found wide applications in terrestrial observation. Specifically, free and open access to the Landsat archive from 2008 has unlocked the critical bottleneck to the cost and access to imagery, which has changed the way towards exploring the planet (Banskota et al 2014, Wulder et al 2019). Moreover, the dense time series image stacks have paved the way for a seamless detection of vegetation disturbances at an annual time step (Hansen et al 2013). Accordingly, various tools and algorithms have been proposed for monitoring vegetation change based on the Landsat time series, including offline change detection (e.g. Huang et al 2010, Kennedy et al 2010, online change detection (e.g. Verbesselt et al 2010, Zhu and, and ensemble strategy (e.g. Healey et al 2018, Bullock et al 2020. However, most related studies are focused on forest disturbances detection (Hansen et al 2016, Shimizu et al 2017, Senf and Seidl 2020, Ye et al 2021a, with fewer investigations conducted for vegetation disturbances, especially that occurring in mountainous areas. We investigate the spatiotemporal characterization of vegetation disturbance of the TP from 1986-2018 for improving knowledge of the impact of human and natural coupled disturbances on the TP's vegetation. Specifically, the aims of this study were to: (a) detect the spatial-temporal characteristics of multiple vegetation disturbance and the inter-annual change of disturbed area for the study period; (b) analyze disturbance-driven vegetation changes among different vegetation types; (c) investigate the potential agents of vegetation disturbance change.

Materials and methods
The overall workflow of study is delivered in figure S1 (SI appendix), which consists of imagery preprocessing, index selection, time-series segmentation, and disturbance mapping.

Image data sets and preprocessing
A time series surface reflectance collection of Tier 1 Landsat 4, 5, 7, and 8 images during the period from 1984 to 2020 were employed in the study  for the actual detected years because no disturbances were mapped for the years 1984-1985 and 2019-2020, as current approach is unable to detect the disturbances in the beginning and end of the period). We only selected images within the plant growing season (from early May through middle September, SI appendix, figure S2) to minimize the influence of ice and snow and guarantee inter-annual spectral comparability of surface independent of phenological stages (Shen et al 2015b). Pixels with clouds and cloud shadows were masked out from all images utilizing the Fmask algorithm (Zhu and Woodcock 2014). To assure the inter-sensing harmonization and temporal continuity between Landsat sensors, we also conducted spectral transformation between Operational Land Imager and Enhanced Thematic Mapper Plus (ETM+) using the coefficients suggested by Roy et al (2016). Images were merged into annual composites using the medoid composites method, which is robust against extreme values and more efficient at generating value that is representative of the time series than the commonly used maximum composite method (Flood 2013). The preprocessing of the images was achieved with the help of Google Earth Engine.
To obtain the distribution characteristics of vegetation disturbance on the TP, we also used a vegetation map (Hou 2019). According to the vegetation map, we have divided five vegetation categories, including forest, steppe, meadow, shrub, and alpine vegetation. The remaining vegetation types were labeled as 'others' (figure 1).

Index selection
An analysis of six spectral indices was carried out in order to determine which one captures the variation in time series most accurately: the normalized burn ratio (NBR), the normalized difference vegetation index (NDVI), normalized difference moisture index (NDMI), and the tasseled cap transformation indices greenness (TCG), brightness (TCB) and wetness (TCW). The analysis consists of two parts, first, we chose a disturbed area covered with different vegetation (forest, meadow, and shrub) to evaluate the sensitivity of the spectral indices at pre-disturbance, during-disturbance, and post-disturbance based on the annual composites derived from the preprocessing (Runge et al 2022), and it is presented in figure 1 with a size of 0.15 M ha; second, we calculated the disturbance signal-to-noise ratio (DSNR) at pixel-level on the selected area for each spectral index, which is a fit metric used to evaluate spectral index's effectiveness in disturbance detection (Cohen et al 2018). Higher DSNR values indicate higher signal corresponding to noise: where y 1 and y n respectively represent the values of the start and end vertices of the segment; y i represents the fitted spectral value and y i represents the original spectral value for all fitted points between the start vertex and end vertex. The boxplot was used to displayed different discriminating power (pre-disturbance, during-disturbance, and post-disturbance) of indices and the distribution of DSNR values for different indices. Specially, our algorithm and analysis only focused on 'vegetated regions' , which hereinafter are defined as the pixels with an average NDVI larger than 0.1 over all growing seasons.

Time-series segmentation process
The pre-processed time series data were applied as input to the temporal segmentation algorithm, LandTrendr (Kennedy et al 2010). It is a spectraltemporal segmentation algorithm which is used to detect change in a time series of moderate resolution satellite imagery for producing trajectory-based spectral time series information. This generated a set of annual fitted images and a series of vertices bounding line segments (SI appendix, The process of LandTrendr temporal segmentation). According to the segmented trajectories, we further utilized several disturbance attributes (figure 2): timing of disturbance is defined as the first year in which the disturbance is visible; magnitude of disturbance is defined as the difference between the spectral value corresponding to the onset vertex and that corresponding to the ending vertex; duration of disturbance is defined as the time span (in years) of the segment .

Disturbance change mapping
In this study, we applied the standard LandTrendr parameters for segmentation and fitting processes. Because of the concerns towards the forest management unit, much forest-related research set thresholds for both the mapping unit and disturbance magnitude (e.g. Griffiths et al 2012. However, in order to maximize sensitivity in detecting vegetation change, we did not limit the mapping unit and disturbance magnitude for the mapping of vegetation disturbance times in the first part of the analysis (Saura 2002). Furthermore, we defined the disturbances with area greater than 0.9 ha and magnitude over 0.2 as large-scale disturbance referring to forest-related studies (Griffiths et al 2012, Hislop et al 2019. To present the characterization of the largescale vegetation disturbances, we picked the maximum change segment from each pixel's LandTrendr trajectory and implemented a set of spatial filters, including setting a minimum mapping unit of 0.9 ha (10 pixels) and setting a minimum magnitude of disturbance of 0.2. We found that there was an obvious difference in the distribution of early and late stages, and therefore divided the analysis into three categories, including the results of the whole period (i.e. 1986-2018), the results before 2002 (midpoint of the whole period), and the results after 2002. From these segments that LandTrendr detected, our study resulted in a set of metrics which depict the timing, magnitude, and duration of the changes in vegetation.

Accuracy assessment 2.5.1. Sampling procedure
To validate the vegetation disturbance map, we built a stratified random sample set, including 918 pixels for undisturbed vegetation and 432 pixels for disturbance class (SI appendix, figure S3), which was determined following the sampling design suggested by Olofsson et al (2014). We also conducted a field survey towards 100 points (68 for undisturbed area and 32 for disturbed areas) in 2020. All samples were manually interpreted by independent interpreters through TimeSync 7 , combined with very high resolution (VHR) imageries accessible in Google Earth and some other satellite-borne VHR imageries, like QuickBird (0.61 m), WorldView (0.3-0.5 m), and RapidEye (5 m) (Cohen et al 2010). Since some samples cannot be labeled because of data gaps, a total of 1353 samples (496 samples for disturbance and 857 disturbance samples for no disturbance) remained for accuracy assessment.

Accuracy analysis
The accuracies of the final disturbance map were assessed using 1353 pixel-based interpreted samples. The accuracies were evaluated with regard to overall accuracy (OA) according to the sample count-based error matrix following the method suggested by Olofsson et al (2014). The class-specific omission error, commission error and F1 score were also calculated based on the error matrix. To evaluate the accuracy of detected disturbance onset, we compared the detected disturbance years (i.e. the disturbance onset assigned by the algorithm) versus the interpreted disturbance years (i.e. the disturbance onset assigned by the manual interpretation) in the validation data. In addition, we inspected two temporal segmentations for two examples to perform qualitative analysis.

Analysis of the potential agents for the vegetation disturbance change
Considering the spatial distribution of the large-scale disturbances after 2002, we used the points from field survey to explore the potential agents for the large-scale vegetation disturbance and further selected five typical areas to inspect the disturbance agents. The agent was divided into three categories, including anthropogenic factor (e.g. urbanization and road construction), natural factor (e.g. fire and earthquake), and others. To comprehensively examine the historical change and change agents, we employed multiple data resources for disturbance characterization, including the statistical yearbook, local history raw Landsat data, and other available satellite images (i.e. QuickBird, WorldView, and RapidEye).

The best performed index for disturbance mapping
The results suggest that the NBR, NDVI, and NDMI had the highest discriminating power between the three stages in their respective time series (figure 3). The NDVI and NDMI exhibited the most marked differences between the three stages, while the changes in the NBR was relatively minor. Besides, the DSNR distribution presents that the NDVI have a higher distribution of DSNR compared to that of the NDMI, showing a better performance for detecting vegetation disturbance. Although the indices as NBR are broadly regarded as better at describing forest structure than the indices using shorter wavelengths (e.g. the ever-present NDVI) (Hislop et al 2018), the NDVI is less sensitive than other indices to varying sun-sensor geometry, hence making a NDVI a safe and robust option for time-series analyses (Veraverbeke et al 2012, Morton et al 2014; in addition, most vegetation types in the TP belong to the species with relatively low biomass. In the case of low vegetation biomass, NDVI performs better than those shortwave-infrared based indices (e.g. NBR and NDMI) which are sensitive to variations in the forest structure (Cohen and Goward 2004, Cuevas-gonzález et al 2009, Gasparri et al 2010), by avoiding common saturation problem. The NDVI was therefore selected as the input for temporal segmentation.

The accuracy of disturbance mapping
After taking the potential bias derived from the stratified sampling design into consideration, our disturbance mapping obtained an OA of 81.15% (table 1). Disturbance commission and omission errors were 21.10% and 33.67%, respectively. As for the undisturbed class, the commission and omission errors were 17.84% and 10.27%, respectively. The F1 score of accuracy assessment reached 0.72, which indicated that the LandTrendr performed well in capturing the vegetation disturbances. Figure 4 displayed that most samples are distributed on or near the 1:1-line, suggesting that the onsets of disturbances have been correctly assigned. The qualitative analysis demonstrated that the breaks and the trends identified from temporal segmentation align well with what we observed from the ground (SI appendix, figure S4).

Disturbance maps and its spatial and temporal characteristics over the TP
The disturbance distribution and frequency are mapped in figure 5. We identified a total of 75.71 M ha of vegetation disturbance across the TP over the period of 1986-2018, which accounts for 29.34% of TP's area. In all vegetated areas, most have not been disturbed (count = 0, 69.6%, figure 5(E)). Within the disturbed vegetation areas, vegetation disturbance occurred mainly for 1-3 times, and they were evenly distributed in the TP. Also, the areas with a single disturbance were dominated, accounting for a proportion of 22.8% of TP's area. The disturbed area was widely distributed in the vegetated area while the distribution was relatively dense in the following areas: the areas around the cities, like Sining (figure 5(A)), Lhasa, and Nyingchi (figure 5(B)); the southern edges of the TP where the land was covered with dense forests (figure 5(C)) were also suffered from a higher disturbance frequency; the eastern edges of the TP with high forest coverage which was once devastated by a major earthquake in 2008 ( figure 5(D)). Figure 6 presents the area of vegetation disturbed by year, and an average of 2.29 M ha of vegetation has experienced disturbance per year. Apart from the year 1988 with a relatively high disturbed areas detected (3.53 × 10 6 ha), the results for other years did not vary greatly. Still, an obvious increasing trend was observed from 2015 to 2018.
We further explored the characteristics of timing, magnitude, and duration of vegetation disturbance for the first, second, and third times within the study period, respectively (figure 7). As for the first time (figures 7(A), (D) and (G)), the timings of disturbances were relatively scattered throughout the vegetated areas over the TP. Specifically, the first disturbances in the central of the TP were earlier than that  in the northern and southern regions of the TP; the magnitude of disturbances was mainly concentrated between 0 and 0.10, and the disturbances distributed in the southeastern TP was larger than that distributed in the northwestern TP; as for the disturbance duration, most disturbances lasted 1-5 yr distributed across the vegetated area. As for the second time (figures 7(B), (E) and (H)), the disturbance in the central and east of the TP was earlier than that in the northern and southern edges of the TP; the results of the magnitude of the second disturbances were similar with that of the first disturbance; most second disturbances lasted 1-5 yr also while some disturbance lasted more than 16 yr, which was mainly distributed in the central and eastern TP. As for the third time (figures 7(C), (F) and (I)), the spatial pattern of    disturbance was similar to that for the first and second times, but the disturbance pixels were much sparser. Figure 8 displays the spatial and temporal patterns of timing, magnitude, and duration of the largescale vegetation disturbances. The timing of the disturbance of different categories (i.e. 1986-2018, 1986-2001, and 2002-2018) was obviously different between the south and north of the east-west section of the Yarlung Tsangpo river (SI appendix, figure S5) in the southeastern of the TP (figures 8(A)-(C)). Most disturbances in the north occurred early, while most in the south occurred late. Furthermore, the disturbance points were evenly distributed in the southeastern of the TP before 2002 while that were In summary, our results indicate that there was a noticeable difference between the situations before and after 2002 in terms of the large-scale disturbances. Specifically, the disturbances before 2002 were evenly dispersed across the southeastern TP, while the disturbances after 2002 were spatially clustered in the south of the Yarlung Tsangpo. Figure 9 presented the results in two parts including all disturbances (disturbance for the first, second, and third time, respectively) and the large-scale disturbances (the results of 1986-2018, 1986-2001, and 2002-2018, respectively). The disturbance distribution over different land covers are similar for the disturbances of different sequence (the first, second, and third time) (figure 9(A)), dominated by Meadow and Steppe (together over 50% of the disturbed area) and followed by Shrub, Forest, Alpine vegetation, and Others.

The distribution of vegetation disturbances among different vegetation types
As shown in figure 9(B), during the whole study period , most large-scale disturbances happened over meadow (40.09%), followed by the shrub (25.79%) and the forest (20.80%); during the period of 1986-2002, large-scale disturbances dominated over meadow (52.06%) with substantial proportions over shrubland (19.83%) and forest (13.62%); during the period of 2001-2018, on the contrary, the disturbances occurred mostly over shrubland (33.06%) and forest (29.56%), but less over meadow (25.49%).

The potential agents for the vegetation disturbance
Among all points, 56.3% of the points induced by anthropogenic factor, 25% by natural factor, and 18.7% by others. We further chose five typical areas to inspect the agents (figure 10). Point A was within the affected region of the Lhasa-Pondo Water Control Project. The detected disturbances in this area coincide with the beginning year of the project. The detailed disturbance map shows that the timing of the main disturbances was 2014 and the outline of disturbance points is similar to the shape of the river in 2014. Point B is located between the G219 National highway and the Lhasa-Nyingchi Railway, and most disturbances here occurred in 2019. We also found that the road and the railway were under intensive construction around that year. Point C is located in the southwest of Nyingchi city, and the disturbances appeared between 2015-2018. The results suggest substantial land use change from vegetation to other types during the period. Point D is located in the northwestern Nyingchi near the G318 National highway. The main disturbances occurred in 2008, which is attributed to a forest fire. Point E is located at the eastern fringe of the TP belonging to Wenchuan County, Sichuan province. A large area of vegetation was disturbed in 2008, which is reasonable as an earthquake of magnitude 8.0 named the Wenchuan Earthquake occurred in this area.

Discussion
For the overall disturbance, the timings of the first and second disturbances in the central of the TP were earlier than that in the northern and southern regions of the TP. The range of these earlier disturbances is also similar to that of rangelands of the TP (Dong and Sherman 2015). Some previous studies have proved that the rangelands of the central TP have been overexploited by local pastoralists (Dong et al 2011). Therefore, these earlier disturbances may be because large areas of grazed grassland (steppe and meadow) in the central of the TP, and the vegetation have been disturbed due to the overgrazing in the early years.
For the large-scale disturbances, it is intriguing that the disturbance locations were evenly distributed in the southeastern of the TP before the 21st century, while the disturbance points were concentrated in the south of the Yarlung Tsangpo and were sparse in the north of the Yarlung Tsangpo after the 21st century.
The uniform and wide distribution of the disturbance areas over meadow before 21st century ( figure 8(B)) was probably caused by overgrazing. Many previous studies have reported that the TP's grasslands were under degradation in recent decades (Wang et al 2005, Dong et al 2013, which might be largely attributed to a combination of overstocking beyond carrying capacity and unscientific livestock management (Li 1994, Wang 1999. After entering the 21st century, a range of ecological restoration projects (e.g. The Grain for Green Project, Grazing Withdrawal Program, and others) have been enforced across the TP (Li et al 2016). These projects have shown their effectiveness in facilitating the restoration of vegetation cover, coincident with the diminished frequency of the disturbances in the meadow found in this study (figure 9(B)) and the sparse distribution of disturbances in the north of the Yarlung Tsangpo after 2002 ( figure 8(C)).
The distribution of the large-scale vegetation disturbances among different vegetation types indicated that the disturbances occurred more frequently in the meadow in the early period, while in the late period, there is a decreased disturbance in the meadow areas with an increased disturbance observed in forest and shrub areas. The meadow is an important part of the rangeland on the TP, accounting for 44.64% of the area (Jian 2002). Therefore, a possible explanation is that the local economy was highly dependent on grazing before the 21st century; however, after the 21st century, local residents became less dependent on livestock husbandry and started to interfere with the forest and shrub more frequently due to the increasing demand for wood in the context of urbanization. For instance, the Tibet's industrial structure has been changing after entering the 21st century: the ratio of the primary industry has decreased gradually while an increase in the ratio of the secondary industry has been observed, showing a pattern with tertiary in the lead, secondary industry as the middle, and primary industry as the smallest one; along with the development of the economy, many labors from the primary industry have transferred to the other industries (Tibet Bureau of Statistics 2020). We further analyzed five typical areas to inspect the disturbance agents. Main disturbances in five areas were caused by anthropogenic factors (i.e. the construction of hydro power station, the construction of highway and railway, urbanization) and natural hazards (i.e. forest fire and earthquake). Moreover, most disturbance sites are found in the triangle area formed by the Lhasa, Shigatse, and Nyingchi, which is the most economically active region on the TP in the 21st century. In general, economic development is always accompanied by intensive human activities, such as highway or railway construction, urban expansion, and water conservancy project. Roads are one of humankind's most productive linear infrastructures, driving habitat transformation directly (Trombulak and Frissell 2000). Especially, the timings of disturbances in point B were within the period of the construction of the G219 National highway and the Lhasa-Nyingchi Railway, and the vegetation was therefore more likely be disturbed by these projects. Li et al (2013) also demonstrated the areas along the traffic arteries that are exposed to frequent human activities had suffered severe degradation on the TP. Many related studies revealed that the impoundment of the reservoir has a great influence on the growth of riparian vegetation (Li et al 2015, Botelho et al 2017. Also, hydropower dams may also fragment habitat, reduce the distribution of primary vegetation, and simplify the types of vegetation along the rivers (Li et al 2012). These factors may explain the occurrence of concentrative disturbances in point A. As for the urbanization, although it plays a central role in economic growth and social development, it also profoundly alters the physical environment far beyond city range, leading to habitat loss, vegetation degradation (Chen 2007). Over the last two decades, Tibet has been witnessing a dramatic growth of economy like many other cities in China. As table 2 shows, there has been a dramatic rise in the number of several significant indicators of economic development, including hydroelectric power, length of highways, length of rail lines, and urbanization rate. We also found that there is an increasing trend with the area of disturbed vegetation each year from 2015 to 2018 when these indicators showed an apparent increase. These further proved that the growth of vegetation has been affected by anthropogenic activities to varying degrees along with rapid economic growth.
Although a high-resolution map of vegetation disturbance was obtained, there are still some limitations regarding the data and method. Firstly, cloud cover was an inherent issue over the TP, especially in the southeastern TP (figure 3). In contrast to the conterminous United States, the data density of Landsat archive in the TP is so far substantially lower (Wulder et al 2016). These factors considerably limit the acquisition of cloud-free Landsat observations, which further increases the difficulties for extracting accurate disturbance information from time-series datasets. Multisensor analysis may potentially alleviate data density issues but could also introduce uncertainties during the data fusion process . Secondly, we were not able to distinguish the disturbance agents for all individual disturbance patch, because there were not adequate agent-specific training data. Future work is warranted to develop disturbance agent algorithms and map products, which will provide a more comprehensive examination of disturbance drivers for the TP region. Despite these limitations, this work provides important baseline disturbance information, which laid a solid foundation for further investigating the relationship between vegetation growth and anthropogenic activities over the last decades across the TP in future.

Conclusion
In this study, the spatiotemporal patterns of vegetation disturbance across the TP were investigated based on the LandTrendr algorithm using 30 m Landsat archives for the period of 1986-2018. The results demonstrated that the disturbances in vegetated land have been existing in the context of vegetation greening across the TP over the last decades. During 1986During -2018.71 M ha (accounting for 29.34% of TP's area) of vegetated area was disturbed, of which 8.44 M ha area experienced large-scale disturbances. The disturbances were dominated in meadow areas, largely attributed to a combination of overstocking and livestock mismanagement before 2002. Along with the implementation of several policies regarding ecology conservation around 2002, the disturbances in the meadow have diminished. Instead, most vegetation disturbances were gathered in the south of the Yarlung Tsangpo after 2002 which was mainly induced by anthropogenic activities such as urbanization, the construction of roadways and railway, fire hazards, and water control projects. These emerging threats may bring about unexpected changes in TP's ecology in the future. Our work provides critical insights into vegetation variation over the TP during the past 40 yr, which are essential for development planning and environmental management.

Data availability statement
The data that support the findings of this study are available upon reasonable request from the authors.