Analysis of Canopy Gaps of Coastal Broadleaf Forest Plantations in Northeast Taiwan Using UAV Lidar and the Weibull Distribution

Canopy gaps are pivotal for monitoring forest ecosystem dynamics. Conventional field methods are time-consuming and labor intensive, making them impractical for regional mapping and systematic monitoring. Gaps may be delineated using airborne lidar or aerial photographs acquired from a manned aircraft. However, high cost in data acquisition and low flexibility in flight logistics significantly reduce the accessibility of the approaches. To address these issues, this study utilized miniature light detection and ranging (lidar) onboard an unmanned aircraft vehicle (UAVlidar) to map forest canopy gaps of young and mature broadleaf forest plantations along the coast of northeastern Taiwan. This study also used UAV photographs (UAVphoto) for the same task for comparison purposes. The canopy height models were derived from UAVlidar and UAVphoto with the availability of a digital terrain model from UAVlidar. Canopy gap distributions of the forests were modeled with the power-law zeta and Weibull distributions. The performance of UAVlidar was found to be superior to UAVphoto in delineating the gap distribution through ground observation, mainly due to lidar’s ability to detect small canopy gaps. There were apparent differences of the power-law zeta distributions for the young and mature forest stands with the exponents λ of 1.36 (1.45) and 1.71 (1.61) for UAVlidar and UAVphoto, respectively, suggesting that larger canopy gaps were present within the younger stands. The canopy layer of mature forest stands was homogeneous, and the size distributions of both sensors and methods were insensitive to the spatial extent of the monitored area. Contrarily, the young forests were heterogeneous, but only UAVlidar with the Weibull distribution responded to the change of spatial extent. This study demonstrates that using the Weibull distribution to analyze canopy gap from high-spatial resolution UAVlidar may provide detailed information of regional forest canopy of coastal broadleaf forests.


Introduction
Canopy gaps (or "gaps", used interchangeably hereafter), caused by large branch losses, individual tree death or several treefalls due to natural processes or disturbance, play a pivotal role in regulating forest ecosystems [1][2][3][4][5][6][7]. Canopy gap dynamics are directly associated with regeneration and succession [8,9], structures [10,11], management [12,13] and disturbance [2,14]. The metrics commonly utilized to quantify canopy gaps are the gap depth and size. Gap depth is defined as openings in the forest canopy extending down to an average height ≤ 2 m aboveground [15]. High variation in canopy gap size may be attributed to biotic and abiotic factors such as the number of trees that have fallen, died or been removed. The range of gap sizes of hardwood and broadleaf stands typically vary Figure 1. (a) The study region, consisting of coastal young and mature broadleaf forest plantations (b) located in the subtropical region of northeast Taiwan (the star). The background image of (a) acquired in November 2015 was from Google Earth. (c) Monthly precipitation (gray-colored bars, the left y-axis) and air temperatures (the black line, the secondary y-axis) records of 1990-2021 for the study site during the observation period.

UAV Data Acquisition
UAVlidar and UAVphoto data were acquired on July 20, 2020 using a LiAir V scanning system (GreenValley International, Berkeley, CA, USA) and Sony α7 RII (Sony Group Corporation, Tokyo, Japan), respectively, carried by a Matrice M600 pro UAV (Da-Jiang Innovations, Shenzhen, China) at an altitude of 180 m a.g.l. The average point density was 11.2 points m −2 . The UAV lidar first return and ground points were gridded with a spatial Figure 1. (a) The study region, consisting of coastal young and mature broadleaf forest plantations (b) located in the subtropical region of northeast Taiwan (the star). The background image of (a) acquired in November 2015 was from Google Earth. (c) Monthly precipitation (gray-colored bars, the left y-axis) and air temperatures (the black line, the secondary y-axis) records of 1990-2021 for the study site during the observation period.

UAV Data Acquisition
UAV lidar and UAV photo data were acquired on 20 July 2020 using a LiAir V scanning system (GreenValley International, Berkeley, CA, USA) and Sony α7 RII (Sony Group Corporation, Tokyo, Japan), respectively, carried by a Matrice M600 pro UAV (Da-Jiang Innovations, Shenzhen, China) at an altitude of 180 m a.g.l. The average point density was 11.2 points m −2 . The UAV lidar first return and ground points were gridded with a spatial resolution of 1 m using adaptive kriging (SCOP++, Department of Geodesy and Geoinformation, Vienna, Austria) to generate the DSM and DTM, respectively. We then generated the CHM (CHM lidar ) by subtracting the DTM (from UAV lidar ) from the lidar DSM. For spatial continuity of the data, 80% forward and side overlaps were set before the flight [45]; the data were georegistered by referencing to nine ground control points using a handheld real-time kinematic (10 cm accuracy in the real-time kinematic (RTK) mode under the open-sky condition) global positioning system (GPS) (GeoXH, Trimble Inc., Sunnyvale, CA, USA). The UAV photo data were acquired using 80% forward and side overlaps. Point-cloud data from the UAV photo were interpolated and resampled to generate a 10 cm DSM using Pix4D (Pix4D SA, Lausanne, Switzerland). Mean x, y and z errors were established to be ≤3 cm by referring to the ground control points. The DSM derived from the UAV photo was further aggregated to a 1 m spatial resolution by a nearest neighbor interpolation method using a geographical information system (QGIS v. 3.4.4, http://www.qgis.org, the last accessed date 28 April 2021), and the CHM photo was also derived by subtracting the DTM from the DSM derived from the UAV photo . Figure 2 shows the flowchart of the study. resolution of 1 m using adaptive kriging (SCOP++, Department of Geodesy and Geoinformation, Vienna, Austria) to generate the DSM and DTM, respectively. We then generated the CHM (CHMlidar) by subtracting the DTM (from UAVlidar) from the lidar DSM. For spatial continuity of the data, 80% forward and side overlaps were set before the flight [45]; the data were georegistered by referencing to nine ground control points using a handheld real-time kinematic (10 cm accuracy in the real-time kinematic (RTK) mode under the open-sky condition) global positioning system (GPS) (GeoXH, Trimble Inc., CA, USA). The UAVphoto data were acquired using 80% forward and side overlaps. Point-cloud data from the UAVphoto were interpolated and resampled to generate a 10 cm DSM using Pix4D (Pix4D SA, Lausanne, Switzerland). Mean x, y and z errors were established to be ≤3 cm by referring to the ground control points. The DSM derived from the UAVphoto was further aggregated to a 1 m spatial resolution by a nearest neighbor interpolation method using a geographical information system (QGIS v. 3.4.4, http://www.qgis.org, the last accessed date 28 April 2021), and the CHMphoto was also derived by subtracting the DTM from the DSM derived from the UAVphoto. Figure 2 shows the flowchart of the study.

Gap Detection
Gaps may be defined as canopy openings reaching within 2 m of the ground [15], with height thresholds being relative to the height of the canopy surrounding a gap. To determine canopy gaps, we defined a height class and took a horizontal cross-section of the CHM smaller than that height threshold, and then recorded agglomerations of empty pixels surrounded by the full pixels [46]. We repeated the analysis for a range of height thresholds up to the maximal canopy height with incremental 2 m intervals. We removed gap size < 5 m 2 and > 2 ha by referring to White et al. [6]. Finally, we applied both criteria To assess gap detection performance

Gap Detection
Gaps may be defined as canopy openings reaching within 2 m of the ground [15], with height thresholds being relative to the height of the canopy surrounding a gap. To determine canopy gaps, we defined a height class and took a horizontal cross-section of the CHM smaller than that height threshold, and then recorded agglomerations of empty pixels surrounded by the full pixels [46]. We repeated the analysis for a range of height thresholds up to the maximal canopy height with incremental 2 m intervals. We removed gap size < 5 m 2 and >2 ha by referring to White et al. [6]. Finally, we applied both criteria on CHM lidar and CHM photo . We note that the primary objective of this study is to only develop ideal toolsets for forest gap quantification. Therefore, only a standard gap delineation approach was applied in this study to make a reasonable comparison.
To assess gap detection performance, we randomly sampled 30 gaps each for young and mature stands (n = 60) in July 2021 by referring to both GAP lidar and GAP photo , and georeferenced those gaps by using a handheld GPS. In addition, we estimated the size of each gap by referring to Runkle [47] and Yao et al. [19], assuming the shape of the gaps was elliptical, which should be appropriate for the application [48]. We then calculated the area (Equation (1)) by measuring the longest axis (L) and the one (W) perpendicular to L.
Gap area = πLW/4 (1) We compared the size differences (e.g., root mean squared error (RMSE)) of field observation and GAP lidar and GAP photo . We note that, to our knowledge, there were no disturbances (e.g., typhoons or forest management) between the times of UAV and field data acquisition. In addition, both the UAV and field campaigns were conducted in July. Therefore, the ramifications of the time difference between UAV and field sampling should be minimal.

Modeling Canopy Gap Distribution
We utilized the zeta and Weibull distributions to model the gap characteristics of coastal forests in the subtropical zone of northeastern Taiwan. The zeta distribution provides a summary of the frequency that is suitable for characterizing the distribution of landscape-level gap area [46,49] (Equation (2)): where the denominator is the Pareto distribution in a discrete power law probability density function. We employed the maximum likelihood to estimate λ for the zeta distribution [50]. The relationship (Equation (2)) becomes linear with a negative slope λ after transforming it to log-log space. The λ values usually vary from 1.0 to 3.0 for forests, and a greater value (e.g., >2) indicates more small gaps present in a forest (high-growth-low-mortality dynamics) and vice versa (mortality of large canopy) [4,36,51]. We derived λs for gaps derived from GAP lidar and GAP photo . We calculated each area of the gap to have their canopy size and frequency, and then used derived parameters λ and k to fit a zeta distribution (Equation (2)) by referring to Asner et al. [4]. The Weibull distribution function has been commonly applied for fitting multishape distributions because of its flexibility in characterizing data profiles [52][53][54][55]. This function can also be used to model the probability of an increasing, decreasing or stable trend. The two-parameter Weibull distribution is suitable for modeling phenomena with a monotonic decrease trend, with its probability density function for gap size given by: where f (g) is the decrease probability trend of gap size, and β and θ are the shape and scale parameters of the distribution with positive values. The β parameter, known as the shape parameter (the slope of the Weibull probability plot), determines the shape form of the Weibull family of distributions that best fits the data. The θ parameter is the characteristic gap size, which is also known as the scale parameter. We specifically studied Weibull p50 (p50 hereafter), the probability for which the gap size is the median of the Weibull probability density function. This metric was selected since p50 indicates that the median gap size is the area at which half of the amount is smaller than the median. Finally, we assessed the areal size effect of gap distribution quantifications for the zeta and Weibull distributions by randomly selecting 1-10 ha areas with a 1 ha increment for 30 times for the gaps (GAP lidar and GAP photo ) of young and mature forest stands and investigating the variation of λ and p50 (also termed the sensitivity analysis). We note that the descriptive statistics (e.g., mean, SD, median, min, max, skewness, kurtosis) were utilized to describe or summarize the characteristics of canopy gap data distribution. We applied the Shapiro-Wilk normality test (W) to investigate normality of the datasets. If the dataset was rejected by the Shapiro-Wilk test (not normally distributed), we then used Dunn's test, which is a nonparametric pairwise multiple comparison procedure based on rank sums.

Canopy Height Model Characteristics
The means (±SD) of young forest stand CHMs were 1.2 ± 0.9 m for CHM lidar and 1.4 ± 1.1 m for CHM photo ; those of mature forest stands were 4.0 ± 2.4 m for CHM lidar and 4.6 ± 2.2 m for CHM photo (for examples, see Figure 3). None of the CHMs were normally distributed (p < 0.001) according to the Shapiro-Wilk normality test. According to Dunn's test for multiple comparisons, a nonparametric pairwise multiple comparisons method [56], there were significant differences (p ≤ 0.001) of median CHM lidar and CHM photo for both young and mature stands, but not for sensors of the same forest type (p = 0.97 and 0.58 for young and mature stands, respectively) ( Table 1). In addition, Dunn's test for multiple comparisons (Table 1) also demonstrated significant differences (p ≤ 0.001) of median CHM lidar and CHM photo for both forest types but not for sensors of the same forest type (p ≥ 0.12). (GAPlidar and GAPphoto) of young and mature forest stands and investigating the variation of λ and p50 (also termed the sensitivity analysis). We note that the descriptive statistics (e.g., mean, SD, median, min, max, skewness, kurtosis) were utilized to describe or summarize the characteristics of canopy gap data distribution. We applied the Shapiro-Wilk normality test (W) to investigate normality of the datasets. If the dataset was rejected by the Shapiro-Wilk test (not normally distributed), we then used Dunn's test, which is a nonparametric pairwise multiple comparison procedure based on rank sums.

Canopy Height Model Characteristics
The means (±SD) of young forest stand CHMs were 1.2 ± 0.9 m for CHMlidar and 1.4 ± 1.1 m for CHMphoto; those of mature forest stands were 4.0 ± 2.4 m for CHMlidar and 4.6 ± 2.2 m for CHMphoto (for examples, see Figure 3). None of the CHMs were normally distributed (p < 0.001) according to the Shapiro-Wilk normality test. According to Dunn's test for multiple comparisons, a nonparametric pairwise multiple comparisons method [56], there were significant differences (p ≤ 0.001) of median CHMlidar and CHMphoto for both young and mature stands, but not for sensors of the same forest type (p = 0.97 and 0.58 for young and mature stands, respectively) ( Table 1). In addition, Dunn's test for multiple comparisons (Table 1) also demonstrated significant differences (p ≤ 0.001) of median CHMlidar and CHMphoto for both forest types but not for sensors of the same forest type (p ≥ 0.12).   Our field observation showed that 25 (16.7% not-gap rate with gap sizes ≤ 5 m 2 ) and 26 (13.3% not-gap rate with gap sizes ≤ 5 m 2 ) gaps were found in young and mature stands. The mean (±SD) of L in Equation (1) of young and mature stands was 15.2 ± 25.3 m and 17.3 ± 18.1 m, respectively; the mean (±SD) W of young and mature stands was 6.7 ± 12.6 m and 6.7 ± 6.5 m, respectively. With these parameters, we calculated the gap sizes for young (mean ± SD = 220.1 ± 915.3 m 2 ) and mature (95.1 ± 402.1 m 2 ) stands ( Figure 4). Estimated errors (RMSE) were 145.2 m for GAP lidar and 256.2 m for GAP photo for the young stands, and 87.2 m for GAP lidar and 218.5 m for GAP photo for the old stands ( Figure 5). There were strong agreements (R 2 ≥ 0.94, p < 0.001) between ground and UAV measurements (Table 2), and the performance of UAV lidar was superior to UAV photo . Our field observation showed that 25 (16.7% not-gap rate with gap sizes ≤ 5 m 26 (13.3% not-gap rate with gap sizes ≤ 5 m 2 ) gaps were found in young and mature The mean (±SD) of L in Equation (1) of young and mature stands was 15.2 ± 25.3 17.3 ± 18.1 m, respectively; the mean (±SD) W of young and mature stands was 6.7 m and 6.7 ± 6.5 m, respectively. With these parameters, we calculated the gap si young (mean ±SD = 220.1 ± 915.3 m 2 ) and mature (95.1 ± 402.1 m 2 ) stands ( Figure 4 mated errors (RMSE) were 145.2 m for GAPlidar and 256.2 m for GAPphoto for the stands, and 87.2 m for GAPlidar and 218.5 m for GAPphoto for the old stands ( Figure 5) were strong agreements (R 2 ≥ 0.94, p < 0.001) between ground and UAV measure ( Table 2), and the performance of UAVlidar was superior to UAVphoto. Table 2. Comparisons ( Figure 5) of ground and UAV observations (GAPground = b0 + b1 GAPse models are significant (p < 0.001).

Sensor
Forest

Gap Characteristics
Total numbers of GAPlidar and GAPphoto derived from CHMlidar and CHMphoto v markedly (n = 154 for GAPlidar and 128 for GAPphoto for young stands; n = 748 for G and 165 for GAPphoto for mature stands) ( Table 3). The mean (±SD) sizes of GAPli young and mature stands were 1392.9 ± 4298.8 m 2 and 74.0 ± 311.9 m 2 , respectively; photo for young and mature stands was 491.3 ± 1778.7 m 2 and 65.9 ± 99.2 m 2 , respect Both GAPlidar and GAPphoto for young and mature stands were not normally distribu < 0.001, the Shapiro-Wilk normality test) ( Table 3).

Zeta and Weibull Distributions
The results of zeta distribution for each dataset types were similar, and λs of G for young and mature forest stands were 1.36 and 1.71, respectively; those of GAPph young and mature forest stands were 1.45 and 1.61 (Figure 6), respectively. The resu fitted Weibull distributions depicted that the shape parameters of young stands we for GAPlidar and 0.5 for GAPphoto, and the scale parameters were 426.3 m 2 for GAPlid 181.4 m 2 for GAPphoto. The shape parameters were 0.6 for GAPlidar and 0.9 for GAPphot the scale parameters were 41.0 m 2 for GAPlidar and 51.3 m 2 for GAPphoto in mature s   Figure 5) of ground and UAV observations (GAP ground = b 0 + b 1 GAP sensor ). All models are significant (p < 0.001).

Sensor
Forest

Gap Characteristics
Total numbers of GAP lidar and GAP photo derived from CHM lidar and CHM photo varied markedly (n = 154 for GAP lidar and 128 for GAP photo for young stands; n = 748 for GAP lidar and 165 for GAP photo for mature stands) ( Table 3). The mean (±SD) sizes of GAP lidar for young and mature stands were 1392.9 ± 4298.8 m 2 and 74.0 ± 311.9 m 2 , respectively; GAP photo for young and mature stands was 491.3 ± 1778.7 m 2 and 65.9 ± 99.2 m 2 , respectively. Both GAP lidar and GAP photo for young and mature stands were not normally distributed (p < 0.001, the Shapiro-Wilk normality test) ( Table 3).

Zeta and Weibull Distributions
The results of zeta distribution for each dataset types were similar, and λs of GAP lidar for young and mature forest stands were 1.36 and 1.71, respectively; those of GAP photo for young and mature forest stands were 1.45 and 1.61 (Figure 6), respectively. The results of fitted Weibull distributions depicted that the shape parameters of young stands were 0.3 for GAP lidar and 0.5 for GAP photo , and the scale parameters were 426.3 m 2 for GAP lidar and 181.4 m 2 for GAP photo . The shape parameters were 0.6 for GAP lidar and 0.9 for GAP photo , and the scale parameters were 41.0 m 2 for GAP lidar and 51.3 m 2 for GAP photo in mature stands (Table 4). In the young stands, the gap sizes of p50 were 185.7 m 2 for GAP lidar and 85.4 m 2 for GAP photo . In the mature stands, the gap sizes of p50 were from 23.1 m 2 for GAP lidar and 39.6 m 2 for GAP photo (Table 4).
x FOR PEER REVIEW 9 of 15

Gap Size Distributions across Spatial Extents
To assess the influence of gap size distribution and the analyzed area, we changed the spatial extents from 1 to 10 ha with a 1 ha increment to detect gap characters by fitting the zeta and Weibull distributions. The locations of each spatial extent (1-10 ha) were randomly selected for 30 times. The values of λ were from 1.59 (1 ha) to 1.35 (10 ha) in young stands derived from GAPlidar, and from 1.84 (1 ha) to 1.41 (10 ha) in young stands derived from GAPphoto (Figure 7a). The analysis of zeta distribution show that λs were from 1.83 (1 ha) to 1.71 (10 ha) in mature stand derived from GAPlidar and from 1.95 (1 ha) to 1.71 (10 ha) in mature stands derived from GAPphoto (Figure 7b). In the young stands, the gap sizes of p50 were from 394.7 m 2 (1 ha) to 185.4 m 2 (10 ha) for GAPlidar, and from 259.6 m 2 (1 ha) to 85.7 m 2 (10 ha) for GAPphoto (Figure 7c). In the mature stands, the gap sizes of p50 were from 29.1 m 2 (1 ha) to 23.3 m 2 (10 ha) for GAPlidar, and from 40.5 m 2 (1 ha) to 39.8 m 2 (10 ha) for GAPphoto (Figure 7d).  Table 4. Summary of the fitted Weibull distributions of canopy gaps detected by UAV lidar and photographs. Note that p50 is the median of the Weibull probability density function (median gap size is the area at which half of the amount is smaller than the median).

Sensor
Forest

Gap Size Distributions across Spatial Extents
To assess the influence of gap size distribution and the analyzed area, we changed the spatial extents from 1 to 10 ha with a 1 ha increment to detect gap characters by fitting the zeta and Weibull distributions. The locations of each spatial extent (1-10 ha) were randomly selected for 30 times. The values of λ were from 1.59 (1 ha) to 1.35 (10 ha) in young stands derived from GAP lidar , and from 1.84 (1 ha) to 1.41 (10 ha) in young stands derived from GAP photo (Figure 7a). The analysis of zeta distribution show that λs were from 1.83 (1 ha) to 1.71 (10 ha) in mature stand derived from GAP lidar and from 1.95 (1 ha) to 1.71 (10 ha) in mature stands derived from GAP photo (Figure 7b). In the young stands, the gap sizes of p50 were from 394.7 m 2 (1 ha) to 185.4 m 2 (10 ha) for GAP lidar , and from 259.6 m 2 (1 ha) to 85.7 m 2 (10 ha) for GAP photo (Figure 7c). In the mature stands, the gap sizes of p50 were from 29.1 m 2 (1 ha) to 23.3 m 2 (10 ha) for GAP lidar , and from 40.5 m 2 (1 ha) to 39.8 m 2 (10 ha) for GAP photo (Figure 7d).

Discussion
Canopy gap dynamics are pivotal metrics, which may indicate the conditions o ests from different perspectives including ecology, carbon sequestration and ma ment. In this study, we verified that UAVlidar with the Weibull distribution may be a timal approach to characterize canopy gaps of young and mature broadleaf plan forests in a coastal region of northeastern Taiwan. To our knowledge, this UAVlid analysis method has not been previously published. In this section, we deliberate th sibility of using UAVlidar and UAVphoto to quantify canopy gaps, and demonstrate th plications on canopy gap monitoring with the availability of a proper analytical too

Canopy Gap Delineation Using UAVlidar and UAVphoto
The performance of remotely sensed canopy gap detection is sensitive to s types, such as high spatial resolution satellite optical imagery [57], point cloud data airborne laser scanning [4,22] and terrestrial laser scanning [5]. In this study, we ut UAVlidar and UAVphoto to quantify and analyze canopy gaps in coastal young and m

Discussion
Canopy gap dynamics are pivotal metrics, which may indicate the conditions of forests from different perspectives including ecology, carbon sequestration and management. In this study, we verified that UAV lidar with the Weibull distribution may be an optimal approach to characterize canopy gaps of young and mature broadleaf plantation forests in a coastal region of northeastern Taiwan. To our knowledge, this UAV lidar gap analysis method has not been previously published. In this section, we deliberate the feasibility of using UAV lidar and UAV photo to quantify canopy gaps, and demonstrate the applications on canopy gap monitoring with the availability of a proper analytical tool.

Canopy Gap Delineation Using UAV lidar and UAV photo
The performance of remotely sensed canopy gap detection is sensitive to sensor types, such as high spatial resolution satellite optical imagery [57], point cloud data from airborne laser scanning [4,22] and terrestrial laser scanning [5]. In this study, we utilized UAV lidar and UAV photo to quantify and analyze canopy gaps in coastal young and mature broadleaf forest plantations in the subtropical zone of northeast Taiwan. A pronounced discrepancy in the detection of the number of canopy gaps was discovered; a similar result was also observed in White et al. [6]. Both methods detected a similar number of canopy gaps in young forest stands, and UAV lidar observed 22.4% more gaps (n = 37) in 23.3 ha (Table 3). Since there was a strong agreement between field, GAP lidar and GAP photo (Figure 5), we conclude that both UAV lidar and UAV photo are suitable in detecting canopy gaps in young coastal broadleaf plantations.
On the other hand, UAV lidar can detect almost five times more gaps than UAV photo (Table 3), especially for small gaps in mature forest stands. There was strong agreement between GAP lidar and field observation (RMSE = 87.2 m 2 ), but not for GAP photo (RMSE = 218.5 m 2 ). This suggests that CHM photo derived from UAV photo using a standard approach may not be effective to map gaps in coastal mature broadleaf forests (also see [6]). The efficacy of UAV photo is strongly related to the intensity of the ambient light in the visual region [58], which makes it challenging to delineate canopy gaps of certain UAV photo view angles surrounded by shadow. Structure-from-motion (SfM) photogrammetry point-cloudderived DSM for dense canopy might lead to a continuous surface between several tree canopies. Due to the presence of dense canopies in the mature forest plantations, the DSM from UAV photo was unable to provide sufficient vertical points, while the DSM of UAV lidar did. Acquiring points in shaded canopy areas may be difficult for UAV photo , which may significantly hinder its ability to detect small canopy gaps [6,59].
Lidar technology is known for being effective in delineating a multilayer canopy structure, especially for mature forests with dense canopies (23,26) (Figure 3). In most cases, point cloud data acquired from UAV lidar were greater than those from airborne laser scanning (e.g., 10+ pts m −2 vs. 1-5 pts m −2 ) [43,44,[60][61][62], making it an ideal tool for mapping gaps in dense forests. Although the UAV photo is known to be cost effective for forest mapping [26] (in this case, costs for UAV lidar and UAV photo were USD 7000 and 3000, respectively), it may not be feasible to delineate the vertical profile in coastal mature forests with dense canopies.

Canopy Gap Structure Status
Gap distribution may reflect the condition of a forest [63]; the power-law zeta distributions (λs) have a narrow range of values across different sites in forests [4,6,35,37]. In this study, we utilized λs to analyze GAP lidar and GAP photo of young and mature coastal forest plantations. We found that λs fell into a narrow range (1.36-1.71) regardless of the forest types ( Figure 6). According to the synthesis by Jucker [37], λs follow the same U-shaped pattern with canopy height and converge on relatively similar minimum values at multiple sites, therefore limiting the use of the zeta method in characterizing gap-size frequency distributions. These λs indicate that both young and mature forests were dominated by large gaps, perhaps due to high forest mortality [4]. This high mortality may be attributed to the periodic disturbance caused by summer tropical cyclones in the region (e.g., from June to October) [64]. This is also in agreement with Fisher et al. [36]. We also found that λs were insensitive to forest maturity (young vs. mature forest plantations) even with significant differences of GAP lidar and GAP photo ( Figure 6). The λ appears to converge on a narrow range of values across differences in forest structure, climate and disturbance history, which may limit its use for inferring the characteristics that shape the canopy structure dynamics of forests [37]. Therefore, we conclude that a power-law zeta distribution may not be feasible to monitor canopy gap variation of coastal plantation forests.
The other approach that we utilized to analyze canopy gaps was the Weibull distribution, forming the distributions of GAP lidar and GAP photo with the shape and scale parameters ( Table 4). The shape parameters were all <1 (an exponential distribution), ranging from 0.3 to 0.9 (Table 4), indicating the decreasing probability increasing with gap sizes [65]. Our result showed that a greater gap size (the young stands) may yield a small shape parameter. The values of the scale parameter derived from GAP lidar were greater than those of GAP photo . In the young stands, the values of the Weibull distribution scale parameters derived from GAP lidar and GAP photo were 426.3 and 181.4 m 2 , respectively; both indicate the presence of large gaps in the young broadleaf forest plantations. Therefore, the scale parameter may clearly distinguish the difference of GAP lidar and GAP photo with the Weibull distribution.

Effects of Detected Areas
Stability of λs of power-law zeta distribution across spatial scales has been rarely investigated. In this study, we found that λs were stabilized after the spatial extent was >2 ha (3 ha) for both GAP lidar and GAP photo in the young (mature) stands. In general, λs were more stable in the young stands than the mature ones due to the presence of several large gaps with relatively few small gaps (Figure 7a,b).
For the young stands, values of p50 GAP lidar consistently decreased until the spatial extent was ≥8 ha. On the other hand, those of GAP photo stabilized when the spatial extents were ≥3 ha (Figure 7c). GAP lidar contained both small and large canopy gaps with greater variation (Table 3), as a consequence of high Weibull scale parameters through a range of spatial extents (Figure 7). Our results indicated that the domination numbers of gap size will affect the distribution, causing an increase of detected areas. The number of small gap size will result in the instability of p50. The method of UAV photo is stabilized in small detected areas due to a lack of detection on the small gaps (Table 3), possibly due to the potential errors caused by shadow and data overlapping. Contrarily, UAV lidar was able to detect small canopy gaps and is therefore stabilized at a larger spatial extent. For mature stands, the values of p50 for both GAP lidar and GAP photo were insensitive to spatial extents (Figure 7d), and UAV lidar was able to detect small canopy gaps for different spatial extents. The value of p50 was stabilized since small canopy gaps were dominant in mature forest stands. The results demonstrate that p50 may be applicable to assess characteristics of gaps in mature forest plantations regardless of spatial extents of the monitored region. Finally, on a side note, the sensitivity analysis ( Figure 7) implies that the spatial extent (60.7 ha) of this study is sufficient for the application since all canopy gap metrics were stabilized before reaching the areal size of 10 ha.

Conclusions
This study assessed the feasibility of UAV lidar and UAV photo to map canopy gaps of young and mature broadleaf plantation forests in a coastal region of northeastern Taiwan, and tested the feasibility of using different mathematical functions (power-law zeta and Weibull distributions) to characterize canopy gaps of forest stands. We found that both UAV lidar and UAV photo may be able to quantify gaps in young plantations. However, only UAV lidar is able to thoroughly delineate gaps in mature plantations with a dense canopy layer. Lidar is able to detect small canopy gaps due to the physical nature of the instrument for better quantification of forest vertical profiles and insensitive to canopy shadow. By referring to the canopy gap analysis conducted in this study, we conclude that the Weibull distribution is a robust tool for coastal canopy gap monitoring. The proposed approach (UAV lidar with the Weibull distribution) may permit frequent monitoring of forest structure dynamics, which is particularly crucial in the era of climate change.

Informed Consent Statement: Not applicable.
Data Availability Statement: Not applicable.