Processing

Two vegetation indices (VIs), NDVI (Normalized Difference vegetation index, Rouse et al., 1974) and NDII (Normalized Difference Infrared Index, Hardisky et al., 1983) were used for the overall and disturbance analysis using QGIS 3.4.2. NDVI was

18

chosen because it is the most widely used VI for various studies around the globe (Kerr

& Ostrovsky, 2003; Pettorelli et al., 2005; Huete, 2012) and has been shown to have a close relationship with leaf area index (LAI, Carlson & Ripley, 1997), and therefore forest productivity (Clark et al., 2008). The NDII was chosen because SWIR-NIR-based indices such as NDII can track defoliation better than does NDVI (Townsend et al., 2012). The topographic position index (TPI), slope steepness and aspect were calculated with the terrain function from R’s ‘raster’ package (Hijmans, 2019). Table 2.2 summarizes the meaning and calculation of topographical and vegetation indices.

Changes of VIs following each disturbance event, ΔVI, were calculated for each VI and disturbance events as: elevation model. R: red band; NIR: near infrared; and SWIR1: short-wave infrared.

Index, variable Meaning Calculation Method

Topographical

Elevation Altitude above sea level -

Slope aspect Surface orientation (e.g., facing North) Orientation of each cell

Slope steepness Angle to the horizontal (°)

Comparison of cell value with surrounding cells

TPI

Topographic Position Index. Landforms, positive values associated with ridges, negative values with

valleys¹ quantifies vegetation by measuring the difference between near-infrared (which vegetation strongly reflects) and red light (which vegetation absorbs), with high values indicating dense vegetation cover.

− +

NDII³

Normalized Difference Infrared Index. It is sensible to vegetation water content and vegetation changes associated with cyclone disturbance⁴ or defoliation⁵.

− 1

+ 1

1 Weiss (2001); ² Rouse et al. (1974); ³ Hardisky et al. (1983); ⁴ Wang et al. (2010); ⁵ Townsend et al. (2012).

19 2.2.5. Analysis of FFDP-FEF Representation

The analyses of the general FFDP-FEF representation (one Landsat image) and the representation for each typhoon disturbance (two images per event, before and after) were carried out similarly. First, mean values for VIs, ΔVIs (for typhoon disturbances only), and the three topographical variables of FEF and FFDP were statistically compared using the bootstrapped difference (via 5000 iterations, n = 289) as:

! "# $% _& − ! "# $% _' . (2.2) There was no significant difference between the FEF and the FFDP if the 95%

confidence interval of their mean difference (CI) includes 0. Slope aspects of FEF and FFDP were compared with a χ² test. Spearman’s ρ was used to examine correlations between VIs, and between VIs and topography of the FFDP and FEF.

For each bootstrapped comparison of means, we calculated the coefficient of variation (CV). The post- and pre-disturbance CVs were compared between the FEF and FFDP as:

() _*

) ₊) _& − () _*

) ₊) _' (2.3)

where ) * is pre-typhoon CV and ) + is post-typhoon CV. A bootstrapped comparison of means was used again to assess if typhoon-induced changes of heterogeneity of vegetation indices were different between the FFDP and the FEF.

Two thresholds were used to evaluate if a section of the FEF and FFDP was disturbed by each of the five typhoons:

ΔVI < 0 (2.4)

ΔVI < ΔVI _& − 0.5 × 2 _& (2.5)

Although a cell with ΔVI < 0 could be considered as a disturbed cell, a lower threshold (i.e., more negative ΔVI) helps minimize mis-identification (e.g., due to differences in image quality). Thus, the zone model threshold, Equation 2.5, was also used to define disturbed cells (Jin et al., 2013). Spatial variation in typhoon frequency, defined as the proportion of pixels experiencing different frequencies of typhoon disturbance (0 to 5), in the FEF and FFDP were compared with a χ² test.

20

The Euclidean distance measure (ED) has been used as a metric to estimate representativeness of eddy flux towers (Hargrove et al., 2003; He et al., 2015) and sampling networks (Hoffman et al., 2013) at the continental scale, and observation plots at the landscape scale (Langford et al., 2016). We used ED to pinpoint well or under-represented areas of the FEF by the FFDP. The VIs (NDVI, NDII) and topographical variables (elevation, slope, TPI) were normalized (i.e., scaled to 0–1)—hence leading to a distance metric akin to the Mahalanobis distance—and then used to compute multi-dimensional distances between all FEF cells and FFDP cells. Several EDs were calculated for each reserve cell because there was more than one cell within the FFDP.

Thus, the minimal values among these distances—the minimal ED (minED)—was kept for each of the reserve cells on the basis that one cell of the FFDP is more representative of a given FEF cell than the others (Figure 2.2). Two types of minED maps were produced: VI-based ED (a combination of NDVI and NDII) and topography-based ED (elevation, slope steepness, TPI). Correlations between topography-based and VI-based minEDs were conducted using Spearman’s ρ.

21

Figure 2.2. Overview of the steps for minimal Euclidean distances (minED) calculation based on vegetation indices (VIs) and topographical variables. This is a modification of regular Euclidean distance (ED) computation in which there is more than one reference cell (the plot contains x cells), where only the minimal value is kept. Steps leading to minED maps are described here: After grouping several data layers (step 1), Euclidean distances were computed between each plot cell and each reserve cell based on variable values standardized on a 0–1 scale (step 2). Then, step 3 involved keeping only the minimal Euclidean distance value computed for each reserve cell (among the n values, where n = the number of cells within the plot), on the basis that there was one plot cell that offered the best representation of the considered reserve cell. Finally, step 4 led to the construction of georeferenced rasters with minED values computed either with vegetation indices or topographical variables. MinEDs show how well the reserve is represented based on the inserted variables.

22 2.2.6. Testing Alternative Forest Plot Designs

Because vegetation cover varies with elevation and the FFDP spans a limited subset of the FEF elevational range, four alternative plot designs were tested to evaluate if increasing plot elevational range and splitting the large plot to several smaller plots, without changing the total plot area, would improve the plot’s landscape representativeness. The four plot designs we tested are: one large-square plot (one 500 m × 500 m, ≥30% of FEF elevational range), four square plots (four 250 m × 250 m,

≥50% of FEF elevational range), two rectangular plots (two 250 m × 500 m transects,

≥50% of FEF elevational range), and four rectangular plots (four 125 m × 500 m, ≥50%

of FEF elevational range). The first strategy is like the FFDP, except that it covers a broader elevational range of the FEF, whereas the other three strategies use subplots that are dispersed across the landscape and cover even greater elevational ranges of the FEF. Ten plot replicates per plot design were created by randomly generating plot locations and orientations within the FEF. Plot locations were selected based on the pre-defined elevation range threshold and by limiting the number of cells with no data to <10%. Because pre-typhoon and post-typhoon images for each typhoon event were less than 30 days apart from the typhoon event (except for Dujuan), it was difficult to get disturbance-related images that lack cloud cover as it rains more than 220 days annually at the FEF (Chang et al., 2017). After removing clouds, maintaining an elevation range >50% was not achievable. Thus, in the analysis of disturbance representativeness, the elevation range threshold was set to ≥30% of the FEF range;

however, the analysis for Typhoon Aere and for the one-square strategy were excluded because clouds limited the image area available for the analyses. Ten to twelve replicates were used for each strategy. The same spatial areas were used to study the four typhoons even though they had varying cloud cover. This permitted us to compare typhoons effects in the exact same locations, and thus to maintain integrity in comparing alternative plot designs by keeping the other variables (such as slope steepness or other unmeasured factors) fixed across different typhoons. Some of the

23

replicates of the same strategy overlap partially in areas, however subplots of each replicate do not overlap (Figure C1).

The analysis of the landscape level representativeness of the four strategies was done similar to the FFDP-FEF comparisons. Elevation was excluded from the minED calculation because this variable was used to create alternative plots. Finally, VI-based minED values for each strategy were compared through a pairwise Wilcoxon test, with Bonferroni adjustment, for each typhoon and the overall FEF. They were then ranked on a 1 to 4 scale for the overall analysis, or 1 to 3 for the disturbance analysis based on the mean minimal Euclidean distance calculated across all replicates, where the strategy with the lowest mean gets first rank.

2.3. Results

2.3.1. Overall Plot Differences between FEF and FFDP

The difference between the FEF and FFDP was significant for both VIs as well as slope and elevation (Figure 2.3). Vegetation cover was greater in the FFDP, with both mean NDII and NDVI being significantly higher for pixels with the FFDP boundary (NDII = 0.25 ± 0.01, NDVI = 0.80 ± 0.009) than the greater FEF (NDII = 0.21 ± 0.03 and NDVI = 0.78 ± 0.033). Mean elevation and slope steepness were significantly higher in the reserve, whereas there was no difference for TPI between the FEF and FFDP. The proportion of the cells with different slope aspects were also significantly different between the FFDP and the FEF (χ² = 41.86, df = 7, p < 0.001, Figure B2), with a greater proportion of NW-N slope aspects and a lower proportion of SW-S aspect in the FFDP than in the FEF. The range of values represented by the FFDP in relation to the FEF was greater than 66% for the TPI, 58% for slope steepness, but only 10% for elevation, 9% for NDVI, and 23% for NDII. The interquartile range (IQR i.e., Q1–Q3) of the FFDP represented only 6.3% of the NDVI, 0% of the NDII, and 9.5% of the elevation, but 33.7% (slope steepness) and 81.1% (TPI) of the FEF’s respective IQRs.

24

Figure 2.3. Boxplots of vegetation indices and topographical variables in the Fushan Forest Dynamics Plot (FFDP) and the Fushan Experimental Forest (FEF) based on Landsat-8 and JAXA digital elevation model (30 m spatial resolution). Significant differences between FFDP and FEF are shown with an asterisk (*), based on bootstrapped comparison on means.

In both the FFDP and FEF, the NDII and NDVI were moderately correlated (ρ > 0.6, p < 0.001, Table B1). However, the topographical variables were only weakly correlated (ρ < 0.2), except for TPI and elevation in the FFDP (ρ = 0.5, p < 0.001). There were no significant correlations between VIs and topographical variables except for elevation, which was negatively correlated with NDVI (ρFEF = −0.56, ρFFDP = −0.34, p < 0.001) and NDII (ρ^FEF = −0.32, ρ^FFDP = −0.21, p < 0.001) in both the FFDP and FEF.

The minimal Euclidean distance varied across the reserve (Figure 2.4). There was a significant, but very weak relationship between VIs-based minED and topographical variable-based minED (ρ = 0.05, p < 0.001).

25

Figure 2.4. Minimal Euclidean distance between the Fushan Forest Dynamics Plot (FFDP) and the Fushan Experimental Forest (FEF) calculated from (a) VIs (left) and topographical variables (right); (b) and the disturbance effect associated with Typhoon Nari (ΔVI-based). White surface within the FEF was either obscured by cloud cover or is non-forested area. The FFDP boundary is shown via the red square.

2.3.2. Disturbances Effect and Representativeness

The mean ΔVIs were mostly negative for the five disturbances, indicating that all typhoons led to vegetation loss (Table B2). For Soudelor, Herb, and Nari, there was less vegetation loss in the FFDP than in the reserve (Figure 2.5, Table B2). There was, however, no significant difference for ΔNDVI with Dujuan and ΔNDII with Nari between the FFDP and the FEF, but there was a greater vegetation loss within the FFPD than the FEF for Typhoon Aere (positive 95% CI, Table B2). Among the five typhoon disturbances, ΔNDVI of the FFDP encompassed 12.5% (Soudelor) to 47.7% (Aere) of the ΔNDVI range of the FEF, and ΔNDII of the FFDP covers 17.6% (Soudelor) to 35.4%

(Dujuan) of the ΔNDII range of the FEF. However, looking at IQRs, the FFDP contained between 31.4% (Soudelor) and 75.9% (Dujuan) of the reserve’s IQR for ΔNDII, and between 34.2% (Aere) and 71.0% (Dujuan) for ΔNDVI.

26

Figure 2.5. ΔNDII and ΔNDVI (calculated as the difference between post- and pre-disturbance values) for Typhoon Nari in the Fushan Forest Dynamics Plot (FFDP) and the Fushan Experimental Forest (FEF). Negative values indicate loss of vegetation.

Significant differences between FFDP and FEF based on bootstrap comparison on means are shown with an asterisk (*).

ΔNDVI and ΔNDII were moderately correlated (0.4 < ρ < 0.7) for Typhoons Herb, Nari, and Soudelor for the FEF, were strongly correlated (0.7 < ρ < 0.9) for Typhoon Soudelor within the FFDP and for Typhoon Dujuan for the FEF and were very-strongly correlated (ρ > 0.9) for Typhoon Dujuan for both the FFDP and the FEF. However, for Typhoon Aere, the correlation between the two VIs was weak (ρ = 0.05 and ρ = 0.28 in the FFDP and the FEF, respectively). Topographical variables were weakly, if ever, correlated with ΔVIs with absolute Spearman’s ρ always below 0.2 (e.g., for elevation).

NDVI-based CV ratios were different in FFDP and FEF except for Typhoons Herb and Nari (Table B3). For both the FFDP and FEF, the ratio was below 1 (i.e., greater post- than before-typhoon CVs), indicating that typhoons caused an increase in the variation of vegetation cover. In addition, the CV ratios of the FEF were significantly closer to 1 than those of the FFDP, signifying that there was a greater difference between pre- and post-disturbance vegetation states within the FFDP than for the

27

greater FEF (Table B3). For NDII-based CVs, significant differences between CV ratios in the FEF and FFDP were seen for Typhoons Nari, Dujuan, and Soudelor. The magnitude of change in NDII-based CVs relative to zero (i.e., no change) varied in both FFDP and FEF, with either increasing, decreasing, or equal values depending on the typhoon considered (Table B3).

The χ² tests indicated different patterns in damage frequencies between FFDP and FEF based on ΔVI < 0 threshold (χ² = 22.48, df = 3, p < 0.001 for NDVI; χ² = 16.26, df = 4, p = 0.003 for NDII). For NDVI, the FFDP had more cells within the intermediate damage frequency classes (2, 3) than the FEF, whereas the proportion of higher frequency (4, 5) and no damage (0) were lower in the FFDP than in the FEF (Figure 2.6). The NDII also showed damage of higher intermediate frequency (3) within the FFDP than within the FEF (Table B4). Using the ΔVI < mean − 0.5 × SD threshold, the differences between the FFDP and the FEF remained (χ² = 49.00, df = 5, p < 0.001 for NDVI; χ² = 33.88, df = 5, p < 0.001 for NDII), even though less disturbance frequency was detected, and there was an increased proportion of undisturbed areas (class 0), especially in the FFDP.

28

Figure 2.6. (a) Spatial representations of disturbance occurrences for NDVI (top panels) and NDII (bottom panels) for 5 typhoons (listed in Table 2.1) based on ΔVI < 0 (left) and ΔVI < mean – 0.5 × SD (right) thresholds. Damage class ranks from 0, being undisturbed for all five typhoons, to 5, having been disturbed by all typhoons. The Fushan Forest Dynamics Plot (FFDP) is shown in red, while the Fushan Experimental Forest (FEF) boundary is in black. White areas were either obscured by clouds or are outside FEF-FFDP limits. (b) Proportions of each damage frequency class within the FEF and the FFDP for both vegetation indices and based on the two thresholds (ΔVI < 0, top; ΔVI < mean – 0.5 × SD, bottom). The percentage of cells for each damage class for both thresholds are shown in Table B4.

2.3.3. Alternative Plot Designs

The replicates of the same plot design had significantly different vegetation and topography from the FEF (Table C1). The reserve had greater VIs (denser vegetation) than half of the alternative plot designs (e.g., 22 out of 40 replicates based on NDVI, Table C1). For topographical variables, most of the alternative designs led to significant differences with the FEF for slope steepness (24 out of 40 replicates, Table C1). However, the aspects proportions of the alternative plots were significantly different between the FEF and every alternative plot design replicate (p < 0.01), except one replicate of the four-rectangular plot design (χ² test, χ² = 12.33, df = 7, p = 0.09).

Despite numerous differences between the FEF and the alternative plots for NDVI and NDII, one plot of the two-rectangular plot design strategy, one of the four-squares

29

strategy, one of the four-rectangular plot, and one of the one-square plot strategy showed no significant differences in NDVI and NDII from those of the FEF (Table C1), despite having much greater elevation ranges.

The topography-based minED of alternative plot designs was significantly greater than that of the FFDP for all but two replicates (rectangles design, plot 8 and four-squares design, plot 3). Figure 2.7 shows the minED values calculated from the FFDP and the four-rectangles design replicates (ranked first in plot designs overall comparison, Table 2.3), whereas other plot designs are shown in Figure C2. However, VI-based minEDs showed the opposite trend, as most alternative plots led to improved representativeness of NDVI and NDII (i.e., for 30 out of 40 replicates the VI-based minED was lower than that of the FFDP, Figure 2.7 and Figure C2).

Figure 2.7. Minimal Euclidean distances (minED) measured between the Fushan Experimental Forest (FEF) and the Fushan Forest Dynamics Plot (FFDP) for the four-rectangular plot design strategy (all designs shown in Figure C2). Two types of minED were calculated: based on topographical variables (left panel), and based on vegetation indices (VI, right panel). Significant differences between the minED calculated for the FFDP and the alternative plots shown with an asterisk (*, see Table C1).

30

For Typhoon Nari, the comparison between alternative plot designs and the greater FEF showed that most of the alternative plot designs were less exposed to disturbances (i.e., ΔVI345678359:6 > ΔVI_;<;, Table C2) than is the FEF. Higher damages within the alternative plot designs were found for 1 of 11 replicates for the four-rectangular plot design, 4 out of 12 for the four-squares plot strategy, and 5 out of 10 replicates for the two-rectangular plot strategy. However, in the four-rectangular plot design, three replicates were not significantly different from the reserve for ΔNDII (Table C2).

The FFDP had a lower mean minED than most of the alternative plot design strategies for Typhoons Dujuan, Herb, and Soudelor (Table C3). However, it was not the case for Typhoon Nari, wherein several alternative plot designs had smaller minEDs than the FFDP. See Figure 2.8 for Typhoons Herb and Nari for the four-rectangular design (other typhoons and plot designs shown in Figure C3).

Figure 2.8. Minimal Euclidean distances (minED) for changes of vegetation indices (ΔVIs) associated with Typhoons Herb and Nari between the Fushan Forest Dynamics Plot (FFDP) and the four-rectangular plot design strategy within the Fushan Experimental Forest (FEF). Data for other plot designs and typhoons are shown in Figure C3. Significant differences between the FFDP and the alternative plots are shown with an asterisk (*) based on 95% CIs (see Table C3).

31

The mean minED of plot designs strategies were significantly different in each typhoon event as well as for the overall vegetation cover (Table C4). The two-rectangular plot design had the best ranking (i.e., smallest mean minED, meaning best FEF representation) among the three alternative plot design strategies across the four typhoon events, but it did not consistently have the smallest minED for all of the four events (Table 2.3). For the overall vegetation cover, the one-square plot design replicates led to the highest mean minEDs (i.e., least representativeness) among the four alternative plot design strategies (Table 2.3).

Table 2.3. Ranks of mean minimal Euclidean distances (ED) of four alternative plot design strategies based on vegetation indices for representativeness of difference effects (i.e., ΔVIs) and of the overall vegetation cover for the Fushan Experimental Forest. Significant differences between means of different plot design strategies were detected in all cases with a Wilcoxon test (p < 0.05, with Bonferroni adjustment), p values shown in Table C4.

Strategies Typhoon Disturbances

Overall Dujuan Herb Nari Soudelor Total

two rectangles 3 1 1 2 7 2

four squares 2 2 3 1 8 3

four rectangles 1 3 2 3 9 1

one large square - - - - - 4

2.4. Discussion

2.4.1. Overall Representativeness

The greater vegetation cover as indicated by both VIs, the overall gentler slopes, and lower altitude of the FFDP in relation to the FEF, and the significant differences in slope aspect between the FFDP and the FEF (Figure 2.3) all suggest that the plot does not adequately represent the greater landscape, both in terms of vegetation and topographic heterogeneity. Given the much smaller area of the FFDP (25 ha) compared to the FEF (1097 ha) and the rough mountain topography of the reserve, it is not surprising that the FFDP cannot adequately represent the topography of FEF.

However, although the FFDP covers only 10% of the elevation range of the FEF, it covers more than 50% of the landform variation (TPI) and slope steepness of the FEF.

32

The high coverage of the landform variability and slope steepness, however, does not lead to a high degree representativeness in vegetation as the range of VIs of the FFDP covers only 9% of the NDVI and 23% of the NDII ranges of the FEF (Figure 2.3). The lack of significant correlation between the topographical variables and vegetation

The high coverage of the landform variability and slope steepness, however, does not lead to a high degree representativeness in vegetation as the range of VIs of the FFDP covers only 9% of the NDVI and 23% of the NDII ranges of the FEF (Figure 2.3). The lack of significant correlation between the topographical variables and vegetation