# An Object- and Shapelet-Based Method for Mapping Planted Forest Dynamics from Landsat Time Series
## 2 Materials and Methods
### Methods
#### 2.2.1 Data and Processing  
Free-access Landsat data were used as the primary data source in this study because of their long-term observations and relatively high spatial resolution [35; 36]. The Google Earth Engine (GEE) is a cloud computational platform that enables users to access and process remote sensing data for planetary-scale geospatial analysis ([[https://earthengine.google.org/](https://earthengine.google.org/)]([https://earthengine.google.org/](https://earthengine.google.org/)) (accessed on 16 April 2022)). To composite 30 cloud-free NHAR images for each year from 1991 to 2020, we collected Tier 1 surface reflectance satellite image products that had been atmospherically and topographically corrected from the Landsat-5 Thematic Mapper (TM; 1982-2013), Landsat-7 Enhanced Thematic Mapper (ETM+; 1999-present), and Landsat-8 Operational Land Imager (OLI; 2013-present) on the GEE platform. All images were selected during the forests' growing season (from June to September) to reduce the impact of phenological factors [37]. An overview of all the Landsat data used in this study is shown in Figure 3. We then masked clouds, cloud shadows, and snow pixels from all the available images using CFMASK [38]. Finally, annual image composites were generated using a medoid approach, which is the value for a given band that is numerically closest to the median of all corresponding pixels among the images considered (all images within a provided annual data range) for a given image pixel. Time-series layer stacks of the NDVI images were subsequently processed from the collected annual Landsat composite images. The auxiliary data comprised the Shuttle Radar Topography Mission (SRTM) digital elevation model (DEM) from the GEE, and high-resolution images in Google Earth (Google Inc., Mountainview, CA, USA).  
The forest map in 2020 served as a mask to remove non-forest pixels from the annual Landsat-NDVI time series. Here, we used the RF classifier to extract the forest (including natural and planted forests) distribution of the NHAR. Sample datasets are composed of 379 forest samples (including 247 planted forest samples) acquired via field study and 1259 points marked manually that rely on high-resolution images from Google Earth. In total, there are 613 forest samples. In addition to the spectral signatures (i.e., original spectral bands and spectral indices) of the Landsat composite images, textural features calculated using the gray-level co-occurrence matrix (GLCM) [39] and terrain characteristics are involved in the classification. To acquire features with higher sensitivity in the classification of land-cover types, the top 15 variables (Table 1) were selected according to their importance in the initial classification. Consequently, the overall accuracy (OA) of the final classification result was 89.1%. The forest area, including trees and shrubs, is 0.83 million ha, accounting for 16.0% of the total area of the study region; this is generally coincident with the statistical data that indicate the forest coverage of the NHAR in 2020 was 15.8%, as reported by the Ningxia Statistical Yearbook ([[http://tj.NHAR.gov.cn/](http://tj.NHAR.gov.cn/)]([http://tj.NHAR.gov.cn/](http://tj.NHAR.gov.cn/)) (accessed on 10 May 2022)).  
Figure 3: Number of Landsat-5 TM, Landsat-7 ETM+, and Landsat-8 OLI scenes analyzed in this study, available annually from 1991 to 2020.