A Long-Term Vegetation Recovery Estimation for Mt . JouJou Using Multi-Date SPOT 1 , 2 , and 4 Images

Vegetation recovery monitoring is critical for assessing denudation areas after landslides have occurred. A long-term and broad area investigation using remote sensing techniques is an efficient and cost-effective approach incorporating the consideration of radiometric correction and seasonality variations across multi-date satellite images. This paper investigates long-term vegetation recovery using 14 SPOT satellite images spanning from 1999 to 2011 over the landslide area of Mt. Jou-Jou in central Taiwan, which was caused by the Chi-Chi earthquake in 1999. The vegetation status was evaluated by the Normalized Difference Vegetation Index (NDVI) with radiometric correction between multi-date images based on pseudoinvariant features, and subsequently a vegetation recovery rate (VRR) model was empirically established after seasonality adjustment was performed on the multi-date NDVI images. An increasing tendency of the vegetation recovery in the landslide area of Mt. Jou-Jou appeared based on the NDVI value rising to 0.367 in March 2011 from −0.044 right after the catastrophic earthquake. The vegetation recovery rate with seasonality adjustment approached 81.5% for the total area and 81.3% for the landslide area through 12 years succession. The seasonality adjustment also enhanced the VRR model with a determination coefficient that increased from 0.883 to 0.916 for the landslide area and from 0.584 to 0.915 for the total area, highlighting the necessity of seasonality adjustment in multi-date vegetation observations using satellite images. Furthermore, the association between precipitation and NDVI was discussed, and the inverse relationship with the reoccurrence of high-intensity short-duration rainfall and yearly heavy rainfall was observed, in agreement with the on-site investigation.


Introduction
Located in a sub-tropical and seismic area, Taiwan often suffers from the impacts of earthquakes, typhoons, and torrential rains, which induce a high probability occurrence of landslides.On 21 September 1999, the Chi-Chi earthquake with ML 7.3, which is the most serious natural catastrophe in Taiwan of the past century, shocked central Taiwan and significantly changed the geographical features of the area.Earthquake-induced landslides represent one of the most hazardous impacts after severe seismic events [1][2][3][4][5][6][7][8][9].Thousands of landslide spots induced by the Chi-Chi earthquake produced a large amount of soil and stones that could turn into debris flow and threaten the residents and agricultural activities [10][11][12][13].The subsequent intense rainfall following the Chi-Chi earthquake exacerbated the situation by causing an unexpected increase in the collapse of the terrain slopes because of the extent of the bare land area in the landslide area.Mt.Jou-Jou was one of the most seriously damaged areas.This broad area of the upper hillsides was bare of vegetation, and dramatically increased the suspended solids in the air causing a serious air quality problem.Thus, Mt.Jou-Jou was established as a Nature Reserve Area by the Taiwan Forest Bureau in 2000, while several ecological engineering attempts were implemented to enhance the vegetation recovery and ecological restoration to stabilize the landslide area.
To evaluate the effect of these ecological restoration projects, a long-term ground observation program was initiated in Mt.Jou-Jou.The distribution of vegetation, soil composition, water level, and the condition of spillways and hillsides were observed in field investigations [14].However, on-site monitoring is labor-intensive and costly for long-term monitoring over such a broad mountainous region.To determine the efficacy of ecological restoration projects, an estimate of the vegetation recovery rate using remote sensing data is more efficient and cost-effective for broad area monitoring [15][16][17][18][19].
Remote sensing techniques have been recognized as a useful tool for the processing and display of spatial data and have been variously employed to provide an effective assessment of remedial measures and aided in decision analyses [20][21][22].Remote sensing data from the Advanced Very High Resolution Radiometer (AVHRR) sensor on the National Oceanic and Atmospheric Administration's (NOAA's) operational series of meteorological satellites were previously used to monitor vegetation dynamics for Africa over a 19-month period [23].As the AVHRR sensor continues operating, Harris et al. (2014) [24] used 25 years of GIMMS (Global Inventory Modelling and Mapping Studies) AVHRR-NDVI data to investigate vegetation cover persistence across southern Africa.Digital remotely sensed images acquired from a wide series of sensors have proven to be effective sources of information to explore and monitor the earth's surface over the last two decades [25][26][27][28][29].The combination of the multi-spectral optical SPOT-4 and the airborne CARABAS-II very-high frequency (VHF) synthetic aperture radar (SAR) has also proved useful for standwise stem volume estimations in forestry applications [30].Lin et al. (2002) [31] used remote sensing and GIS (geographic information system) techniques to assess landslides and debris flows.An extensive bibliography of works is available on the use of differential SAR interferometry (DInSAR) and multi-interferograms SAR interferometry (A-DInSAR) for landslide monitoring and integrating with in situ monitoring instrumentation [32][33][34].A series of Landsat satellites have provided numerous opportunities in earth surface observations [35].For instance, Landsat-5 thematic mapper (TM) imagery was applied to study post-fire vegetation recovery in North-East Victoria [36].Landsat-7 enhanced thematic mapper plus (ETM+) satellite imagery was used to identify spectrally anomalous locations with a 98% overall accuracy for rangeland management [37].A combination of the Landsat 8 satellite and SAR imageries change detection method has been proposed for landslide mapping in vegetated areas [38].Remotely sensed data extracted using the fuzzy c-mean algorithm and GIS, coupled with sediment yield estimation and vegetated buffer strip models, were used to evaluate the efficiency of erosion control for land use changes from cultivated to reforested land [39].One other significant advantage of remote sensing and GIS techniques over the traditional numerical methods is their capability of visually displaying the results of spatial assessments and with relative ease.However, for a long-term remote sensing observation, image sources from different satellite sensors may cause a systematic error on multi-date images.A series of SPOT satellite sensors could be a possible solution to provide a more consistent imaging system for a study period of more than a decade.
Based on a couple of previous reports using satellite images for observation purposes, the vegetation of Mt.Jou-Jou was claimed to recover its vegetation by 47% in 2000 [40], 59% in 2001 [17], 66% in 2002 [41], and 86-89% in 2005 [19,42].However, the on-site vegetation recovery of Mt.Jou-Jou was observed to be less than the projected vegetation recovery rate (VRR) based on previous studies due to several possible reasons.First of all, a very limited number of satellite images used for evaluating vegetation status may statistically mislead a regression model with deficient representativeness.Secondly, the projected VRR of previous studies established based on short-term observations could be insufficient for a long-term trajectory prediction.Moreover, natural seasonal effects affecting vegetation variation should be taken into account when evaluating the vegetation recovering status.Therefore, in order to precisely assess the vegetation restoration condition of Mt.Jou-Jou after the Chi-Chi earthquake, this study applied the Normalized Difference Vegetation Index (NDVI) derived from 14 SPOT satellite images with radiometric correction based on pseudoinvariant features.Additionally, the vegetation status was evaluated for the landslide area and the total area separately across multi-date images spanning from 1999 to 2011.Moreover, a seasonality adjustment was used to accommodate seasonal variations among multi-date images to fairly represent a long-term vegetation recovery.Furthermore, since high-intensity short-duration rainfall events have been considered as the main triggering factors for shallow slope failures [14,43,44], the vegetation recovery rate affected by the corresponding precipitation was analyzed to quantitatively reveal the impact of high-intensity short-duration rainfall on the vegetation variation in Mt.Jou-Jou.

Study Site
The pronunciation of "Jou" is the same as "9" in Mandarin, and 99 means numerous in Chinese; the special landscape of numerous peaks on Mt.Jou-Jou is shown in Figure 1.Mt. Jou-Jou, with a watershed area of 4396 ha, is located on the north shore of the Wu Stream in Nan-Tou County (Figure 1).The gravel layer with high hydraulic conductivity within the area of Mt.Jou-Jou is 1000 m thick.The slope ranges from 50 • to 85 • and over 50% of the slope area lies between 60~70 • [45].The main hillslopes of Mt.Jou-Jou are characterized as raised-straight and double-raised shapes.The geology of Mt.Jou-Jou is dominated by Pliocene and Pleistocene rocks.Mt.Jou-Jou is mostly comprised of Toukoshan formations consisting of three layers: an upper layer of bulk gravel rocks, a middle layer of sand, clay, and, gravel interbeds, and a lower layer of sand, shale, and thin gravel interbeds [14].The surface of Mt.Jou-Jou is rugged, with many deep gullies due to the consolidated gravel layer under dry conditions.The surface geologic formation erodes into the deeper layers during the rainy season.Nearby Mt.Jou-Jou are the Shuang-dong and Chelongpu faults, the latter of which triggered the Chi-Chi earthquake.Mt.Jou-Jou has a mean annual rainfall of about 1800 mm with a non-uniform temporal and spatial distribution.After the Chi-Chi earthquake, a five-year on-site monitoring program was conducted and revealed that 80.2% of the annual rainfall on Mt.Jou-Jou occurs during the rainy season (April-September), compared to 19.8% during the dry season (October-March) [14].Shallow failures frequently occurred on Mt.Jou-Jou in response to high-intensity and short-duration rainfalls, as well as earthquake events.
Remote Sens. 2017, 9, 893 3 of 17 (NDVI) derived from 14 SPOT satellite images with radiometric correction based on pseudoinvariant features.Additionally, the vegetation status was evaluated for the landslide area and the total area separately across multi-date images spanning from 1999 to 2011.Moreover, a seasonality adjustment was used to accommodate seasonal variations among multi-date images to fairly represent a longterm vegetation recovery.Furthermore, since high-intensity short-duration rainfall events have been considered as the main triggering factors for shallow slope failures [14,43,44], the vegetation recovery rate affected by the corresponding precipitation was analyzed to quantitatively reveal the impact of high-intensity short-duration rainfall on the vegetation variation in Mt.Jou-Jou.

Study Site
The pronunciation of "Jou" is the same as "9" in Mandarin, and 99 means numerous in Chinese; the special landscape of numerous peaks on Mt.Jou-Jou is shown in Figure 1.Mt. Jou-Jou, with a watershed area of 4396 ha, is located on the north shore of the Wu Stream in Nan-Tou County (Figure 1).The gravel layer with high hydraulic conductivity within the area of Mt.Jou-Jou is 1000 m thick.The slope ranges from 50° to 85° and over 50% of the slope area lies between 60~70° [45].The main hillslopes of Mt.Jou-Jou are characterized as raised-straight and double-raised shapes.The geology of Mt.Jou-Jou is dominated by Pliocene and Pleistocene rocks.Mt.Jou-Jou is mostly comprised of Toukoshan formations consisting of three layers: an upper layer of bulk gravel rocks, a middle layer of sand, clay, and, gravel interbeds, and a lower layer of sand, shale, and thin gravel interbeds [14].The surface of Mt.Jou-Jou is rugged, with many deep gullies due to the consolidated gravel layer under dry conditions.The surface geologic formation erodes into the deeper layers during the rainy season.Nearby Mt.Jou-Jou are the Shuang-dong and Chelongpu faults, the latter of which triggered the Chi-Chi earthquake.Mt.Jou-Jou has a mean annual rainfall of about 1800 mm with a non-uniform temporal and spatial distribution.After the Chi-Chi earthquake, a five-year on-site monitoring program was conducted and revealed that 80.2% of the annual rainfall on Mt.Jou-Jou occurs during the rainy season (April-September), compared to 19.8% during the dry season (October-March) [14].Shallow failures frequently occurred on Mt.Jou-Jou in response to high-intensity and short-duration rainfalls, as well as earthquake events.In response, the Council of Agriculture in Taiwan designated an 1198-ha area of Mt.Jou-Jou as a nature reservation for a vegetative succession study based on the Taiwan Cultural Heritage Protection Law.The vegetation is comprised of woody plants, such as Trema orientalis, Ficus erecta Thunb, Koelreuteria henryi Dumm, Cyclobalanopsis galuca, Arundo formosana, Pinus taiwanensis Hay, Pinus morrisonicola Hay, and Arundo formosana (the dominant species on Mt.Jou-Jou) [16].

SPOT Satellite Images
Fourteen SPOT 1, 2, and 4 satellite images with a 20-m spatial resolution were collected from 1999 to 2011 (Table 1).The image taken on 1 April 1999 covers the pre-earthquake area of Mt.Jou-Jou and represents the original pre-earthquake condition.The image taken on 27 September 1999 was recorded immediately after the Chi-Chi earthquake and represents the environment at its most damaged condition.The other 13 sequential images demonstrate the subsequent recovery in the vegetation condition after the earthquake until 2011, when the missions of SPOT 1, 2, and 4 satellites terminated in 2002, 2009, and 2013, respectively.Multi-temporal satellite images are often used in long-term programs of the earth's surface monitoring and change detection under the assumption of similar imaging conditions.In addition to geometric corrections, a radiometric correction has a critical influence on the analysis of multi-date satellite imagery.Due to the differences in the imaging environments for the different acquisition dates, radiometric correction is needed to eliminate the effects of these variations.Both deterministic and empirical methods can be used for radiometric correction; however, various imaging conditions must be considered with deterministic radiometric correction, such as the radiometric differences induced by the incidence angle of the sun; solar elevation; atmospheric absorption and scattering; variations in scan angle; and system noise, in addition to the reflectance of the objects.The empirical method is therefore more popularly used, and was adopted in this study.Pseudoinvariant Features (PIFs), which are possible "dark" and "bright" reference points with constant radiance, such as buildings, roads, and water, can be used as temporally-independent reference points on satellite images [21,[46][47][48].In this research, an aerial photograph with a better spatial resolution was used as a ground truth for reference to identify the locations of PIFs.The radiometrical correction was performed band by band using the 10 July 2003 image as the reference satellite image.A total of 93 pseudoinvariant features (bright and dark reference points) were selected on satellite images (Figure 2).A linear regression model can be established between the corrected and reference images to reduce the radiometric differences as: in which DN' and DN are the digital numbers on the corrected and reference images, respectively; m and k are constants for each band of all satellite images, and can be calculated using a least-squares difference of PIFs that appear both on the reference image and the corrected images.Through this equation, the multi-temporal images were then radiometrically corrected by examining the band-to-band scatter-grams in which the pixels display little variation between the imaging dates.
The radiometrically corrected false color images are shown in Figure 3.
Remote Sens. 2017, 9, 893 5 of 17 in which DN' and DN are the digital numbers on the corrected and reference images, respectively; m and k are constants for each band of all satellite images, and can be calculated using a least-squares difference of PIFs that appear both on the reference image and the corrected images.Through this equation, the multi-temporal images were then radiometrically corrected by examining the band-toband scatter-grams in which the pixels display little variation between the imaging dates.The radiometrically corrected false color images are shown in Figure 3.  in which DN' and DN are the digital numbers on the corrected and reference images, respectively; m and k are constants for each band of all satellite images, and can be calculated using a least-squares difference of PIFs that appear both on the reference image and the corrected images.Through this equation, the multi-temporal images were then radiometrically corrected by examining the band-toband scatter-grams in which the pixels display little variation between the imaging dates.The radiometrically corrected false color images are shown in Figure 3.

The Normalized Difference Vegetation Index (NDVI)
The extent and localization of the vegetation recovery were estimated by using the Normalized Difference Vegetation Index (NDVI).Based on the characteristic of chlorophyll absorption decreasing at red wavelengths and increasing in the infrared wavelengths, the spectral response of living vegetation can be measured by the band ratio of reflected red over infrared wavelengths [49].Because chlorophyll exhibits high reflectance in the near-infrared range and low reflectance in the red range, a high NDVI value represents an area covered with healthy vegetation.Recently, NDVI was adopted as an index to identify the inundated terrain and damaged areas where vegetation was seriously destroyed by the Sumatra tsunami [26].With the aid of probability theory, an alternative methodology for assessing VI (vegetation index) efficiency was also developed for vegetation variation [50].To quantify the growth of vegetation, NDVI was analyzed to give an overview of the vegetation before and after the earthquake [51].NDVI can be calculated as: in which NIR is the brightness of the near infrared band, and RED is the brightness of the red band.
Theoretically NDVI ranges between possible values of −1 and 1.

Landslide Area Extraction
In order to evaluate the vegetation status for the same denudation location, this study extracted the landslide area by referencing the maximum denudation area from a series of radiometrically corrected images corresponding with negative NDVI values and excluding areas with a slope of less than 10 degrees to avoid the misclassification of river banks/channels.In addition, total and nonlandslide areas were extracted from satellite images of Mt.Jou-Jou for comparison.The areas of the total, non-landslide, and landslide areas are 30.4km 2 , 24.5 km 2 , and 5.9 km 2 , respectively.

The Normalized Difference Vegetation Index (NDVI)
The extent and localization of the vegetation recovery were estimated by using the Normalized Difference Vegetation Index (NDVI).Based on the characteristic of chlorophyll absorption decreasing at red wavelengths and increasing in the infrared wavelengths, the spectral response of living vegetation can be measured by the band ratio of reflected red over infrared wavelengths [49].Because chlorophyll exhibits high reflectance in the near-infrared range and low reflectance in the red range, a high NDVI value represents an area covered with healthy vegetation.Recently, NDVI was adopted as an index to identify the inundated terrain and damaged areas where vegetation was seriously destroyed by the Sumatra tsunami [26].With the aid of probability theory, an alternative methodology for assessing VI (vegetation index) efficiency was also developed for vegetation variation [50].To quantify the growth of vegetation, NDVI was analyzed to give an overview of the vegetation before and after the earthquake [51].NDVI can be calculated as: in which NIR is the brightness of the near infrared band, and RED is the brightness of the red band.
Theoretically NDVI ranges between possible values of −1 and 1.

Landslide Area Extraction
In order to evaluate the vegetation status for the same denudation location, this study extracted the landslide area by referencing the maximum denudation area from a series of radiometrically corrected images corresponding with negative NDVI values and excluding areas with a slope of less than 10 degrees to avoid the misclassification of river banks/channels.In addition, total and non-landslide areas were extracted from satellite images of Mt.Jou-Jou for comparison.The areas of the total, non-landslide, and landslide areas are 30.4km 2 , 24.5 km 2 , and 5.9 km 2 , respectively.

Seasonality Adjustment
A simple seasonality adjustment was employed to accommodate seasonal variations among multi-date images.The present study assumed that the spatial distribution of natural seasonal influence is uniform across an image.The non-landslide area, which was only affected by seasonality and not the earthquake, was taken as a reference area.Specifically, the NDVI across all pixels of the non-landslide area was extracted from each image and averaged to provide a reference baseline.The differences between the reference baseline and NDVI values of the reference area for each image were calculated and utilized as offset values for seasonality adjustment.After adjusting the landslide and total areas NDVI using the offset values, the adjusted NDVI can be obtained and the modified VRR can be consequently estimated for each image.Thus, the modified VRR can be used to monitor the timing of vegetation recovery and provide insight into the post-earthquake vegetation response.

Estimation of Vegetation Recovery Rate (VRR)
NDVI variation can be used to express the level of vegetation recovery in damaged areas.Identified by the field investigation and aerial photographs, pixels with positive NDVI values were considered as areas covered by vegetation, and were then included in binary vegetation maps.Based on the difference of vegetated area between these thematic maps derived from the pre-earthquake image and the post-earthquake images, a vegetation recovery rate (VRR) can be calculated at any specific time point, as follows: in which N 0 , N 1 , and N 2 represent the number of pixels of the area covered by vegetation as derived from the satellite images acquired before the earthquake, right after the earthquake, and at a specific time point after the earthquake, respectively.A positive VRR means that the vegetation condition is recovering and reaches full recovery as the VRR approaches 100%.

Rainfall Data
Daily rainfall data were collected from Shuang-dong station, the nearest rainfall station to Mt. Jou-Jou, from 1999 to 2011.Based on the heavy rainfall definition established by the Taiwan Central Weather Bureau in 2004, daily rainfall exceeding 130 mm, 200 mm, and 350 mm can be considered as extremely heavy rain, torrential rain, and extremely torrential rain, respectively.The corresponding number of days with rainfall exceeding the abovementioned thresholds in Shuang-dong station were extracted and the 24-h rainfalls above 130 mm were aggregated on a yearly basis for the study period.

SPOT Satellite Images
To quantify the vegetation of Mt.Jou-Jou during the consecutive 13 years after the earthquake, the NDVI was calculated to approximate the vegetation on the 14 frames of SPOT images, because the spectral vegetation index is highly correlated with the green leaf biomass or projected green-leaf area.First of all, the radiometric correction was used to eliminate the imaging difference between multi-date satellite images.By taking the 10 July 2003 image with a solar elevation of 68.7 • and an azimuth of 93.3 • as the reference image due to its lowest shadow influence, a total of 93 PIFs of Mt.Jou-Jou were selected on all satellite images, as shown in Figure 2. Linear regression models were established between the reference image and other images for each band to reduce the radiometric differences.The adjusted multi-date satellite images are presented in chromatic consistency in Figure 3.The abundant satellite images provide a longer temporal coverage and more spatial information to better estimate the vegetation status of Mt.Jou-Jou.

NDVI Analysis
The NDVI was calculated on a pixel-by-pixel basis for all radiometrically corrected images (Figure 4).Locations excluding the river channel of the Wu stream with negative NDVI values were identified as landslide areas and are shown in brown in Figure 4. Approximately 5.9 km 2 earthquake-inducing destruction to the vegetation was identified, which was previously a full-grown forest on Mt.Jou-Jou.For comparative purposes, this 5.9 km 2 area was classified as the landslide area.The NDVI values of the total, non-landslide, and landslide areas were extracted and are shown in Figure 5, which reveals an overall increasing tendency in NDVI after the earthquake despite fluctuations due to associated typhoon events.Figure 6 shows that the NDVI value dropped dramatically immediately after the earthquake for the landslide area and fluctuated in the subsequent years.By comparing the thematic maps in Figure 4, one can easily identify the significant vegetation recovery in the landslide areas.One typical site of landslides was observed at close range in the years 2000, 2001, 2004, and 2008 to illustrate the vegetation recovery (Figure 7).
Table 1 illustrates the NDVI values for the total, non-landslide, and landslide areas of Mt.Jou-Jou.The pre-quake vegetation condition in terms of the NDVI value was excellent, with an average NDVI value of 0.523 for the whole landscape, whereas those spots where the landslide occurred after the quake had a slightly lower NDVI value (0.518).The average NDVI value of the landslide area sharply declined to −0.044 after the earthquake.Twelve years later, the NDVI value on 13 March 2011 had risen to 0.431, 0.448, and 0.367 for the total, non-landslide, and landslide areas, respectively.

VRR and Modified VRR
The VRR and the modified VRR, which were calculated using the adjusted NDVI, were estimated for the landslide and total areas in Mt.Jou-Jou (Table 2) and show a relatively quick growth in the early periods after the earthquake and then a gradual improvement in vegetation recovery.The modified VRR of Mt.Jou-Jou in 2011 approached 81.3% and 81.5% for the landslide and total areas, respectively.In general, the modified VRR in the landslide area is higher than that in the total area in the first two years due to artificial aerial afforestation right after the earthquake.
In order to evaluate the VRR of Mt.Jou-Jou, a trend analysis of VRR in terms of the time after the earthquake was performed for the total area (Figure 8a) and landslide area (Figure 8b).Both the VRR without seasonality adjustment and the modified VRR using adjusted NDVI were analyzed.The coefficient of determination, denoted as R 2 , and p-values were reported to provide measurements of goodness-of-fit and significance levels for VRR and modified VRR models, respectively.P-values from the VRR and modified VRR of the total and landslide area indicate that all recovery rates and intercept values are significant.Based on the reported R 2 values, the model for the landslide area performs better than that of the total area.For the total area (Figure 8a), the R 2 value for the modified VRR model is 0.915, which is a significant improvement compared to the R 2 value for the original VRR model (0.584).The modified VRR model of the landslide area has a higher R 2 value of 0.916, while the R 2 value for the VRR model of the landslide area is 0.883 (Figure 8b).
According to the formula of THE modified VRR model for the landslide area, the vegetation of the landslide area at Mt. Jou-Jou is estimated to require approximately 68 years to reach full recovery (Figure 8b).Additionally, for the total area, a relatively low VRR value of 0.031 was observed in 2 January 2001, whereas two relatively high VRR values, 0.715 and 0.772, were observed in 27 February 2003 and 10 July 2003, respectively.After the seasonality adjustment, the relatively low VRR observed in 2 January 2001 was adjusted upward to 0.393, whereas the two relatively high VRR values observed in 2003 were adjusted downward to 0.598 and 0.516.In terms of R 2 and p-values, the modified VRR model outperforms the original VRR model without seasonality adjustment for both areas, indicating that the seasonality adjustment provides a better estimation of the vegetation recovery status.

The Relationship between NDVI and Corresponding Precipitation
In addition to seasonal effects, high-intensity short-duration rainfall is another major factor associated with vegetation recovery and triggering shallow slope failures [14,43,44,52].According to the early-warning regulation in Taiwan, the Central Weather Bureau issues a heavy rainfall alert publically when rainfall exceeds 130 mm in order to minimize possible impacts, such as landslides and debris flows.To illustrate the relationship between high-intensity short-duration rainfall and the NDVI, the values for the number of days per year with 24-h rainfall exceeding 130 mm, 200 mm, and 350 mm were obtained for each year (Figure 9a

The Relationship between NDVI and Corresponding Precipitation
In addition to seasonal effects, high-intensity short-duration rainfall is another major factor associated with vegetation recovery and triggering shallow slope failures [14,43,44,52].According to the early-warning regulation in Taiwan, the Central Weather Bureau issues a heavy rainfall alert publically when rainfall exceeds 130 mm in order to minimize possible impacts, such as landslides and debris flows.To illustrate the relationship between high-intensity short-duration rainfall and the NDVI, the values for the number of days per year with 24-h rainfall exceeding 130 mm, 200 mm, and 350 mm were obtained for each year (Figure 9a).

The Relationship between NDVI and Corresponding Precipitation
In addition to seasonal effects, high-intensity short-duration rainfall is another major factor associated with vegetation recovery and triggering shallow slope failures [14,43,44,52].According to the early-warning regulation in Taiwan, the Central Weather Bureau issues a heavy rainfall alert publically when rainfall exceeds 130 mm in order to minimize possible impacts, such as landslides and debris flows.To illustrate the relationship between high-intensity short-duration rainfall and the NDVI, the values for the number of days per year with 24-h rainfall exceeding 130 mm, 200 mm, and 350 mm were obtained for each year (Figure 9a).Furthermore, the yearly differences in the NDVI, ∆NDVI, were calculated to investigate the relationship between the gain/loss of NDVI for the total and landslide areas and yearly heavy rainfalls, which comprise a total of 24-h rainfall exceeding 130 mm (Figure 9b).In Figure 9b Furthermore, the yearly differences in the NDVI, ΔNDVI, were calculated to investigate the relationship between the gain/loss of NDVI for the total and landslide areas and yearly heavy rainfalls, which comprise a total of 24-h rainfall exceeding 130 mm (Figure 9b).In Figure 9b

Conclusions
In this research, multi-date SPOT 1, 2, and 4 images have proved to be useful for comprehensively monitoring the spatial and temporal variation of the vegetation on Mt.Jou-Jou after the atmospheric correction.By using a series of 14 SPOT images of Mt.Jou-Jou spanning from 1999 to 2011, the NDVI was analyzed and the corresponding VRR model was empirically established.The VRR was analyzed for the consistent total and landslide areas across multi-date images.To eliminate the influence of seasonal variations, a seasonality adjustment was performed on the satellite images and the association between precipitation and NDVI was discussed.
This study supports the notion that radiometric correction is critical for the long-term detection of change using multi-date satellite images.Relative radiometric correction based on PIFs provides a simplified method to efficiently adjust the differences in the imaging environments for different acquisition dates to ensure a comparable intensity reflected by an object with a constant spectrum.Over 12 years of vegetation succession of Mt.Jou-Jou after the catastrophic earthquake, the NDVI value had risen from 0.278 to 0.431 and from −0.044 to 0.367 for the total and landslide areas, respectively, and the modified VRR approached 81.5% for the total area and 81.3% for the landslide area.The fact that the correlation coefficient of the modified VRR is higher than the original VRR highlights the necessity of seasonality adjustment in multi-date vegetation observations using satellite images.For the total area, the modified VRR regression model has an R 2 value of 0.915, with a significant improvement in the VRR model with an R 2 value of 0.584.For the landslide area, the seasonality adjustment enhanced the VRR with an R 2 value that increased from 0.883 to 0.916.Furthermore, the association between precipitation and NDVI was discussed, and the inverse relationship with the reoccurrence of high-intensity short-duration rainfall and the yearly heavy rainfalls was observed in agreement with the on-site investigation.Further work should be conducted to constantly monitor the vegetation recovery by utilizing more up-to-date satellite imagery and different series of sensors for establishing a reliable regression model of vegetation recovery for Mt.Jou-Jou by considering seasonal change and environmental impacts.

Figure 7 .
Figure 7. Close-range photographs showing the vegetation recovery of a typical landslide spot at Mt. Jou-Jou shot in (a) 2000; (b) 2001; (c) 2004; (d) 2008 and (e) corresponding landslide area NDVI value from satellite images.

Figure 7 .
Figure 7. Close-range photographs showing the vegetation recovery of a typical landslide spot at Mt. Jou-Jou shot in (a) 2000; (b) 2001; (c) 2004; (d) 2008 and (e) corresponding landslide area NDVI value from satellite images.

Figure 7 .
Figure 7. Close-range photographs showing the vegetation recovery of a typical landslide spot at Mt. Jou-Jou shot in (a) 2000; (b) 2001; (c) 2004; (d) 2008 and (e) corresponding landslide area NDVI value from satellite images.
The NDVI values in 2001 and 2003 were estimated by averaging two images in 2001 (2 January 2001 and 10 November 2001) and 2003 (27 February 2003 and 10 July 2003), respectively.In Figure 9a, no rainfall record exceeds 130 mm in 2002 and 2003, and a large increase in NDVI is observed.Likewise, in the years 2010 and 2011, the total area NDVIs increase as no heavy rainfall occurs.In 2004, as the number of days of heavy rainfall exceeding 130 mm and 200 mm increases to four and two, respectively, a large drop in the NDVI is observed.Three, two, and one days in 2008 had rainfall records above 130 mm, 200 mm, and 350 mm, respectively, and the corresponding total area NDVI decreased.

Figure 8 .
Figure 8. VRR and modified VRR for (a) total area and (b) landslide area.

Figure 8 .
Figure 8. VRR and modified VRR for (a) total area and (b) landslide area.
, during years when the total rainfall ≥130 mm increases, ∆NDVI generally decreases, and vice versa.Likewise, in the years 2002, 2003, 2010, and 2011, no rainfall record exceeds 130 mm and positive ∆NDVIs are observed.The largest drop in the NDVI is observed in 2004 with a correspondingly large negative ∆NDVI value, which was associated with the impacts of Typhoon Mindulle.According to the on-site observations by Chen and Wu (2006), Typhoon Mindulle brought almost 600 mm rainfall during 2 July 2004 to 5 July 2004 to Mt. Jou-Jou.A peak discharge of 36.9 m 3 /s with a rainfall intensity of 67.0 mm/h was observed in the major stream in Mt.Jou-Jou on 2 July 2004 [14].Such heavy rainfall caused serious landslides and degraded the vegetation in Mt.Jou-Jou.Based on the abovementioned observations, the relationship between rainfall and NDVI infers an inverse association.Additionally, based on the research results of Hsu et al. [53], 2002 to 2004 was identified as the most severe period of drought in the 21st century in Taiwan with yearly total rainfall and rainfall days far below normal.In the present study, no rainfall records exceeding 130 mm were observed during 2002 to 2003, which is in good agreement with Hsu et al. and Tsai and Yang[53,54].Furthermore, in Figure9b, the difference between the NDVI values of the total area and landslide area generally reduces with time, indicating that the vegetation condition of the landslide area is gradually recovering.Remote Sens. 2017, 9, 893 13 of 17 , during years when the total rainfall ≥130 mm increases, ΔNDVI generally decreases, and vice versa.Likewise, in the years 2002, 2003, 2010, and 2011, no rainfall record exceeds 130 mm and positive ΔNDVIs are observed.The largest drop in the NDVI is observed in 2004 with a correspondingly large negative ΔNDVI value, which was associated with the impacts of Typhoon Mindulle.According to the on-site observations by Chen and Wu (2006), Typhoon Mindulle brought almost 600 mm rainfall during 2 July 2004 to 5 July 2004 to Mt. Jou-Jou.A peak discharge of 36.9 m 3 /s with a rainfall intensity of 67.0 mm/h was observed in the major stream in Mt.Jou-Jou on 2 July 2004 [14].Such heavy rainfall caused serious landslides and degraded the vegetation in Mt.Jou-Jou.Based on the abovementioned observations, the relationship between rainfall and NDVI infers an inverse association.Additionally, based on the research results of Hsu et al. [53], 2002 to 2004 was identified as the most severe period of drought in the 21st century in Taiwan with yearly total rainfall and rainfall days far below normal.In the present study, no rainfall records exceeding 130 mm were observed during 2002 to 2003, which is in good agreement with Hsu et al. and Tsai and Yang[53,54].Furthermore, in Figure9b, the difference between the NDVI values of the total area and landslide area generally reduces with time, indicating that the vegetation condition of the landslide area is gradually recovering.

Figure 9 .
Figure 9. NDVI and ΔNDVI values varying with (a) number of days with precipitation over rainfall warning thresholds and (b) yearly total of 24-h rainfall exceeding 130 mm.

Figure 9 .
Figure 9. NDVI and ∆NDVI values varying with (a) number of days with precipitation over rainfall warning thresholds and (b) yearly total of 24-h rainfall exceeding 130 mm.

Table 1 .
NDVI values derived from SPOT satellite images.Image date is provided in yyyy/mm/dd format.
Note: NDVI values of pre-and post-Chi-Chi earthquake were marked in bold.

Table 2 .
VRR and modified VRR of the total and landslide areas.Image date is provided in yyyy/mm/dd format.Values mark in bold correspond to the discussion in Section 3.3., and one days in 2008 had rainfall records above 130 mm, 200 mm, and 350 mm, respectively, and the corresponding total area NDVI decreased. two

Table 2 .
VRR and modified VRR of the total and landslide areas.Image date is provided in yyyy/mm/dd format.
Note: Values mark in bold correspond to the discussion in Section 3.3