Study areas
This study was conducted in two areas: the Heilihe town of Inner Mongolia and the Dahebei town of Liaoning Province (Fig. 1). Due to the different introduction times of RTB from Shanxi Province to the above two regions, the initial occurrence times of RTB in the Dahebei area and Heilihe area were 2014 and 2017, respectively. According to the unmanned aerial vehicle (UAV) data of these two areas and the four grading criteria proposed by White et al. (2005), the Heilihe area was determined to be in the early stage of RTB outbreak with 1%-5% trees having a red crown, while the Dahebei area was in the middle stage of RTB outbreak with 5%-20% of trees showing a red crown.
The Heilihe town covers approximately 531 km2, rising from 750 m to 1200 m above sea level. The mean annual precipitation and temperature of this area are 470 mm and 6 °C, respectively (Heilihe Forestry Station). The Dahebei town covers an area of 176 km2 and its elevation ranges from 428 to 1018 m above sea level. Similarly, the mean annual precipitation and temperature of this area are 450.9 mm and 8.6 °C, respectively (China Meteorological Data Service Center, 2021). The two areas are dominated by Chinese pine pure forests. Besides, larch (Larix principis-rupprechtii) (not a host tree of RTB) and some broad-leaved tree species also grow in the study areas. The remaining landscapes are characterized as grasslands, agricultural lands, and urban areas.
Stand selection and field sampling
Area in the early stage of RTB invasion
Once forest stands are attacked by RTB, indicators such as pitch tubes and boring materials can be used to inspect whether the trees have been infected. In early August of 2018 and 2019, 79 sample plots (30 m × 30 m), in which the trees showed a continuum of damage degree caused by RTB, were randomly selected in the study area (Fig. 2). The sample plots were located at least 500 m from each other to reduce spatial autocorrelation of landscape elements. The coordinates of these sample plots were recorded using a GPS receiver (Garmin eTrex 309x, Beijing, China) with a precision of < 3 m. In each sample plot, we recorded the diameter at breast height (DBH), species name, and status (infected or uninfected) for each tree with a DBH > 8 cm. We used the canopy projection method to calculate canopy closure (Bunnell & Vales, 2011). Mineral soil (20 cm depth) was collected and composited at the center and four corners of the sample plots. Then, soil organic matter and total nitrogen contents were calculated by the potassium dichromate oxidation and Kjeldahl methods, respectively (Bao & Jiang, 2013). As for topographical characteristics, we recorded the elevation, slope, and aspect for each sample plot. In the subsequent data processing, we transformed the aspect to a south-west-ness index to express the stand environment in sun exposure or dryness (Beers et al., 1966). We recorded the number of RTB entrance holes in the sample plots to represent the damage degree, which was used as the response variable in model construction statistics.
Area in the middle stage of RTB outbreak
In September 2018, we observed that all of the stands in the Dahebei research area were suffered large-scale RTB outbreak through field investigation. The damage in the middle stage of RTB outbreak can be directly detected using UAV images. Therefore, we used the DJI Inspire 2 drone (DJI, Shenzhen, China) to collect RGB images of the whole stand. We obtained 24 synthetic UAV images, each covering more than 20 hectares. We delineated 32 (30 m × 30 m) sample plots (Fig. 2) on UAV images and the distance between the sample plots was more than 500 m as mentioned earlier.
The mean DBH and canopy density of the sample plots were obtained from the National Forest Resources Intelligent Management Platform (2021) and corresponding UAV images. We used digital elevation model (DEM) (Aster GDEM 30 m resolution, 2021) to obtain the elevation, aspect, and slope data of the sample plots. We counted the numbers of damaged Chinese pines (yellow, red, and grey crowns in UAV images) in these sample plots to represent different degrees of RTB damage and took them as the response variable for areas in the middle stage of RTB outbreak.
The plots showing different stages of RTB outbreak were selected from areas without external interferences, such as fire and pest management that can promote or inhibit the occurrence of bark beetles (Agne et al., 2016; Mezei et al., 2017). The variables related to the stand characteristics, topographical characteristics, and soil properties of sample plots were treated as stand-level variables (Mezei et al., 2014). In addition, through remote sensing imagery of Gaofen-2 (China Centre for Resources Satellite Data and Application, 2020) and field investigation, the forest coverage of the buffer with a radius of 1000 m around the sample plots ranges from 10% to 75% to ensure a reasonable distribution of the landscape-level variables (Wang et al., 2019).
Landscape-level variables
We selected the landscape metrics that are most likely to affect the bark beetle outbreak according to the previous studies (Bone et al., 2013; Simard et al., 2012; Wang et al., 2019). We measured landscape composition using the proportion of Chinese pines (PLAND) and Shannon's diversity index (SHDI) calculated from all land cover types. We used mean shape index (SHAPE_MN) and cohesion index (COHESION) to quantify landscape configuration.
Landscape-level variables were quantified in the area around each sample plot at three radii (250, 500, and 1000 m) based on high‐resolution remote sensing imagery, which were acquired from the Gaofen-2 satellite covering all study areas in September 2017 and June 2018. Land cover types were classified into five categories: Chinese pines, larches, broad-leaved trees, grasslands, farmlands, and residence communities. Additional details such as the classification method and map of the remote sensing images are available in Supporting Information (Fig S1, Table S1)
Model selection
Prior to model selection, we measured multicollinearity among all variables, and COHESION and soil organic matter were removed (see methods in Supplementary S2). The final set of variables is listed in Table 1. Next, to evaluate the response of RTB outbreak to stand- and landscape-level variables, we constructed generalized linear mixed models (GLMMs) for the area in the early stage of RTB outbreak and generalized linear models (GLMs) for the area in the middle stage of RTB outbreak. The year was used as random effect. Since overdispersion was observed in our data, we used a negative binomial error distribution that would provide better parameter estimates than the Poisson distribution (Militino, 2010).
Models were ranked based on corrected Akaike's information criterion (AICc) adjusted for small-sized samples (Krasnov et al., 2019). We selected the top model set with a ∆AICc (i.e., AICc - AICcmin) < 2, which is considered as effective. Model averaging was performed to produce Akaike weights (ω) and model-averaged partial regression coefficients for each variable. The relative importance value of each variable was quantified by the sum of the Akaike weights for each model in which the variable appeared. To approximate a normal distribution and enhance model stability, variables were log(x+1)-transformed as needed prior to model fitting (Schmiedel et al., 2015). Finally, we calculated the R2 of the best model (with the minimum AICc value) to observe the goodness of fit.
Variables found to be the most important in the early and middle stages of RTB outbreak (aspect, canopy density, and PLAND) were grouped if they were in the same direction (for aspect) or at equidistant intervals (for canopy density and PLAND), and a nonparametric Kruskal–Wallis test with post hoc was conducted for each variable. While the model located the variables that significantly affected RTB outbreak, the purpose of this analysis was to elucidate which groups of variables differed from the others, and to help forest managers effectively monitor and control RTB outbreak. All statistical analyses were performed in the R statistical programming environment (R Core Team, 2020) including the packages arm (Gelman, 2008), lme4 (Krasnov et al., 2019), MuMIn (Bartoń, 2019; Lukacs et al., 2009), piecewiseSEM (Nakagawa et al., 2017), and glmmTMB (Brooks et al., 2017).
Table 1 List of explanatory variables used in model analysis of areas in different stages of RTB invasion
Variable category
|
Variable
|
Scale of
measurement
|
Data source
|
Variable description
|
Area in the early stage of RTB outbreak
|
Area in the middle stage of RTB outbreak
|
Stand-level factors
|
DBH (cm)
|
Stand
|
FI
|
NFRD
|
Stand mean DBH
|
Canopy density (%)
|
Stand
|
FI
|
NFRD
|
Degree of canopy closure
|
Slope (°)
|
Stand
|
FI
|
DEM 30 m
|
Relief slope
|
Aspect (°)
|
Stand
|
FI
|
DEM 30 m
|
Stand orientation
|
Elevation (m)
|
Stand
|
FI
|
DEM 30 m
|
Stand elevation
|
Total N (g/kg)
|
Stand
|
LM
|
—
|
Total nitrogen of soil
|
Landscape-level factors
|
SHAPE_MN
|
Landscape
|
FC
|
FC
|
Shape complexity of host patches
|
PLAND
|
Landscape
|
FC
|
FC
|
Connectivity of host patches
|
SHDI
|
Landscape
|
FC
|
FC
|
Diversity of landscape types
|
FI, field investigation; NFRD, national forest resource; LM, laboratory measurement; FC, forest classification; — means there is no data.