Impact Assessments of Typhoon Lekima on Forest Damages in Subtropical China Using Machine Learning Methods and Landsat 8 OLI Imagery

Abstract

Wind damage is one of the major factors affecting forest ecosystem sustainability, especially in the coastal region. Typhoon Lekima is among the top five most devastating typhoons in China and caused economic losses totaling over USD 8 billion in Zhejiang Province alone during 9–12 August 2019. However, there still is no assessment of its impacts on forests. Here we detected forest damage and its spatial distribution caused by Typhoon Lekima by classifying Landsat 8 OLI images using the random forest (RF) machine learning algorithm and the univariate image differencing (UID) method on the Google Earth Engine (GEE) platform. The accuracy assessment indicated a high overall accuracy (>87%) and kappa coefficient (>0.75) for forest-damage detection, as evaluated against field-investigated plot data, with better performance using the RF method. The total affected forest area by Lekima was 4598.87 km², accounting for 8.44% of the total forest area in Zhejiang Province. The light-, moderate- and severe-damage forest areas were 2106.29 km², 2024.26 km² and 469.76 km², respectively. Considering the damage severity, the net forest canopy loss fraction was 2.57%. The affected forest area and damage severity exhibited large spatial variations, which were affected by elevation, slope, precipitation and forest type. Our study indicated a larger uncertainty for affected forest area and a smaller uncertainty for the proportion of damage severity, based on multiple assessment approaches. This is among the first studies on forest damage due to typhoons at a regional scale in China, and the methods can be extended to examine the impacts of other super-strong typhoons on forests. Our study results on damage severity, spatial distribution and controlling factors could help local governments, the forest sector and forest landowners make decision on tree-planting planning and sustainable management after typhoon strikes and could also raise public and governmental awareness of typhoons’ damage on China’s inland forests.

1. Introduction

Typhoons, also called hurricanes, are a rotating, organized system of clouds and thunderstorms that originates from tropical or subtropical waters. Typhoons/hurricanes are usually accompanied by torrential rain, severe inundation and extremely strong wind velocity, causing damages to both human and natural systems [1]. Typhoons/hurricanes have caused about 2000 people’s death and USD 7 billion in losses worldwide every year. The recorded strongest hurricanes, the Great Miami and Galveston hurricanes and Hurricane Katrina, caused economic losses of USD 208.3, 130.2 and 118.8 billion, respectively, in the continental United States [2]. China is one of the countries most frequently and severely affected by typhoons. About seven typhoons pass over China’s land area every year, and the total accounts for about 30% of all the tropical cyclones in the western Pacific Ocean [3]. From 2000 to the present, over 100 strong typhoons (about five per year) have made landfall in China. Typhoons mainly affect the southeastern coastal areas in China, such as Zhejiang, Fujian, Guangdong, Hainan, and Guangxi provinces. There have been plenty of reports on addressing the effects of typhoons/hurricanes on forest damage; however, few studies have been conducted in China. One reason is due to too few field-observation data for typhoon damage monitoring. Another reason, possibly the most important one, is that few scientists in China are concerned about damage to forests. Governmental or public reports on economic losses due to typhoon damage generally account for losses in human death, living places, crop production and facilities but either no estimates or only a “rough” estimate of forest damage.

Typhoons/hurricanes can significantly affect the sustainable development of terrestrial ecosystems, especially forests, due to damage from strong wind, heavy rainfall and land erosion [4,5,6]. The estimation of large-scale forest damage is difficult because most typhoon/hurricane events only cause fractional loss of forest canopy [7]. Many approaches, including field investigation and satellite observations, have been developed to estimate forest damage caused by typhoons/hurricanes. Compared with field investigation, satellite-based estimations are more widely applied at a large scale. Various spatial-resolution satellites and onboard sensors, such as Landsat, ASTER, Sentinel, AVHRR and MODIS, are often chosen (e.g., [7,8,9,10,11,12]). In addition to the satellite platforms, the choice of damage-detection algorithm is also important. Many disturbance-detection algorithms have been developed and applied, such as LandTrendr, univariate image differencing (UID), selective principal component analysis (PCA), change vector analysis (CVA), post-classification comparison (PCC), vegetation change tracker (VCT) and breaks for additive season and trend monitor (BFAST) [13,14,15,16]. Several studies (e.g., [12,17,18,19,20]) applied the UID method to assess forest damage due to Hurricanes Andrew, George and Yolanda. Wang and Xu [13] compared the accuracy of UID, PCA, CVA and PCC in detecting forest damages due to Hurricane Katrina and concluded that the performance of UID and CVA was similar and that the PCA method was not stable, while a PCC was the most accurate method. Generally, forest damage was estimated based on these algorithms by comparing the changes in various vegetation indices (VIs) after and before a typhoon event. The most frequently used VIs include the normalized difference vegetation index (NDVI), EVI (enhanced vegetation index), soil-adjusted vegetation index (SAVI), normalized difference water/moisture index (NDWI/NDMI), as well as the tasseled cap greenness (TCG), brightness (TCB) and wetness (TCW) [9,10,20,21,22]. There have been many different arguments about which indices were the most robust. For example, Wang et al. [23] assessed forest damage after Hurricane Katrina using the NDVI, EVI and normalized difference infrared index (NDII) derived from the 1 km MODIS images and forest inventory data and indicated that NDII had higher disturbance-detection accuracy than either the NDVI or EVI because the NDII is more sensitive to canopy water content. In contrast, based on the BFAST method, Schultz et al. [24] compared the performance of various vegetation indices in detecting deforestation and concluded that wetness-related VIs (e.g., the NDWI and TCW) were spatially more accurate than greenness-related VIs (the NDVI and TCG).

In recent years, machine learning classification/regression models (e.g., random forest (RF), convolutional neural network (CNN) and artificial neural network (ANN)) have emerged and been widely applied in detecting disturbance events and their impacts [12,17,18]. These methods can combine the advantages of various VIs through integrating multiple VIs, training and validation datasets into different machine learning methods. Kathrin et al. [25] applied the RF method and further developed a two-step method to detect forest damages. Their study proved that the RF classifier could effectively detect forest-damage severity with robustness over 90%. Zhang et al. [10] assessed the performance of different machine learning methods, including support vector machines (SVM), RF, k-nearest neighbors (K-NN) and ANN, in detecting hurricane-caused forest damage and also indicated that RF has the highest robustness. In addition, Cohen et al. [14,26] proposed a new approach by integrating a two-period change-detection algorithm with machine learning methods in order to improve the detection accuracy of fractional forest disturbance. These studies suggested that the RF machine learning method could be adequately applied to estimating forest damage caused by strong winds. Unlike the conventional regression method, no data-distribution assumptions are needed for RF regression, and the RF regressor can detect interactions and higher-order relationships between predictors during the modeling process without the modeler specifying these terms [27].

Typhoon Lekima landed on the coastal region of Zhejiang Province, in subtropical China, during 9–12 August 2019. It was a super-strong typhoon with a maximum wind speed of 52 m/s when landing, which makes it a Category 16 typhoon (strongest) according to the Chinese GBT 19201-2006 standard or a Category 3 hurricane (devastating degree) according to the Saffir–Simpson Hurricane Scale. It is among the top five strongest typhoons in the history of China (http://www.cma.gov.cn/, accessed on 20 January 2021), but its economic loss is the highest recorded, reaching about USD 8 billion [28]. There have been some studies to explore the causes of high economic losses and the spatial patterns of rainfall and its impacts on soil erosion, flooding intensity and crop yield [28,29,30,31,32], but no studies have yet reported the losses in forest sectors because of the lack of accurate data on forest damage. In this study, based on plot-level training and validation data from field investigation and Landsat 8 OLI images, we applied the RF classification/regression method and a LandTrendr-based single vegetation index approach on the Google Earth Engine (GEE) platform to assess forest damages due to Typhoon Lekima in Zhejiang Province in subtropical China. Our objectives were to: (1) identify the affected forest area and its spatial distribution by Typhoon Lekima, (2) assess forest-damage severity and its spatial-distribution pattern, and (3) explore the impacts of multiple factors on forest-damage severity. This study is among the first estimations of forest damage due to typhoons at a regional scale in China. Our results will not only provide accurate forest-damage data for further assessing the impacts of typhoons on economic losses in forests but also help guide forest recovery management after typhoons.

2. Methods

2.1. Study Region

Zhejiang Province is located on the southeastern coast of China (118°01′–123°10′ E, 27°06′–31°11′ N) (Figure 1). Zhejiang Province includes 11 prefectural-level cities, including Huzhou, Jiaxing, Hangzhou, Shaoxing, Ningbo, Zhoushan, Jinhua, Quzhou, Taizhou, Lishui and Wenzhou City along the north to south gradient (Figure 2). It belongs to the subtropical monsoon climate zone with a mean annual precipitation of 980–2000 mm and a mean annual temperature of 15–18 °C. The terrain varies from mountains (average altitude of 800 m) in the southwest, to hills in the central region, to plains in the northeast [33]. The mountainous and hilly regions account for about 70.4% of the total land area, and most forests are located in these areas. The major forest types are subtropical evergreen needleleaf forest, evergreen broadleaf forest and mixed needleleaf and broadleaf forest [34]. In addition, it has 0.9 million ha bamboo area, accounting for 14.89% of total forest area and ranking the third-highest in bamboo area compared to other provinces in China. The current forest cover is 60.91% (6.05 million ha) in Zhejiang Province according to the ninth national forest resource inventory (NFRI) (http://forest.ckcest.cn/sd/si/zgslzy.html, accessed on 20 January 2021). The main soil types are red and yellow soils.

2.2. Typhoon Lekima

Typhoon events have occurred frequently and, in total, 43 typhoons have made landfall in Zhejiang Province during 1949–2018. Typhoons caused about USD 4.6 billion economic losses per year [35]. Typhoon Lekima was a super-strong typhoon that was ranked the third-strongest typhoon in Zhejiang Province since 1949 (http://www.cma.gov.cn/, accessed on 20 January 2021). Typhoon Lekima first made landfall at the coastal area of Wenling City, in southeastern Zhejiang Province, in the early morning of 10 August 2019 and moved out to neighboring Jiangsu Province at around 23:00 p.m. [30,31]. The central track of Lekima mainly passed across Taizhou, Shaoxing and Jiaxing cities in the east, reducing wind velocity from 40–50 m/s in Taizhou to 20–30 m/s in Jiaxing. It left China’s mainland from the coastal region of Shandong Province, lasting for about 44 h, which ranked as the sixth-longest retention time in China’s history since 1949. Although the main track crossed over eastern Zhejiang, strong wind outbreaks were reported throughout the entire province. About 7.57 million people were affected, and 45 people died in Zhejiang Province [30]. The mean rainfall during the landfall period (10–14 August) was 165 mm over the whole region, with higher rainfall, over 400 mm, in some areas along the eastern coast of Taizhou City and Ningbo City.

2.3. The Work Flows

To achieve its objectives, this study was conducted with the workflows in Figure 3. The main steps include: (1) data collection, including training and validation data, Landsat 8 OLI images and other assisting data; (2) Landsat images preprocessing; (3) masking out the nonforest area; (4) calculating VIs and other variables; (5) implementation of the UID method and RF model to detect forest damage; (6) validations for the derived forest damage; and (7) producing and analyzing the final products.

2.4. Data Sources

2.4.1. Landsat 8 OLI Images

The Landsat 8 OLI images were obtained from GEE platform (https://earthengine.google.com/, accessed on 20 January 2021). The Landsat 8 satellite has a 16 day repeating cycle, but, due to cloud cover and other noise during the rainy season in Zhejiang Province, it is difficult to collect a full coverage of clear images for the entire province right before and after Typhoon Lekima. Therefore, this study extended the collection time from 1 June to 8 August 2019 to represent the forest status before Lekima and from 12 August to the end of October 2019 to represent the forest status after Lekima. To smooth the noise caused by cloud, data transmission errors, ice/snow reflectance and others, we also collected the time series of Landsat 5, 7 and 8 images from 1986 to 2019 for the above two periods. The obtained Landsat images have already been atmospherically corrected using LEDAPS or LaSRC (surface reflectance code). Both algorithms used CFMASK as a built-in cloud-detection method and generated QA (quality assessment) bands that identified the pixels with adverse instrument, atmospheric or surficial conditions [36] We removed the pixels marked as cloud or other improper objects and filled them with clear pixels from other neighboring months to ensure that only cloud-free images were generated. Finally, annual mosaicked images were produced for the two periods (before and after Lekima) during 1986–2019.

2.4.2. Training and Validation Sample Data

To train the RF algorithms and validate the derived results, we collected 640 sample plots, in total, from sources including field-investigated plots (FIP) and Google-Earth- and Worldview-2-based visual interpretation plots (Figure 1). The field investigation was conducted from 7 October to 12 October 2019. The investigation started in the first landfall region (Wenling City) and went through the predesigned cities and plots. The high spatial-resolution (≤1 m) composite remote-sensing images from Baidu Map (https://map.baidu.com/, accessed on 20 January 2021) were used as base maps for visually identifying field trips and plots. The investigated plot area is 30 m × 30 m by matching the 30 m Landsat image grid cells. The recorded information included main forest types, the forest fraction in each plot, the fraction of forest damage and damage types (tree fall, branch fall or canopy fall). Obvious damages caused by other disturbances were excluded. In total, 200 plots affected by Typhoon Lekima were investigated (Figure 1). These plots were scattered in Taizhou City, Linhai City, Zhuji City and Ningbo City. Among them, 120 plots were randomly selected for RF model training, which left 80 plots to be used for validation.

To increase the representativeness, in addition to the above-mentioned samples, we also obtained sample plots based on available high-resolution Google Earth Pro and WorldView2 composite RGB images (data source: https://discover.digitalglobe.com/, accessed on 20 January 2021) using the visual interpretation (GEVI) method. Due to dense cloud cover in the study region, the WorldView2 image scenes after 12 August and before the end of 2019 were used to represent the post-Lekima forest status. The same methods for obtaining training and validation data have been applied in [37]. These GEVI sample plots were generally selected in obviously disturbed areas with affected areas larger than 20 × 20 m² that can be easily identified directly by eyes based on our long-term classification experience. In total, 440 sample plots were obtained through the GEVI method. Among them, 280 plots were randomly selected for model training, which left 160 plots to be used for validation. Including the FIP plots, the training plot numbers totaled 400, and the validation plots totaled 240.

2.4.3. Topographic and Climatic Data

Topographic factors such as elevation, slope and aspect were reported to significantly affect the forest damage caused by typhoons/hurricanes [1,6,23,38,39]. These data can be derived from the digital elevation model (DEM) data. The high spatial-resolution DEM data (30 m), originally sourced from the NASA/USGS Shuttle Radar Topography Mission (SRTM), were collected from the GEE platform using the ee.Terrain.products() function. The slope and aspect were calculated based on the DEM data.

In addition, many previous studies have also indicated that climatic factors, such as wind velocity and rainfall, affect the impact of typhoons/hurricanes on forest damages [1,6]. Therefore, we used wind velocity and rainfall as input variables in the RF model to monitor forest damage. Both datasets were collected from the CLDASV2.0 climatic products in China’s Meteorological Administration (CMA) website for the period 9–12 August 2019 (http://data.cma.cn/, accessed on 20 January 2021). These climate data have a spatial resolution of 0.0625° at an hourly time step. The data were further downscaled to a 30 m spatial resolution using the inverse distance weighted (IDW) interpolation method. The spatial distributions of elevation and total rainfall were shown in Figure 4.

2.5. Mask of Forest Areas

To eliminate the noise caused by nonforest areas, we first mask out the forest areas before the classification of Lekima damage on forests. The RF model with a classifier mode (ee.classifier.smileRandomForest() function in GEE) was applied for forest classification. Our collected field-investigation plots (200 plots in Figure 1) were selected and used as training (140 plots) and validation (60 plots) for forest-area classification. The NDWI, NDBI (normalized difference built-up index) and NDVI indices (Table 1) were used as input variables in the RF model. The classified forest images were validated, and the overall classification accuracy was 0.95, indicating a strong robustness for screening out forests. We assumed no land use/cover change before and after Typhoon Lekima (from 1 June to the end of October 2019).

2.6. The Calculations of VIs and Other Input Variables

In previous studies, VIs, including NDVI, EVI, NDII, CVI, SAVI and NDMI, were widely used to represent forest damage in response to typhoon/hurricanes ([6,9,10,21,22,23,37,40]). Therefore, we also used these VIs as inputs to the RF algorithm. The calculation methods for these VIs are shown in Table 1.

In addition to VIs, image-texture information can be used as input variables in the RF model to assist decision-making. Seven texture variables were used, including CONTRAST, IDM (inverse difference moment), CORR (correlation), DISS (dissimilarity), VAR (variance), ASM (angular second moment) and ENT (entropy). CONTRAST measures the local contrast of an image; IDM measures the homogeneity; ENT measures the randomness of a gray-level distribution; VAR measures how spread out the distribution of gray levels is; CORR measures the correlation between pairs of pixels; ASM measures the number of repeated pairs. These texture variables were widely applied in land-use/cover-change and disturbance-detection studies. We applied the gray-level co-occurrence matrix (GLCM) method to calculate these texture variables on the GEE platform.

2.7. Mapping Forest Damage Using Univariate Image Differencing of VIs

As indicated from many previous studies, changes in the NDVI, EVI, NDII, CVI and SAVI can be used to represent forest-damage rates due to typhoons/hurricanes. Therefore, we applied the univariate image differencing (UID) method to calculate the differences of these five VIs before and after Typhoon Lekima.

Before the calculation of difference in VI (ΔVI, represents the change rate), we applied the LandTrendr despiking algorithm to smooth the noise caused by unreasonable fluctuations in reflectance and then fit the VI trajectory over time. The LandTrendr trajectory-fitting and segmentation methods have been widely applied for disturbance detections [14,26,41,42]. Although this algorithm was designed for detecting long-term large-scale forest loss/gains, it has superior capability to smooth noise and fit VI trajectories. The despiking and trajectory-fitting algorithms can smooth noise over trends but faithfully capture abrupt change. The trajectories of VIs were fitted after a standard F-test using the p-value as the indicator of the goodness of fit. The available time-series of VIs from 1986 to 2019 for the period 1 June to 8 August and from 12 August to the end of October were used as inputs in the LandTrendr algorithm on the GEE platform. After smoothing, the ΔVI was calculated as:

(1)$\Delta \mathrm{VI}=({\mathrm{VI}}_{\mathrm{after}}-{\mathrm{VI}}_{\mathrm{before}})/{\mathrm{VI}}_{\mathrm{before}}$

where VI_after denotes the VI obtained from the mosaic images before Lekima (represented by the images from 1 June to 8 August 2019) and after Lekima (represented by the images from 12 August to the end of October 2019).

The change in the VI (ΔVI) represented the fractional loss (0–100%) of forest canopy, which was further grouped into four damage-severity categories, including no (<5% loss), light (5–20%), moderate (20–50%) and severe (>50%) damage. We speculated, with high probability, that severity between 0–5% was caused by system and image errors; therefore, the category with <5% loss was considered to be the forest area that was not affected by Typhoon Lekima.

2.8. Mapping Affected and Damaged Forest Area Using RF Classifier/Regressor Algorithms

The RF model with classifier and regressor modes was used for forest-damage detections. The collected 400 training plots (Figure 1) were all selected for RF model training and regression model cross-validation. We wrote a script and implemented the RF algorithm (numberOfTrees = 1000) on the GEE platform (ee.classifier.smileRandomForest()), and other parameters were set to the default. We set RF into regressor mode with the code “setOutputMode(‘regression’)”. Our purpose was to more accurately find out which forest pixels experienced damage and what the loss fractions of canopy (damage severity) were. The input variables, including image texture, topographic factors (elevation, slope and aspect), reflectance changes of individual bands (from ΔBand 1 to ΔBand 7) and changes in VIs (ΔVI) were used as RF model inputs. The same smoothing methods for VIs were applied to calculating the reflectance difference of bands before and after Lekima. During the implementation of RF, the z-score was calculated and the importance of each factor was computed and ranked (Table 2). The NDVI, B4, B2, SAVI, EVI, B1, B3 and NDMI had the highest z-score and thus had more obvious impacts on the detection of forest damages.

2.9. Accuracy Assessment

The affected forest area and damage-severity classifications, based on the RF and UID methods, were validated against validation sample data. The randomly selected validation plots in Figure 1 were used for validation. Accuracy assessment was conducted using an error matrix, which is a table that shows the performance of a method wherein each row of the matrix represents the predicted category, while the columns represent the actual class (or vice versa) [37]. Based on the error matrix, we can calculate different assessment indicators, including overall accuracy, user accuracy, producer accuracy and the kappa coefficient. Higher overall accuracy reflects stronger robustness of the results. The kappa coefficient is a measure of how the classification results compare to values assigned by chance. If the kappa coefficient equals 0, there is no agreement between the classified image and the reference image. If the kappa coefficient equals 1, then the classified image and the ground truth data are totally identical [43].

3. Results

3.1. Accurancy Assessment

The 240 field-investigated plot and GEVI plot data for validation listed in Figure 1 were used to evaluate the performance of classification results for the affected area based on the RF method and the single VI-based UID method. The consumer, producer and overall accuracy, as well as the kappa coefficient, were used as evaluation criteria (Table 3). The results indicated that all methods showed a high overall accuracy (>87%) and kappa coefficient (>0.75). The RF method had the highest overall accuracy (0.93) and kappa coefficient (0.86), indicating a stronger robustness. At a spatial scale, we also compared the detected and observed affected area and damage severity (Figure 5). The observed affected forest plots were mostly located within the detected damage area using the RF method. The light, moderate and severe damage severities reflected by the smoothed NDVI time-series data matched well with the detected damage categories and magnitudes using the RF method.

To compare the uncertainties from different methods, we further applied the ordinary least squares method to evaluate the relationship between the observed and estimated damage severities from the RFR and UID methods using a single VI and calculated the root mean square error (RMSE) and regression coefficients. The evaluation results are shown in Figure 5. The regression coefficients (R²) for RFR, EVI, NDII, NDVI and SAVI were 0.82, 0.68, 0.63, 0.62, 0.56 and 0.29, respectively (Figure 6). This suggested that RFR could accurately estimate the Typhoon-Lekima-caused forest-damage severity.

From the above analyses, we finally selected the RF method to estimate forest-damage area and severity. The UID method using a single VI was used to compare to the RF method and provide uncertainty ranges.

3.2. Spatial and Regional Patterns of Affected Forest Area

According to the RF algorithm, Typhoon Lekima has, in total, affected 4598.87 km² of forest, accounting for 8.44% of the total forest area in Zhejiang Province (Figure 7). Although the wind velocity declined along the typhoon track from the coastal to the inland area, the affected forest area did not show an increasing trend with increasing wind velocity. Instead, the affected forest areas were mostly located at the fringes of forest landscapes, where there is less forest cover and the area is more easily affected by strong wind and rainfall. In addition, less forest area was affected by Lekima in the western and southwestern regions because Lekima mainly passed through eastern Zhejiang Province. The affected proportions (area) for evergreen needleleaf forest, bamboo, mixed and evergreen broadleaf forests were 11.74% (1140.15 km²), 6.59% (1006.30 km²), 7.40% (1575.18 km²) and 3.84% (189.00 km²). As shown in Figure 1, most of the pure evergreen needleleaf forests are located at the edge of the forest landscape, which causes a relatively higher affected area. In addition, bamboo canopies are generally more susceptible to wind damage due to their taller and thinner stems; however, most bamboo forests are located in the lower-elevation and easy-access parts of Zhejiang Province and the top crown of most managed bamboo forests are trimmed to increase more product-allocation to bamboo shoots and to prevent damages from snowstorms and strong winds. Therefore, the affected proportion of bamboo forests was relatively lower than that of evergreen needleleaf forests.

At a prefecture-level city scale, we found that the largest affected forest areas were located in Hangzhou City (1076.18 km²) and Lishui City (610.81 km²), while the least were located in Jiaxing City and Zhoushan City, due to those cities having the lowest forest area (Table 4). Zhoushan City had the highest damaged forest fractions (18.74%), while Lishui City had the lowest affected forest fraction (4.97%). Lishui is in the west of Zhejiang, far from the typhoon track, and experienced the lowest rainfall and wind velocity during the typhoon period (Figure 2 and Figure 4). Zhoushan is closer to the typhoon track and the landing region, and most importantly, the forest distribution pattern there was more scattered and thus, had a larger portion of fringes. In addition, this city received the heaviest rainfall and higher wind velocity. The first landing coastal region, Taizhou City, had less damaged forest area, which may be because few forest areas are distributed along the coastline, and forests there are more adapted to the frequent typhoon outbreaks.

3.3. Spatial and Regional Patterns of Forest Damage Severity

Based on the RF regression model, the damaged forest area and its proportions with four categories of severity, including no (<5% loss of canopy), light (5–20%), moderate (20–50%) and severe (>50%) damage, was estimated. By accounting for damage severity, the net-loss fraction of forest canopy was 2.57% caused by Lekima in Zhejiang Province (Table 4; Figure 8). The area of light, moderate and severe damage severity accounted for 45.79%, 44.00% and 10.21% of the affected forest area. The greater damage severity generally occurred at the edge of forest landscapes, since these areas are more susceptible to wind, intense rainfall or floods. Forest-damage severity did not show an obvious tendency along the wind-velocity gradient from the south to the north of Zhejiang Province.

At a prefecture-level city scale, Taizhou (17.60%), Quzhou (16.54%) and Wenzhou (13.15%) had higher fractions of severe damage rates, while Huzhou, Hangzhou and Shaoxing had the lowest fractions (about 7.00%; Table 4). Quzhou and Taizhou also had relatively higher fractions of moderate damage severity, suggesting forests in these two cities were more susceptible to Lekima’s damages. Taizhou was the first landfall region of Typhoon Lekima. Although forests in this coastal city are generally adapted to frequent outbreaks of typhoons, the stronger wind velocity and more intense rainfall caused by Lekima still resulted in more damage in this region.

3.4. The Impact Factors on Forest Damage Caused by Lekima

The relationships between major environmental factors and forest-damage severity were analyzed to discover the driving factors to the unique spatial and regional patterns in forest damage caused by Lekima. The effects of slope on forest-damage severity were more complex but generally showed a higher forest damage severity when slopes were greater than 20° and lower than 60°, while the lowest damage severity was found in areas with slopes between 10–20° and >60° (Figure 9). Forest-damage severity showed an increasing trend with increasing elevation when elevation was lower than 1400 m; this trend then declined when elevation was greater than 1400 m. Forest-damage severity was not significantly affected by rainfall when rainfall was lower than 250 mm but greatly increased when rainfall was higher than 250 mm. This indicated that severely strong rainfall could significantly increase forest-damage probability. Evergreen needleleaf forests and mixed forests had greater damage severity compared to evergreen broadleaf forests and bamboo forests. As mentioned before, evergreen needleleaf forests were mostly distributed at the edge of the forest landscape, making them more easily damaged by strong wind and devastating flooding.

3.5. The Uncertainty of Affected Forest Area and Damage Severity

Many individual VIs have been applied to estimate forest damage due to hurricanes/typhoons and have proven their robustness; therefore, we also ran the UID approach using the smoothed ΔNDVI, ΔEVI, ΔCVI, ΔNDII and ΔSAVI. The spatial distribution patterns for affected forest area and damage severity were similar among these methods (Figure 10). However, significant differences existed for the estimated affected forest area and severity (Figure 11). The estimated affected forest area ranged between 4598.87 to 7726.78 km², with the lowest from the RF method and the highest from the ΔNDII method. All of the UID methods showed significantly larger affected forest area (ranging from 6128.94 km² to 7726.78 km²) than that estimated by the RF classifier/regressor method, which may because the RF method considered the overall changes of all these ΔVIs and other factors. For damage severity, the proportions for light, moderate and severe damage were 45.79–53.92%, 35.66–44.00% and 7.06–10.43%, respectively, with smaller uncertainty ranges compared with the uncertainty ranges of the affected area.

4. Discussion

4.1. Impact Factors on Typhoon-Caused Forest Damages

Many factors were found to affect hurricane/typhoon damage severity [1,6,38,39]. Wind velocity is one of the major factors affecting wind damage [6,44]. Martin and Ogden [44] even reviewed the impact of wind velocity on wind damage on a national scale. According to the Saffir–Simpson Hurricane Scale and maximum wind velocity (52 m/s), Typhoon Lekima was a Category 3 typhoon (H3; devastating degree). Our results indicated the net-loss fraction of canopy was 2.57%, which was significantly lower than those for Hurricanes Rita (H3; 12.37% loss), Katrina (H3; 12.38% loss) and Ivan (H3; 11.97% loss) occurred in the southern United States [6]. The similar wind velocity but different damage severity implied that, except for wind velocity, other factors could significantly affect forest-damage severity. This explained why our detected forest-damage severity did not show a declining tendency with wind velocity along the hurricane track from the south (52 m/s) to the north (23 m/s) (Figure 2). Indeed, our study found that topographic factors (slope and elevation) could significantly affect severity through secondary modifications to wind velocity and direction. We found that slopes ranging between 30° and 60° had the highest damage severity, and damage severity increased with elevation when elevation was lower than 1400 m. Most of the land area (74%) in Zhejiang Province is located in hilly or mountainous terrain. About 45% of the land area has a slope greater than 25°, and most forests are located in these areas. Many previous studies (e.g., [6,39,44,45]) have proven the effects of topography. Precipitation was a unique factor influencing damage severity during the Lekima landfall period in our study region. Heavy rain was accompanied with strong wind, with regional-average total rainfall of 165 mm and the highest rainfall over 400 mm during 9–12 August 2019. Many trees were actually uprooted by flood-caused soil erosion, especially in steep mountainous areas [28,29,30]. This further suggested that soil and tree conditions could greatly affect damage severity [39,44]. Our study further found that evergreen needleleaf forests had the highest affected area and damage severity. A reason is that most needleleaf forests are located at the edge of forest landscapes, where they are more susceptible to wind strikes. Another reason is that the stems and crowns of pure needleleaf forest are more vulnerable to wind damage, which was indicated in [6,39,44].

4.2. Uncertainties, Implications, and Outlooks

Although Landsat has proven successful in monitoring hurricane impacts, the acquisition of cloud-free Landsat data immediately after Typhoon Lekima (9–12 August 2019) is a challenge in subtropical China due to the rainy season. Our study collected all images within two to three months before (June to August) and after (August to October) Lekima to cover the entire study region. These long time-spans could either introduce noise from other disturbance events or the damaged forests could have recovered during this period. Therefore, our study may slightly over- or underestimate Lekima’s impacts on forest damage. In addition, it is challenging to determine what fraction of trees is dead in satellite images at 30 m resolution [7]; therefore, our study used the canopy-loss rates to represent damage severity, and our results for damage severity may not be equivalent to forest-mortality rates.

Besides image quality, the change-detection method may cause additional uncertainty. Many algorithms (e.g., BFAST, LandTrendr, machine/deep learning and VCT) and vegetation indices (e.g., NDVI, EVI, NDMI and TCW) have been applied and proven robust in detecting the impact of typhoons/hurricanes. We chose the RF algorithms, which used machine-learning methods to generate relationships between forest damage and to input variables trained with the field observational data. This method has previously proven more effective than some other approaches [10,25]. To look at the uncertainty ranges from different methods and VIs, we also applied the UID method, which is based on individual vegetation indices to fit VIs’ trajectories to determine the percentage of canopy changes caused by disturbance. Based on the accuracy assessments, our study also indicated that the RF method was more accurate than the single VI-based UID method to detect typhoon impacts as evaluated against the field-investigated data. The comparisons also suggested that RF methods estimated a significantly lower affected forest area compared with the UID method, but the forest fractions of different severity categories were not significantly different. Cohen et al. [14] proposed a combination with LandTrendr and RF methods to detect forest loss/gains after disturbance. They argued that this combined method is more accurate in detecting forest gain/loss due to large-scale disturbance events. Based on our previous study [46], we also found that this combination method can be effectively applied to detect the impacts of disturbance. In the future, we will apply this combined method to estimate forest damage in response to other devastating typhoon events in China.

Ground-truthing sample plot data is key to training the RF model to detect the effects of typhoons on forest damage and thus, reduce uncertainty. However, field observational data are often lacking in China. High spatial (generally <1 m) and temporal resolution images are now available through some online cloud systems, such as the Google Earth Engine platform, digitalglobe, USGS data warehouse and other commercial cloud platforms. Therefore, although there were few available field observations for forest damage after typhoons, we can still apply these high-resolution satellite images to visually interpret the affected forest area and damage severity and to provide training for the RF detection algorithms and validation for the detected results at a larger scale. However, visual interpretation generally uses personal expertise knowledge or changes in reflectance values and, thus, may bring in some uncertainties in the training and validation data.

Although there are some uncertainties, our study is among the first attempts to estimate the impacts of typhoons on forest damage based on satellite data at a regional scale in China. Our estimation can help assess forest ecosystem function and structure changes in response to typhoons and further help forest managers and local governments evaluate economic losses in the forest sector, as well as guide sustainable management practices for forest recovery after typhoon strikes. Our results indicated that evergreen needleleaf forests are more vulnerable to wind damage; wind damage severity increased with elevation when elevation < 1400 m, and damage severity significantly increased when rainfall was >250 mm during typhoon outbreaks. The pure needleleaf forests in the study region were mostly planted forests rather than natural forests. Therefore, we suggest avoid planting pure needleleaf forests in places with high elevation and slope. Forest managers and the forestry sector should pay more attention to forest areas with a high probability of soil erosion when typhoon events are accompanied by heavy rainfall within a short period. Recovery management after typhoons should give priority to these areas. Based on our estimated mean damage severity (2.57%) and NFRI-reported total forest volume (2.81 × 10⁸ m³ during 2014–2018), we can roughly estimate a loss of forest volume of 7.23 × 10⁶ m³ caused by Typhoon Lekima. Certainly, an accurate estimation of biomass and economic loss should consider more variables.

The coastal area in China has suffered from frequent and serious impacts from typhoons [35]. Although the vegetation in the affected area has adapted to typhoons’ impacts, super-strong typhoons, such as Lekima, Rammasun and Saomai, can result in severe damage to forest ecosystem function, such as biomass, carbon storage, soil erosion and productivity. Our study indicated that the RF machine learning method accompanied by field- and visually-interpreted sample plots could adequately detect forest-damage area and severity. In the future, we will extend this method of estimating the impacts of previous and coming super-strong typhoon events on forest structure and function. In addition, we will also estimate the economic losses and carbon-emission rates due to forest damage, which will fill in the data gaps for economic and ecological loss in the forest sector caused by typhoons, and thus, increase public and governmental concerns regarding typhoon-caused forest damage. Compared with other countries worldwide, China’s government and researchers pay less attention to wind damage on forest ecosystems.

5. Conclusions

This study provided a comprehensive assessment of the forest damage caused by Typhoon Lekima in Zhejiang Province in the subtropical China, based on the random forest classifier/regressor and UID algorithms using the Landsat 8 OLI images. Our results revealed that this typhoon disturbed 8.44% of the forests and resulted in 2.57% forest canopy loss. Most of the affected forest areas experienced light and moderate damage severity. The complex spatial variation patterns in affected forest area and damage severity generally resulted from the impacts of topographic, climatic and biotic factors. As reflected by the accuracy assessments, the random forest approach was proven superior to the UID method in wind-damage detection. Although some uncertainties existed, our study can provide a high spatial-resolution dataset for further assessing forest biomass and productivity loss and resultant economic loss; in addition, our study can also inform local governments and the forest sector, as well as forest landowners about damage severity, its spatial distribution patterns and major controlling factors, and thus, help rationally manage forest ecosystems to reduce typhoons’ damage.

Author Contributions

Conceptualization, G.C.; methodology, X.Z. and G.C.; validation, X.Z., L.C., H.J., J.H. and X.L.; formal analysis, X.Z.; investigation, X.Z., J.H. and L.C.; data curation, X.Z., H.J., X.L., J.H. and L.C.; writing—original draft preparation, X.Z. and G.C.; writing—review and editing, G.C. and X.W.; supervision, G.C. and X.W.; project administration: G.C.; funding acquisition: G.C. and X.W. All authors have read and agreed to the published version of the manuscript.

Funding

This study was funded by the Natural Science Foundation of Zhejiang Province (Grant number LY20C030001), the Scientific Research Foundation of Zhejiang A&F University (Grant number 2034020080) and the 111 Project (Grant number D18008).

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The data that support the findings of this study are available from the corresponding author (G.C.) upon justifiable request.

Acknowledgments

Authors are grateful to the field-observation data from Zhenming Shen in the Bureau of Agriculture and Rural Affairs of Lin’An District, Hangzhou, Zhejiang Province.

Conflicts of Interest

The authors declare no conflict of interest.

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Figures and Tables

Figure 1. The study region, forest types and sample plots for Typhoon Lekima’s impacts in Zhejiang Province. Note: FIP, field-investigated plots (200 plots); GEVI, Google Earth and Worldview-2 based visual interpretation plot (440 plots); Red circles or triangles: plots for training, black circles or triangles: plots for validation. Forest type source: [34].

Figure 2. The administrative boundaries of the prefectural-level cities and the path and wind velocity of Typhoon Lekima during 9–12 August 2019 in Zhejiang Province.

Figure 3. The workflow and study logic.

Figure 4. The digital elevation model (DEM) (m; Left) and cumulative precipitation (mm; Right) during 9–12 August 2019.

Figure 5. The selected Landsat 8 OLI images before and after Typhoon Lekima (Left two panels), the detected forest damage severities (light, moderate, and severe categories) and area based on the random forest algorithm and the NDVI time series for the period 8 August to the end of October during 2000–2019. Note: the green dots are the ground-truth sample plots.

Figure 6. The comparisons between detected and observed forest-damage severity based on different change-detection methods. (a): UID method based on ΔEVI; (b): ΔNDII; (c): ΔNDVI; (d): ΔCVI; (e): ΔSAVI; (f): random forest method.

Figure 6. The comparisons between detected and observed forest-damage severity based on different change-detection methods. (a): UID method based on ΔEVI; (b): ΔNDII; (c): ΔNDVI; (d): ΔCVI; (e): ΔSAVI; (f): random forest method.

Figure 7. The spatial distribution of forest areas affected by Typhoon Lekima in Zhejiang Province.

Figure 8. The spatial distribution of forest-damage rates (%) caused by Typhoon Lekima.

Figure 9. Impacts of elevation (m), slope (°), forest type and precipitation (mm) on damage severity for the affected forest areas by Typhoon Lekima. Note: forest type 1 is evergreen broadleaf forest; forest type 2 is bamboo forest; forest type 3 is evergreen needleleaf forest; forest type 4 is mixed forest.

Figure 9. Impacts of elevation (m), slope (°), forest type and precipitation (mm) on damage severity for the affected forest areas by Typhoon Lekima. Note: forest type 1 is evergreen broadleaf forest; forest type 2 is bamboo forest; forest type 3 is evergreen needleleaf forest; forest type 4 is mixed forest.

Figure 10. The detected affected forest area and damage severity using UID methods. (a) NDVI, (b) EVI; (c) NDII, (d) SAVI.

Figure 10. The detected affected forest area and damage severity using UID methods. (a) NDVI, (b) EVI; (c) NDII, (d) SAVI.

Figure 11. The comparisons of the detected affected forest area with different damage severity using the random forest regressor (RFR) method and the UID method based on single vegetation indices including NDVI, EVI, CVI, NDII, SAVI and CVI.

Table 1

The selected vegetation indices and their calculation methods.

Vegetation Indices	Formulas
NDVI	$\mathrm{NDVI}=\left(\mathrm{NIR}-\mathrm{Red}\right)/\left(\mathrm{NIR}+\mathrm{Red}\right)$
EVI	$\mathrm{EVI}=2.5\times \left(\mathrm{NIR}-\mathrm{Red}\right)/\left(\mathrm{NIR}+6\times \mathrm{Red}-7.5\times \mathrm{Blue}+1\right)$
NDII	$\mathrm{NDII}=\left(\mathrm{NIR}-\mathrm{SWIR}\right)/\left(\mathrm{NIR}+\mathrm{SWIR}\right)$
NDMI	$\mathrm{NDMI}=\left(\mathrm{NIR}-\mathrm{MIR}\right)/\left(\mathrm{NIR}+\mathrm{MIR}\right)$
SAVI	$\mathrm{SAVI}=\left(\mathrm{NIR}-\mathrm{Red}\right)/\left(\mathrm{NIR}+\mathrm{Red}+\mathrm{L}\right)\times \left(1+\mathrm{L}\right)$
NDWI	$\mathrm{NDWI}=\left(\mathrm{Green}-\mathrm{NIR}\right)/\left(\mathrm{Green}+\mathrm{NIR}\right)$
NDBI	$\mathrm{NDBI}=\left(\mathrm{MIR}-\mathrm{NIR}\right)/\left(\mathrm{MIR}+\mathrm{NIR}\right)$
CVI	$\mathrm{CVI}=\mathrm{NIR}\times \left(\mathrm{RED}/{\mathrm{GREEN}}^2\right)$

Table 2

The importance (z-score) and rank of input factors in the random forest model.

Factors	Importance	Rank
ΔNDVI	11.55	1
ΔB4	11.47	2
ΔB2	9.50	3
ΔSAVI	9.12	4
ΔEVI	8.76	5
ΔB1	7.98	6
ΔB3	6.96	7
ΔNDMI	6.44	8
ΔB7	6.29	9
ΔNDII	5.65	10
ΔCVI	4.87	11
Elevation	4.76	12
Precipitation	3.91	13
ΔB5	3.62	14
ΔB6	3.35	15
contrast	1.13	16
idm	1.02	17
corr	0.89	18
diss	0.87	19
var	0.80	20
asm	0.76	21
ent	0.51	22

Table 3

The performance evaluations for the classified affected area by Typhoon Lekima based on different methods against the observed sample plot data.

Methods *	Consumer Accuracy	Producers Accuracy	Overall Accuracy	Kappa
RF	0.92	0.94	0.93	0.86
NDVI	0.90	0.91	0.89	0.79
EVI	0.90	0.91	0.90	0.79
NDII	0.88	0.90	0.88	0.75
CVI	0.88	0.89	0.87	0.75
SAVI	0.90	0.90	0.89	0.79

* RF: random forest classifier/regressor method; NDVI, EVI, NDII, CVI and SAVI are selected single vegetation indices in the UID method.

Table 4

The forest areas and proportions of different damage severities affected by Typhoon Lekima in different prefecture-level cities of Zhejiang Province.

City	Forest Area (km2)	Affected Area (km2)	Fraction (%)	Light (%)	Moderate (%)	Severe (%)
Wenzhou	5647.32	370.43	6.56	45.38	41.47	13.15
Hangzhou	11,030.66	1076.18	9.76	48.46	44.52	7.03
Zhoushan	270.01	50.62	18.75	55.80	32.09	12.11
Huzhou	1891.40	245.30	12.97	63.43	29.54	7.03
Jiaxing	15.56	2.14	13.75	52.29	38.36	9.35
Jinhua	5802.89	578.78	9.97	39.39	51.22	9.39
Lishui	12,295.74	610.81	4.97	56.43	35.65	7.92
Ningbo	3458.70	464.73	13.44	57.45	34.59	7.96
Quzhou	5430.95	366.48	6.75	34.95	48.51	16.54
Shaoxing	4205.59	495.69	11.79	48.87	44.06	7.07
Taizhou	4464.09	337.71	7.57	32.98	49.42	17.60
Entire Zhejiang	54,512.91	4598.87	8.44	45.79	44.00	10.21