The research area is located in the Maoer Mountain Experimental Forest Farm in Shangzhi City, Heilongjiang province, China (45°15′–45°29′ N, 127°23′–127°43′ E), with a total area of 2.6 × 103 hm2 (Fig. 1). The original forest in this area began to suffer serious damage during the construction of the Middle East Railway in 1906. Bald and barren hills, shrubs, swamps and damaged remnant communities were predominant in this area until 1949. After the founding of New China, all kinds of secondary forest communities recovered and developed naturally under the influence of reasonable management,, forming a typical natural secondary forest in the east of northeastern China. At the time of writing, a typical secondary forest has formed in this area. The forest types are hard broad-leaved forest, soft broad-leaved forest and mixed coniferous and broad-leaved forest. The main broad-leaved species include Mongolian Oak (Quercus mongolica Fisch), Aspen (Populous davidiana Dode) and Manchurian walnut (Juglans mandshurica Maxim).
LiDAR data were acquired using a LiteMapper 5600 LiDAR system (Riegl, Horn, Austria), adopting an LMS-Q560 (Riegl, Horn, Austria) as the laser scanner and a Yun-10 aircraft (Harbin Aircraft Manufacturing Company, Haerbin, China) as the airborne platform. Data collections were made on 14 and 15 September 2015 (deciduous season). When collecting data, the weather is clear and cloudless, which has no impact on lidar data collection. The LiteMapper 5600 system integrates laser ranging, global positioning system (GPS) and inertial navigation system (INS), including a single narrow band laser and a receiving system. The working wavelength of the laser was 1550 nm, the divergence angle of the laser beam was 0.5 mrad the waveform data recording interval was 1ns. The flight altitude was 1.2km, and the single point density was about 2 point/m2. After taking into account multiple echoes and repeated coverage, the maximum point cloud density was found to be more than 10 point/m2 and the average point cloud density to be about 3.7 point/m2. The spatial resolution of the digital orthophoto map (DOM) acquired by the CCD (Charge-Coupled Device) camera simultaneously is 20 cm.
Using the adaptive TIN model filtering of TerraScan software (Terrasolid, Heisinki, Finland), the original LiDAR point cloud data were divided into ground point cloud and non-ground point cloud. The ordinary Kriging interpolation method was used to interpolate the ground point cloud to generate the digital elevation model (DEM). The inverse distance weighting method was used to interpolate the first echo point to generate the digital surface model (DSM) (Guo et al. 2010). Compared with the measured elevation of the differential GPS, the obtained LiDAR data has an elevation accuracy better than 0.3 m and a plane accuracy better than 0.5 m. The CHM was obtained from the reduction of DSM and DEM data (Fig. 2). The spatial resolution was 1 m, and the data were the floating point type (32 bits). DOM and CHM data were in TIFF format, using the Xi’an 80 geographic coordinate system and Gauss Kruger 3-degree belt projection coordinate system.
According to Katarzyna et al. (2016), forest gaps are defined as areas with a height of less than 3m and an area of 5-1000 m2 on the ground. Because the difference between the inner microenvironment of less than 5 m2 and that of under forest is not significant, the inner microenvironment of more than 1000 m2 is similar to that of open space. We adopted a stratified sampling method to select the complete age groups; that is, the types of stand included young forest, middle-aged forest, near-mature forest, mature forest and over-mature forest. There were 75 sub-compartments (Table 1): stands of Aspen (7 sub-compartments), and stands of mixed forest with Mongolian Oak (41 sub-compartments) and Manchurian walnut (27 sub-compartments) as the main tree species. We used aerial images of the 75 sub-compartments as the base map, and manual digitization method was used to extract forest gaps. The area of each gap (m2) was calculated, and gaps less than 5 m2 and greater than 1000 m2 were excluded. ArcMap10.5 software (Redlands, California, America) was used to extract the grid with height of 0-3 m from CHM data and transformed it into vector graphics. Gap polygons were obtained and the geometric data were calculated as the basic parameters for quantification.1343 forest gaps of three forest types (Aspen, Mongolian Oak and Manchurian walnut) were extracted as sample of spatial feature quantification. This study did not adopt the method of full-automatic forest gap recognition, mainly because it is difficult to ensure the accuracy of the area and position of the forest gap needed for the quantification of the spatial characteristics of the forest gap, and there will be a lot of missing errors (Malahlela et al. 2014; Bonnet et al. 2015).
Quantitative method for analyzing spatial characteristics of forest gaps
Spatial characteristics of a single gap
The spatial characteristics of a single gap are quantified by three parameters: the size and shape of the gap and the spatial heterogeneity within the gap (Denslow et al. 1990; Runkle 1991). In this study, the area (m2) and GSCI of a single gap were used to quantify the size and shape of the gap, respectively (Blackburn and Milton 1996). The method of calculating GSCI is shown in Formula 1. The shape complexity index value of a perfectly circular gap is 1. The more complex the shape of the gap is, the higher the value of the shape complexity index is. The spatial heterogeneity in forest gaps is quantified by the GHD index. The Shannon Formula (Formula 2) is used to calculate the GHD in forest gaps (Zenner and Hibbs 2000). In this study, we used Arcmap10.5 software (Redlands, California, America) to reclassify the CHM grids corresponding to forest gaps; that is [0,0.5], (0.5,1], (1,1.5], (1.5,2], (2,2.5], (2.5,3]. The percentage of the number of grids at six height levels was obtained, and the GHD index was calculated to quantify the height diversity of vegetation in the forest gaps (Formulas 1 and 2):
where P is the perimeter of a single forest gap, A is the area of a single forest gap, is the proportion of grid elements in the i-th height class and N is the number of grades into which the data are divided.
Quantification of spatial characteristics between forest gaps
The spatial distribution of forest gaps is the spatial distribution pattern of forest gaps. The boundary of sub-compartments in the study area was irregular. If K function (Ripley 1976, 1977) is used to quantify the pattern distribution, complex calculation will be produced. So in this study, a modified nearest neighbor index method, Clark–Evans (CE), proposed by Fidner (1995), was used to quantify the spatial distribution pattern of forest gaps. The CE index was obtained by dividing the average distance between the center of mass of each gap and its nearest neighbor with the expected value when the individuals in the forest gaps were randomly distributed (Formula 3). When CE = 1, it was considered that the forest gaps were randomly distributed. When CE < 1, distribution was clustered; and, when CE > 1, the gaps were evenly distributed. The deviation degree of and was tested by a normal distribution test. When the test was not significant (|u| < 1.96) (Formula 4), it was considered that the forest gaps in the sub-compartments were randomly distributed (Kint et al. 2000):
where CE is the Clark–Evans index, is the distance between the ith individual and the nearest neighbor; N is the number of forest gaps in the sub-compartments, A is the area of the sub-compartments and P is the perimeter of the sub-compartments.
Through the study of frequency distribution, the overall distribution of gap area, GSCI and GHD in the study area was obtained, and the correlation between area and GSCI, area and GHD was obtained by correlation analysis. One-way ANOVA was used to test whether there were significant differences in gap characteristics among three forest types (Aspen, Mongolian Oak and Manchurian walnut) and five age groups (young forest, middle-aged forest, near-mature forest, mature forest and over-mature forest). All analyses were realized by MATLAB R2019a software (Natick, Massachusetts, America).