INTER-ANNUAL DIFFERENCE OF LAND SURFACE PHENOLOGY (LSP) TRANSACTION DATES USING TIME SERIES LANDSAT IMAGES

Land surface phenology (LSP) is a kind of vital information for land cover classification and vegetation growth monitoring. Time series Landsat images, with the advantages of long observations and high spatial resolution, have been widely used in LSP identification. However, LSP transaction dates, such as start of season (SOS) and end of season (EOS), are highly influenced by the coarse temporal resolution. In this study, we compare the inter-annual difference of LSP SOS from 5 years interval, 10 years interval and all years interval Landsat images, and improve the SOS estimated model by considering the accumulated growing degree-days (AGDD) of soil temperature and soil moisture. Results indicate that LSP SOS can serve as a good proxy for reflecting ground vegetation phenology, especially using 5 years interval Landsat images. Soil temperature and soil moisture have certain influence on SOS estimation, and the R-squared value reached 0.9 after model adjustment. This study can provide guidance for estimating suitable inter-annual LSP transaction dates under different sceneries in the future. * Corresponding author


INTRODUCTION
Land surface phenology (LSP) is a kind of vital information for land cover classification and vegetation growth monitoring, reflecting the changes in terrestrial ecosystems and climate. Vegetation indexes (VIs) derived from remote sensing images can potentially serve as a proxy for responding seasonal and annual changes of ground vegetation phenology. Advanced Very High Resolution Radiometer (AVHRR), Landsat series, Moderate Resolution Imaging Spectroradiometer (MODIS), Visible Infrared Imaging Radiometer Suite (VIIRS) and other new remote sensing images have been widely used in the last few decades. For example, Zhang et al., (2018) created one LSP product derived from the VIIRS remote sensing images at a 500-meter gridded spatial resolution to supplement the LSP products derived from MODIS since 2001. Zhang et al., (2020) proposed a new algorithm of LSP product at a 30-meter gridded spatial resolution fusing the operational harmonized Landsat and Sentinel-2 (HLS) products and VIIRS surface reflectance products during 2016 and 2018.
With the increasing demand of the long term and high precision LSP products, Landsat series are becoming more and more important. However, due to the coarse temporal resolution (16-day revisit cycles), yearly fulfil requirement images are often insufficient to support the LSP identification such as start of season (SOS) and end of season (EOS), especially under cloud, rain or other bad weather conditions. To overcome this limitation, several approaches have been development to improve the temporal resolution and keep the high spatial resolution at the same time. Multi-source remote sensing image fusion is a common way to enrich the dataset pool. Spatial and Temporal Adaptive Reflectance Fusion Model (STARFM) (Gao et al., 2006), and Improved Flexible Spatiotemporal Data Fusion (IFSDAF) (Liu et al., 2019) have been successfully applied in the fusion work of coarse spatial resolution (MODIS) and Landsat series. In addition, empirical models such as general liner regression are also useful to fill up the missing observation (Pouliot et al., 2018). However, due to the obvious difference in spectral information among these remote sensing images, and the lack of long term observations, both fusion and empirical models have certain limitations.
Currently, a two-step model is often applied to (1) calculate the average LSP transaction dates by accumulating Landsat images over years, (2) and evaluate inter-annual variability according to the seasonal similarity of vegetation growth (Li et al., 2019). It cannot be ignored that even slight inter-annual environmental changes will have a certain impact on the VIs. Therefore, in order to quantify the difference caused by environmental conditions and acquire high accuracy estimated SOS results, this study firstly collects Landsat images during 1995 and 2016 from 1000 sample sites in New York State under different land types. Then, the inter-annual difference of LSP SOS in 5 years interval, 10 years interval and all years interval Landsat images are compared. Finally, ordinary linear regression (OLR) model is built with environmental factors to improve the estimation accuracy of LSP SOS.
The remainder of this paper is organized as follows. The study area and datasets are detailedly described in Section 2. The methodology and quantification indexes are described in Section 3. Results and discussion are presented in Section 4. Finally, conclusions are listed in Section 5.

Study Area
We randomly sampled 1000 sample sites in New York State. These sample sites were evenly distributed on different latitude and land types. In addition, we collected PhenoCam sites (totally 31 sites, 195 years results) with continuous observations in northeastern United States (Figure 1).

Datasets
We calculated Normalized Difference Vegetation Index (NDVI) value for all available Landsat images at the 1000 sample sites on the Google Earth Engine (GEE) platform (Gorelick et al., 2017). These datasets contains L1T-level data products from the Landsat5 TM (Thematic Mapper, during 1995 to 2011), Landsat7 ETM+ (Enhanced Thematic Mapper Plus, during 1999 to 2016), and Landsat8 OLI (Operational Land Imager, during 2013 to 2016). We used the surface reflectance datasets which have already processed radiation and atmospheric correction (Masek, et al., 2006) and removed clouds and shadows (Zhu and Woodcock, 2012). Due to the substantial difference between Landsat8 OLI images and Landsat7 ETM+ images, we used linear regression model to calibrate the Landsat8 NDVI value (Li, et al., 2014).
PhenoCam server (http://klima.sr.unh.edu/) provides digital images with red, green, and blue (RGB) three colourful channels for every 30 minutes. In this study, we used green chromatic coordinate (GCC) index, which describes the proportion of green channel brightness relative to the total image brightness to derive ground (near-surface) phenological transaction dates.
Temperature and precipitation are considered to be the important environmental factors affecting vegetation phenology (Dannenberg et al., 2014). In this study, we summarized the accumulated growing degree-days (AGDD) information from daily soil temperature (ST) and soil moisture (SM) datasets during January to April using the ERA5 reanalysis products

LSP Transaction Dates
The LSP transaction dates are retrieved from time series Landsat NDVI value using double logistic regression model (DLRM) (Jönsson et al., 2010).
where is the time in the day of year (DOY)， 1 , 2 , 1 , 2 are the parameters of the fitted model. and are corresponding to the minimum NDVI and maximum NDVI value in over years observation.
The start of season (SOS) is defined as the date when NDVI increases to 35% percentage of the difference between the maximum and minimum fitness NDVI value during growing period, and end of season (EOS) is the date where the NDVI decreases to 35% percentage of the difference between the maximum and minimum fitness NDVI value during senescence period. The time period before NDVI reaching the maximum value is defined as growing period, and the rest of time period is defined as senescence period (Figure 2).
Normalizing the LSP transaction dates under different time intervals can be used to quantify the annual variations (Li et al., 2019).
where and are the mean and standard deviation of LSP SOS in 5 years interval , 10 years interval or all years interval.
is the start of season (SOS) for every single year. is the total number of years.
Accumulated growing degree-days (AGDD) is accumulated from daily GDD over the growing season.

Estimated Method
Previous studies have proved that AGDD and NDVI results are in high correlation (Beurs et al., 2004). In our study, we assume that the annual transaction dates would be affected by the variability of yearly soil temperature and soil moisture. Therefore, we utilize environmental factors to improve the estimated accuracy of ordinary linear regression model (OLR).
= 0 + 1 1 + 2 2 +. . . + where is the difference of SOS between LSP SOS and ground SOS, could be LSP SOS, AGDD of soil temperature, soil moisture and other factors.

RESULTS AND DISCUSSION
Affected by the image quality and weather conditions, the average fulfil requirement Landsat images is less than half of the total counts (Table1). In addition, it might lead to obvious errors when calculating the average LSP SOS with time interval less than 5 years.  Figure 3).
RMSE results verified the validation of our assumptions that collecting Landsat images at different time intervals could have impacts on the LSP SOS ( Figure 4). RMSE results in 10-all years interval had the minimum value, followed by RMSE results in 5-10 years interval. RMSE results had the maximum value in 5-all years interval due to the biggest time gaps. RMSE results among different time intervals generally increased in the higher latitude.   The LSP SOS in different time intervals was not significantly different comparing with the ground SOS. Difference between LSP SOS and ground SOS also existed among latitude. Difference of SOS during 41° and 42°, 44° and 45° was smaller than difference during 42° and 43 °, 43° and 44°. However, difference of SOS between 41° and 42° was larger than other latitude with the median value more than 10 days. SOS between 44° and 45° had the smallest different with the median value around 2 days ( Figure 5).
LSP SOS and ground SOS had high correlation. There were 0.759 R-squared value between ground SOS and LSP SOS with 5 years interval, 0.711 R-squared value between ground SOS and LSP SOS with 10 years interval and 0.697 R-squared value between ground SOS and LSP SOS with all years interval. Based on our perspectives, the difference between LSP SOS and ground SOS should have correlations with the AGDD of soil temperature and soil moisture in the previous months before SOS. Through comparison, we found that the AGDD results of soil temperature and soil moisture in the first previous month before SOS had the highest correlation with the difference of SOS, and increasing soil temperature and soil moisture could both promote the earlier SOS. After model adjustment by adding environmental factors into the original OLR model, the R-squared value between ground SOS and 5 years interval estimated LSP SOS reached to 0.913, while R-squared value between ground SOS and 10 years interval estimated LSP SOS reached to 0.902 and R-squared value between ground SOS and all years interval estimated LSP SOS reached to 0.898.

CONCLUSIONS
(a) (b) Figure 6. Correlation of SOS among LSP SOS, ground SOS and estimated LSP SOS ((a) correlation between LSP SOS and ground SOS, (b) adjusted correlation of SOS between estimated LSP SPS and ground SOS) In this study, we compared the difference of SOS in different time intervals and explored the relationship between ground SOS and LSP SOS. First, we used a two-step model based on Landsat datasets over years to identify the LSP SOS based on 1000 sample sites in New York State. Then, in order to quantify the difference of SOS under different time intervals, we compared the inter-annual difference of LSP SOS in 5 years interval, 10 years interval and all years interval. Finally, based on the possible impact of soil temperature (ST) and soil moisture (SM) on SOS, we processed the AGDD of soil temperature and soil moisture from ERA5 reanalysis products and improved the estimation accuracy when taking these factors into OLR model.
Our results indicate that LSP SOS can serve as a good proxy for reflecting ground SOS, especially in 5 years interval. Many sample sites have the close mean value of SOS among different time intervals in the long term (1995-2016) within 9 days. Due to the small changes affected by ST and SM in the urban built-up area over years, LSP SOS in different time intervals has high consistency. Croplands can be greatly affected by the other external environmental condition and human reaction, so the LSP SOS among time intervals is different. Moreover, annual variability of SOS agrees well in different time intervals, especially in low latitude areas. The AGDD of soil temperature and soil moisture have a certain influence on the difference of SOS between LSP SOS and ground SOS.