Assessing Spatio-Temporal Variation and Associated Factors of Forest Fragmentation from Morphological Spatial Pattern Analysis and Geo-Detector Analyses: A Case Study of Xinyu City, Jiangxi Province of Eastern China

Abstract

In the context of economic boom and climate change, monitoring the spatio-temporal dynamics of forest fragmentation induced by disturbances and understanding its corresponding associated factors are critical for developing informed forest management strategies. In this study, based on multi-temporal Landsat images acquired from 1999 to 2020, a SVM classifier was first applied to produce high-accuracy land cover maps in Xinyu City. Next, morphological spatial pattern analysis (MSPA) was implemented to characterize the spatio-temporal patterns of forest fragmentation by producing maps of seven fragmentation components, including the core, islet, perforation, edge, bridge, loop, and branch. Then, both natural and human factors responsible for the observed forest fragmentation dynamics were analyzed using the geo-detector model (GDM). The results showed that over the past two decades, Xinyu City experienced a process of significant forest area loss and exacerbating forest fragmentation. The forest area decreased from 1597.35 km² in 1999 to 1372.05 km² in 2020. The areal ratio of core patches decreased by 8.49%, and the areal ratio of edge patches increased by 5.98%. Spatially, the trend of forest fragmentation exhibited a progressive increase from the southern and northern regions towards the central and eastern areas. Large-scale forest core patches were primarily concentrated in the northwestern and southwestern regions, while smaller core patches were found in the eastern and central areas. Notably, human activities, such as distance from the roads and land use diversity, were identified as significantly associated with forest fragmentation. The interaction effect of these factors had a greater impact on forest fragmentation than their individual contributions. In conclusion, Xinyu City possesses the potential to further alleviate forest fragmentation by employing the regional differentiation development strategies: (1) intensive development in the northwest and southern regions; (2) high-density development in the western, northwestern, and southern regions, and (3) conservation development in the southwest, northeast, and east-central regions, thus aligning with the path of local social advancement.

1. Introduction

Human activities have been swiftly transforming land cover patterns and their spatial arrangement [1]. Forest ecosystems, in particular, are profoundly affected by processes such as deforestation, excessive development, and land reclamation [2,3], leading to substantial changes in the spatial functionality of forest ecosystems, such as the eroded ecosystem services [4], and escalated threats to biodiversity [5]. Moreover, these alterations potentially hinder human society’s progress and sustainability.

Unsustainable or negative forest exploitation and land use pattern transformation have substantially altered the spatial distribution and composition of existing forests, resulting in a progressive process of once-large and continuous forest patches converting into smaller, isolated forest remnants, which is a worldwide phenomenon known as forest fragmentation [5,6]. Forest fragmentation induces significant consequences, including the spread of diseases [7], abrupt shifts in microclimates at forest edges [8], increased forest fire occurrences [9], diminished habitat for species survival [10], increased greenhouse gas emissions [11], and the loss of biodiversity [12]. Thus, forest fragmentation has emerged as a notable global environmental challenge. For example, research has indicated that roughly 70% of the world’s forests are directly or indirectly impacted by human activities, leading to forest degradation resulting from fragmentation [13]. Consequently, investigating the progression of forest fragmentation and its underlying factors holds immense importance for both forest ecosystems and human society to aid in implementing measures to mitigate forest fragmentation and promote harmonious development between humans and nature.

In previous studies, researchers usually employed traditional landscape indices, such as the number, shape, size, and connectivity of forest patches, to quantify the degree of forest fragmentation [14,15,16] via a landscape analysis software named Fragstats (v4.2) to calculate the values of these indices to assess the temporal trends or spatial pattern differences of forest fragmentation. For instance, Slattery and Fenner employed Fragstats to assess forest fragmentation in seven Brazilian cities, revealing that regions experiencing agricultural expansion exhibited more intricate shapes in their remaining forest patches and significant losses in their core forest patch areas. Conversely, areas affected by commercial deforestation tended to produce disconnected core patches within the forest landscape [17]. Bulut et al. employed spatial land metrics to analyze forest fragmentation on the eastern Black Sea using the number of fragments, average fragment size, minimum and maximum fragment sizes, and the area-weighted average figure index [18]. Zhang et al. used six indices, such as the area-weighted mean patch area and area-weighted shape index, to analyze forest fragmentation in Yingkou City [19]. These endeavors have enabled a quantitative overview of the overall forest fragmentation situation in different regions. However, they are limited by the spatial relationship between forests and non-forests and cannot visually represent forest fragmentation locations and severity [14], which severely limits the usability or practicality of their findings in guiding the development of forest restoration programs and anti-fragmentation management actions [20].

A crucial aspect of describing forest fragmentation involves accurately distinguishing between ‘internal’ and ‘external’ fragmentation [21]. In recent decades, numerous novel methods have been developed to describe forest fragmentation in a spatial explicit manner instead of just a numeric manner. For instance, Riitters et al. employed a moving window analysis of land cover maps to calculate the forest area density and forest connectivity indices to classify the central forest pixel of the window into different fragmentation components, such as the core, perforated, edge, and patch [22,23]. Additionally, Forman et al. considered that landscape fragmentation was a temporal course, involving five distinct and spatially explicit processes: perforation, dissection, subdivision, shrinkage, and attrition [24]. Based on this theory, Li and Yang, and Ren et al. all considered landscape fragmentation processes like perforation, subdivision, shrinkage, and attrition to detect spatially explicit processes of forest fragmentation, but Ren et al. further advanced other landscape restoration processes, including increment and expansion, to monitor changes in the forest landscape and ecological consequences [25,26]. In the meantime, Soille pioneered the method of morphological spatial pattern analysis (MSPA) [27]. This method can effectively classify structure in order to simplify image data by associating morphological features to specific shapes in raster images [14], and represents a clear distribution map of forest fragmentation [28]. The integration of mathematical morphology and image processing has brought about a surge of interest in utilizing morphological spatial pattern analysis (MSPA) methods for assessing forest fragmentation to compensate for the shortcomings of other methods [29,30,31,32]. For example, Yang et al. combined MSPA and the minimum resistance model (MCR) to construct the ecological network of “Ailaoshan-Wuliangshan” in Yunnan to provide suggestions for reducing the impact of national park construction on the surrounding environment [33]. Zhang et al. utilized the vegetation change tracking (VCT) algorithm in conjunction with MSPA to map forest cover and fragmentation patterns in Nanjing from 1987 to 2017 to explore the impact of urban expansion on forest cover and fragmentation dynamics [28]. Compared to the traditional landscape index-based methods, these existing model-based forest fragmentation mapping studies have had varying degrees of success in spatially locating fragmentation components to enable better targeting for protection and remediation efforts; however, they have rarely explored the associated factors of forest fragmentation from a spatial perspective.

The mechanisms driving forest fragmentation are complex and linked to multiple contributing factors [19,34,35]. The geo-detector model holds the capacity to depict the spatial distribution of geographical phenomena, elucidate the driving influence of various factors, and reveal the interplay between these factors. It boasts the added advantage of accommodating diverse data types and is free of linearity assumptions [36]. Currently, the application of the geo-detector model has been expanded into numerous domains, including ecology [37], disaster management [38], pollution prediction [39], and urban expansion [40]. However, it has received little attention in the evaluation of the difference in contribution to forest fragmentation caused by different factors.

Given that 85% of studies on forest fragmentation are located in the United States and Europe [41], and with the forests in southern China having experienced significant forest fragmentation changes between 2000 and 2020 [42], the major objective of the current work was to take Xinyu City of Jiangxi Province as the prototype study area to test the effectiveness of the proposed framework that integrates mathematical morphological analysis and geo-detector analysis in forest fragmentation characterization and spatial detection of different driving forces’ contribution to forest fragmentation. Moreover, we also expected to propose targeted strategies or recommendations to improve forest management to mitigate forest fragmentation and achieve sustainable forest management goals based on the outcomes of the analysis.

2. Materials and Methods

2.1. Study Area

Xinyu City is located in the central region of Jiangxi Province, China (27°33′–28°05′ N; 114°29′-115°24′ E) (Figure 1). Encompassing an area of 3178 km², the city extends 101.9 km from east to west and 65 km from north to south. As of 2020, the GDP of Xinyu city reached 100 billion RMB, with a population of 1,202,100. Xinyu City exhibits a subtropical monsoon-influenced climate characterized by four distinct seasons, an average annual temperature of 17.7 °C, and an annual precipitation of 1600 mm.

Notably, Xinyu has earned acclaim as a renowned forest city in China, abundant in forestry resources. It has 18 national-level forest villages and 31 provincial-level forest villages [43]. During the period of 2016–2020 only, Xinyu City has planted a total of 13,653 ha of forests. The forestry industry in Xinyu City has experienced a rapid development, primarily including industries such as tea oil production, flower and nursery stock cultivation, forest tourism, and timber processing, with a total output value of 12.19 billion RMB in 2020 [44].

2.2. Data Sources

The land cover data used in this analysis were derived from implementing a SVM classification based on the Landsat TM and OLI images acquired in 1999, 2006, 2013, and 2020. These Landsat images were directly downloaded from the USGS EROS Data Center (https://glovis.usgs.gov/app, accessed on 10 December 2021). Before placing our data order, we limited the image acquisition dates to the peak season of vegetation growth, e.g., from mid-June to mid-September for this mid-latitude region, to minimize the influence of vegetation phenology on classification accuracy.

According to the relevant studies [45,46,47], forest fragmentation is associated with multiple factors, including anthropogenic factors, such as agricultural and forestry development, expansion of construction land, irrational forest activities and forest management, property rights, and uncoordinated management policies, and natural factors, such as the climate and geography. Considering the availability and reliability of data, a total of six natural factors were identified, including the altitude (X₁), slope (X₂), aspect (X₃), annual rainfall (X₄), accumulated sunshine hours (X₅), and distance from rivers (X₆). Additionally, five human activity factors were selected, encompassing the distance from roads (X₇), the areal ratio of cropland (X₈), the areal ratio of construction land (X₉), distance from consumption centers (X₁₀), and land use diversity (X₁₁). These factors associated with fragmentation were compiled from various sources, including DEM data, climate data, infrastructure data, land cover data, and POI data. Subsequently, we recorded and categorized these data in accordance with pertinent industry standards and specifications.

Table 1. Description of the associated factors responsible for fragmentation involved in this analysis.

Factor	Resolution	Time	Classification Quantity	Unit	Source
Altitude (X1)	30 m	2018	6	m	https://earthexplorer.usgs.gov/, accessed on 14 December 2021
Slope (X2)	30 m	2018	4	°	https://earthexplorer.usgs.gov/, accessed on 14 December 2021
Aspect (X3)	30 m	2018	9	°	https://earthexplorer.usgs.gov/, accessed on 14 December 2021
Annual rainfall (X4)	30 m	2020	5	mm	https://data.cma.cn/, accessed on 15 December 2021
Accumulated sunshine hours (X5)	30 m	2020	4	hours	https://data.cma.cn/, accessed on 15 December 2021
Distance from rivers (X6)	1:20,000	2020	4	km	https://www.openstreetmap.org/, accessed on 14 December 2021
Distance from roads (X7)	1:20,000	2020	4	km	https://www.openstreetmap.org/, accessed on 14 December 2021
Areal ratio of cropland (X8)	30 m	2020	5	%	Calculated from Landsat images
Areal ratio of construction land (X9)	30 m	2020	5	%	Calculated from Landsat images
Distance from consumption centers (X10)	30 m	2020	4	km	http://www.bigemap.com/, accessed on 14 December 2021
Land use diversity (X11)	30 m	2020	5		Calculated from Landsat images

summarizes the basic information of these associated factors contributing to forest fragmentation in this work.

2.3. Method

2.3.1. Landsat Image Preprocessing

The preprocessing of the downloaded Landsat TM and OLI images mainly involved their atmospheric corrections to generate highly accurate surface reflectance images to support subsequent classifications of land cover. Landsat 8 Surface Reflectance data were generated from the Land Surface Reflectance Code (LaSRC). The LaSRC makes use of the coastal aerosol band to perform aerosol inversion tests, employs auxiliary climate data from MODIS, and uses a unique radiative transfer model. The USGS currently only provides LaSRC-based products within the U.S., but remote sensing desktop (RSD) software (version 3.2.1, Shenyang, China) introduced the core code of the LaSRC, which can be used for the atmospheric correction of Landsat 8 OLI data within China. According to the comparison between RSD products and USGS LaSRC products for the same dataset, the difference was relatively small, so we used RSD software for the atmospheric correction of Landsat 8 OLI images [48]. Landsat 5 TM Surface Reflectance data were generated using the Landsat Ecosystem Disturbance Adaptive Processing System (LEDAPS) algorithm (version 3.4.0, U.S. Geological Survey, Reston, VA, USA), a specialized software tool originally developed through a National Aeronautics and Space Administration (NASA). Making Earth System Data Records for Use in Research Environments (MEaSUREs) grant by the NASA Goddard Space Flight Center (GSFC) and the University of Maryland [49], based on the second simulation of a satellite signal in the solar spectrum (6S) radiative transfer models. The LEDAPS model takes the cloud-free atmosphere as an assumed case and integrates the absorption and scattering effects of water vapor, carbon dioxide, and aerosols. The gain and bias obtained from the header file are used to calculate the radiant brightness and reflectance at the upper boundary of the atmosphere in combination with the atmospheric parameters. Then, we obtained the surface reflectance image by interpolating the generated aerosol optical thickness (AOT) and the collected data, such as water vapor and water pressure. These two techniques are very mature and routinely implemented to generate highly accurate Landsat analysis ready data for surface change monitoring [50].

2.3.2. Land Cover Classification and Validation

We chose the support vector machine (SVM) classifier to produce land cover maps of Xinyu City. Numerous studies have demonstrated that SVMs have the advantages of coping with the problems of high-dimensional data, which means that the input dataset will be projected into a higher-dimensional feature space where the training samples will become linearly separable [51], possessing a strong generalization capability when solving small sample learning problems [52,53], and having high accuracy and robustness in numerous experiments [54] relative to classical statistical algorithms, like the minimum distance classifier (MD) and maximum likelihood classifier (ML), and other machine learning algorithms, such as the random forest classifier (RF) and classification and regression tree classifier (CART) [55,56,57], when dealing with classification issues. Before classifying the images, by referring to the existing land cover classification systems in China [58,59] and the local surveys in Xinyu City, the final classification scheme was defined as five classes, including the construction land, forest, water, cropland, and unused land. Next, based on the temporally corresponding Google Earth high-resolution maps or original Landsat images, we implemented an intensive visual interpretation to produce the training dataset consisting of 800 pixels (

Table 2. The number of the training and validation samples for land cover classification in 2020.

Land Cover Type	Training Samples	Validation Samples
Cropland	200	176
Forest	350	225
Water	45	25
Construction land	70	42
Unused land	35	32

and Figure 2) for the five classes. Table 2 shows the number of training samples for different classes in 2020. Considering the wide adoption of the Kanth–Thomas transform, NDVI, and NDBI for monitoring vegetation and urban development [59,60], after several comparative experiments, we calculated the (1) NDVI, (2) NDBI, and (3) brightness, greenness, and wetness components derived from the K–T transform, and (4) an optimized combination of green, near-infrared, and shortwave infrared bands (bands 2, 4, and 7 for Landsat 5 TM, and bands 3, 5, and 7 for Landsat 8 OLI) [61,62] as the input features for SVM classification. Since the training sample set in this study comprised non-linear low-dimensional space data, the radial basis function (RBF) was utilized to map to the high-dimensional space to obtain the optimal hyperplane. We chose to perform the tuning with a triple cross-validation with grid search, whereby multiple classes were classified for land cover mapping. The RBF function used in this study is formulated as follows:

$$K\left(x_i,x_j\right)=\mathrm{e}\mathrm{x}\mathrm{p}(-y {\parallel x_i-x_j\parallel }^2)$$

(1)

where $x_i$ and $x_j$ are the sample input features, ${\parallel x_i-x_j\parallel }^2$ is the square of the Euclidean distance of the samples, and $y$ is the free parameter.

Finally, 500 random sample points (Table 2 and Figure 2), as the validation dataset, were generated from the classification maps, and their corresponding actual classification attributes were determined by visually interpreting their corresponding Google Earth maps or original Landsat images; this was followed by a one-to-one comparison to derive the accuracy statistics, including the overall accuracy and kappa coefficient [28].

2.3.3. Characterization of Forest Fragmentation

We employed the MSPA method to analyze the spatial pattern of the forest landscape. MSPA relies on set theory and topology for image processing and implements morphological operators, including erosion, dilation, and opening/closing operations. These operations are instrumental in analyzing the geometry and connectivity of binary images [63,64]. Prior to the implementation of MSPA, the land cover map of Xinyu City was aggregated into two categories: forest and non-forest. Then, the Guidos toolbox (version 2.8, European Comission, Ispra, Italy) was used to set the edge width to one pixel and an adjacent eight-pixel neighborhood as the structure of a moving window for MSPA analysis. The Guidos toolbox provides a generic framework that is applicable to image analysis at any scale and for any kind of digital raster data [21]. The MSPA analysis finally classified each forest pixel into one of the seven fragmentation components, including the core, islet, perforation, edge, bridge, loop, and branch [65], as outlined in

Table 3. Forest fragmentation components and their corresponding ecological definitions involved in the MSPA analysis.

Fragmentation Component	Ecological Definition
Core	The larger habitat patches in the foreground pixel that can provide a larger habitat for species, which is of great significance for the protection of biodiversity and is the ecological source in the ecological network.
Islet	Isolated and broken small patches that are not connected to each other. The connectivity between patches is relatively low, and the possibility of internal material and energy exchange and transfer is relatively small.
Perforation	Non-forest land inside the ecological core area; this area does not have ecological benefits.
Edge	The outer edge of the forest pixel. The transition area between the forest core area and the non-forest area, which has an edge effect.
Bridge	This area is an ecological land connected to the core area, such as the regional corridor, which can promote the flow exchange of energy and material within the region.
Loop	The ecological corridor connecting the same core area is small in scale and low in connection with the surrounding natural patches.
Branch	The area with only one end connected to the core patch, which is the channel for species diffusion and energy exchange with the peripheral landscape.

.

2.3.4. Detection of Factors Associated with Forest Fragmentation

In order to accurately analyze the spatial heterogeneity of forest fragmentation, the optimal analysis scale or unit must be determined beforehand. We chose the coefficient of variation (CV) to obtain the optimal analysis unit by comparing five different scale fishnet grids (1 km × 1 km, 2 km × 2 km, 3 km × 3 km, 4 km × 4 km, and 5 km × 5 km) after calculating the degree of variation in forest fragmentation severity, which was defined as the areal ratio of the remaining six fragmentation components, excluding the core component within the analytical unit over the area of the forest (including all the seven fragmentation components defined in Table 3) in the analytical unit here. We extracted the area of different forest fragmentation components in each unit through the ‘zonal statistics as table’ tool in ArcGIS 10.8 and then derived the fragmentation severity in each analytical unit. Next, the mean and standard deviation of forest fragmentation severity for the entire study area were obtained to derive the coefficient of variation for each analytical scale, and the optimal analytical scale or unit, which gave the biggest CV, was accordingly determined. The larger the CV value, the greater the spatial variation characteristics of forest fragmentation severity. When the CV ≤ 10%, it shows a weak variability; when 10% < CV < 100%, it exhibits a moderate variability; and when the CV ≥ 100%, it shows a strong variability [66]. The calculation formula of CV was shown in Equation (1):

$$CV=\raisebox{1ex}{$S$}\!\left/ \!\raisebox{-1ex}{$U$}\right.\times 100\mathrm{\%}$$

(2)

where $CV$ is the coefficient of variation of forest fragmentation severity, $S$ is the standard deviation of forest fragmentation severity, and $U$ is the average of forest fragmentation severity in the study area.

Based on the relevant studies [67,68,69] and the experimental comparisons, we reclassified the forest fragmentation severity of each analytical unit into five grades: low fragmentation, above-low fragmentation, medium fragmentation, sub-high fragmentation, and high fragmentation, by following the criteria of natural breakpoint reclassification.

In this study, the geo-detector model was applied to evaluate the spatially stratified heterogeneity of forest fragmentation and to identify its related factors in Xinyu City [70,71]. The geo-detector model operates on the foundational assumption that when an independent variable (x) shares a similar spatial distribution with a dependent variable (y), the independent variable (x) possesses substantial explanatory power concerning the occurrence and progression of the dependent variable (y) [72]. In simpler terms, if a factor, such as the presence of edge patches in the forest fragmentation process, is influenced by a specific factor, like the distance from roads, then both factors should exhibit congruent spatial distributions. Based on the determined optimal scale or grid, we used the “zonal statistics as table” tool in ArcGIS 10.8 to count the values of the factors in each grid, and the parameter “statistics type” was set to “mean”. Then, we matched the results with the attribute table of the grids. Finally, we classified the data according to the classification quantity mentioned in Table 1 and placed all the classified data into the geo-detector to obtain the final results.

The geo-detector model consists of four primary modules: the risk detector, factor detector, interaction detector, and ecology detector. In this study, the factor detector and interaction detector were primarily utilized to investigate the principal-related factors impacting forest fragmentation in Xinyu City and to evaluate the interactive effects of different factors.

(1)

Factor detector: This module was employed to monitor the spatial heterogeneity of the forest fragmentation process and to gauge the explanatory power (q-value) of its associated factors on its spatial differentiation. By conducting significance tests on the mean value differences, it quantifies the explanatory power of various factors to assess their relative importance. The calculation formula of the q-value was shown in Equation (2):

$$q=1-\frac{{\sum }_{h=1}^L{N}_h{\delta }_h^2}{N{\delta }^2}$$

(3)

where $q$ is the explanatory power of each factor to the forest fragmentation process, and its value ranges from 0 to 1. $h=1$, …, $L$, is the classification number of the factor. $N_h$ and $N$ are the number of the analytical units in stratum $h$ and the whole study area, respectively. ${\delta }_h^2$ is the variance of forest fragmentation severity in stratum $h$, and ${\delta }^2$ is the variance of forest fragmentation severity in the entire study area. Specifically, a larger $q$-value indicates a higher similarity for the spatial distribution between the factors and forest fragmentation, that is, the value of $q$ means that factor $X$ explains $100\times q\%$ of $Y$ variance. In simpler terms, this means that the larger the value of $q$, the greater the impact of the factor ($X$) on forest fragmentation ($Y$).

(2)

Interaction detector: The interaction detector module was used to assess the degree, direction, and linear/non-linear relationship of the interaction between two independent factors potentially associated with forest fragmentation severity. By comparing the q-value of the individual factor’s impact on the forest fragmentation process with the q-value of the interaction between different independent factors, it can be ascertained whether the two independent factors collaborate to enhance or weaken the explanatory power of forest fragmentation or whether their effects on forest fragmentation are mutually independent. For example, in the case of a two-factor enhancement relationship, it holds that $\left\{q\right(X_1\cap X_2)>Max[q\left(X_1\right),q\left(X_2\right)]$}, while in the case of a non-linear enhancement relationship,$\left\{q\left(X_1\cap X_2\right)>\left[q\right(X_1)+q(X_2\left)\right]\right\}$.

3. Results

3.1. Land Cover Classifications and Validations

3.1.1. Accuracy Assessment

Table 4. The accuracy statistics of the land cover classifications derived from an independent validation.

	Year	1999	2006	2013	2020
Statistic
OA (%)		93	92.4	94	93.6
Kappa		0.898	0.877	0.914	0.903

summarizes the accuracy statistics of the land cover classifications during the period from 1999 to 2020. The detailed confusion matrices used to derive these accuracies were appended in Appendix A (

Table A1. Confusion matrix of the classification in 2020.

Land Cover Classification Map in 2020
		Cropland	Forest	Water	Construction Land	Unused Land	Total	Producer Accuracy
	Cropland	159	10	0	0	2	171	92.98%
Reference map	Forest	5	215	0	0	0	220	97.73%
	Water	1	0	23	0	0	24	95.83%
	Construction land	6	0	2	42	1	51	82.35%
	Unused land	5	0	0	0	29	34	85.29%
	Total	176	225	25	42	32	500
	User accuracy	90.34%	95.56%	92.00%	100.00%	90.63%
	Overall accuracy	93.60%		Kappa coefficient	0.903

,

Table A2. Confusion matrix of the classification in 2013.

Land Cover Classification Map in 2013
		Cropland	Forest	Water	Construction Land	Unused Land	Total	Producer Accuracy
	Cropland	150	3	1	2	1	157	95.54%
Reference map	Forest	6	200	0	0	0	206	97.09%
	Water	0	1	34	0	0	35	97.14%
	Construction land	3	0	2	56	2	63	88.89%
	Unused land	5	1	0	3	30	39	76.92%
	Total	164	205	37	61	33	500
	User accuracy	91.46%	97.56%	91.89%	91.80%	90.91%
	Overall accuracy	94.00%		Kappa coefficient	0.914

,

Table A3. Confusion matrix of the classification in 2006.

Land Cover Classification Map in 2006
		Cropland	Forest	Water	Construction Land	Unused Land	Total	Producer Accuracy
	Cropland	149	5	2	0	0	156	95.51%
Reference map	Forest	11	215	6	2	3	237	90.72%
	Water	2	0	40	0	0	42	95.24%
	Construction land	3	1	1	34	0	39	87.18%
	Unused land	1	0	0	1	24	26	92.31%
	Total	166	221	49	37	27	500
	User accuracy	89.75%	97.29%	81.63%	91.89%	88.89%
	Overall accuracy	92.40%		Kappa coefficient	0.887

and

Table A4. Confusion matrix of the classification in 1999.

Land Cover Classification Map in 1999
		Cropland	Forest	Water	Construction Land	Unused Land	Total	Producer Accuracy
	Cropland	164	4	1	1	0	170	96.47%
Reference map	Forest	9	200	2	0	0	211	94.79%
	Water	0	0	33	0	0	33	100.00%
	Construction land	9	0	0	40	3	52	76.92%
	Unused land	3	0	1	2	28	34	82.35%
	Total	186	204	37	43	31	500
	User accuracy	88.65%	98.04%	89.19%	93.02%	90.32%
	Overall accuracy	93.00%		Kappa coefficient	0.898

). The OAs of the classification results were all above 92%, and the Kappa coefficients were all above 0.87, indicating that all the classifications were reliable, and they could act as the input for the subsequent fragmentation analysis. The 2013 classifications had the highest accuracy, with an OA of 0.94 and a Kappa coefficient of 0.914, respectively, while the 2006 classifications gave the lowest accuracy, with an OA of 0.924 and a Kappa coefficient of 0.887, respectively.

3.1.2. Spatio-Temporal Changes in Land Cover Types

Figure 3 shows the spatio-temporal patterns of the land cover classes in Xinyu City from 1999 to 2020. Obviously, forest was always the major land cover type in this region, followed by cropland, water, construction land, and unused land in the four years. Forests were mainly distributed in the northwest and southwest regions of the study area, such as Mengshan Mountain in the north of Yushui District and Dagangshan Mountain in the southwest of Fenyi County. During the study period, there was a notable reduction in forest area, particularly in the northern and eastern regions of Yushui District and the central area of Fenyi County. Construction lands were mainly distributed in the urban area of Yushui District and the county seat of Fenyi County. During the period from 1999 to 2020, the scattered and small patches of construction land gradually increased, resulting in a concentrated distribution of construction land in the central and eastern regions. The distribution of unused land was similar to that of the construction land. Cropland was mainly distributed in the eastern and central regions, as well as around towns and villages. The distribution of cropland has undergone significant changes; for example, in the northern part of the city center, cropland gradually replaced the original forest during the study period. Water bodies were primarily located in the southwestern part of the study area, with a slight change in its area or size. Unused lands were widely distributed in the study area, mainly located between croplands and forests.

Table 5. Area changes of different land cover types in Xinyu City from 1999 to 2020 (km²).

	Year	1999	2006	2013	2020
Land Cover
Cropland (km2)		1251.96	1105.17	1153.11	1194.15
Forest (km2)		1597.35	1622.05	1511.32	1372.05
Water (km2)		98.02	135.40	87.80	94.67
Construction land (km2)		150.83	209.59	282.92	376.52
Unused land (km2)		61.93	87.88	124.94	122.70

shows the change statistics of different land cover types during the period from 1999 to 2020. It can be seen that forest area increased slightly from 1999 to 2006 and rapidly decreased after 2006, and a maximum area of 1622.05 km² in 2006 and a minimum area of 1372.05 km² in 2020 were observed. Construction land and unused land showed a substantial growth, almost doubled. Cropland was reduced from 1251.96 km² in 1999 to 1105.17 km² in 2006, then constantly rebounded slightly to 1194.15 km² in 2020. The overall area of water bodies remained relatively stable during the study period.

3.2. Spatio-Temporal Variations of Forest Fragmentation

Figure 4 shows that large-sized core patches were mainly distributed in the southwest and northwest regions, while the small-sized and medium-sized core patches were distributed in the eastern and central regions. Perforations were mainly distributed in the northwest and southwest regions, while bridges and loops were primarily located between medium-sized and small-sized core patches. The islets were distributed in the central and eastern regions, where the forest area was very small.

The forest fragmentation in the study area showed a trend of worsening between 1999 and 2020. Between 1999 and 2006, the overall characteristics of forest core patches remained relatively stable, but in proximity to urban areas, there was a decline in forest core areas, while the number of bridges and loops between core patches increased. From 2006 to 2013, the phenomenon of forest fragmentation extended from areas near urban centers towards the mountains in the northern and southern regions; there were more perforations connecting to large forest core patches in the southwest region, and more bridges, rings, and islets in the central region. From 2013 to 2020, forest fragmentation continued to spread, with the forest core patches in the central part of the study area having further shrank or even disappeared. The northern and southern parts of the study area suffered a significant decrease in the forest core areas, leading to the separation of core patches.

The temporal variation of forest fragmentation is outlined in

Table 6. Changes in the proportion of forest fragmentation types to forest area in Xinyu City from 1999 to 2020 (%).

	Year	1999	2006	2013	2020
Type
Core (%)		78.29	73.8	73.28	69.8
Islet (%)		1.38	1.82	1.83	1.65
Perforation (%)		3.2	3.85	3.48	3.13
Edge (%)		12.52	12.98	15.07	18.5
Loop (%)		1.38	2.14	1.68	1.71
Bridge (%)		1.08	1.94	1.33	1.47
Branch (%)		2.76	3.25	3.33	3.74
Total (%)		100	100	100	100

. During the past 20 years, cores were the dominant morphological spatial pattern, followed by edges, branches, and perforations, whereas the bridges were the least common. The core area proportion (the core area divided by the total forest area) continuously decreased from 78.29% in 1999 to 73.8% in 2006, then to 73.28% in 2013, and finally to 69.8% in 2020. Conversely, the trends of edges and branch patches differed from those of core patches during the study period, exhibiting a continuous upward trend in the proportion. The areal proportion of edges showed an apparent increase over time from 12.52% in 1999 to 18.5% in 2020. These facts indicated that an accelerated forest fragmentation process was witnessed in Xinyu City during the period from 1999 to 2020. In contrast, from 1999 to 2020, the perforation patches showed fluctuating changes; the ratio of perforation patches to forest areas fluctuated between 3.13% and 3.85%. The proportion of islets, bridges, and loops also fluctuated, but the maximum proportion was only 2.14% for loops in 2006.

3.3. The Spatio-Temporal Differentiation Law of Forest Fragmentation

Figure 5 shows the sensitivity analysis results of forest fragmentation intensity; obviously, the highest sensitivity was obtained at a CV = 60%, when the analysis unit was 3 × 3 km². Based on this analytical unit, a forest fragmentation severity map of each analytical unit within the study area was produced from the natural breakpoint criteria of low-fragmented areas (<6%), above-low-fragmented areas (6%–12%), medium-fragmented areas (12%–18%), sub-high-fragmented areas (19%–24%), and high-fragmented areas (>25%), as shown in Figure 6.

High fragmentation severity mainly existed in the northwestern portion of the study area, while low fragmentation severity was principally observed in the central, eastern, southwestern, and northeastern portions in the four years (Figure 6). Medium fragmentation severity always dominated the region, followed by above-low and sub-high severity. Additionally, more and more medium and sub-high fragmentation areas in Xinyu City were developing towards higher levels of fragmentation over time (Figure 6).

Geo-detection was employed to reveal the impact of human and environmental factors on forest fragmentation, and the detection results are summarized in

Table 7. Detection results of the natural factors responsible for forest fragmentation.

Factor	Altitude	Slope	Aspect	Annual Rainfall	Accumulated Sunshine Hours	Distance from Rivers
q-value	0.061	0.007	0.034	0.401	0.233	0.406

and

Table 8. Detection results of the human factors contributing to forest fragmentation.

Factor	Distance from Roads	Areal Ratio of Cropland	Areal Ratio of Construction Land	Distance from Consumption Centers	Land Use Diversity
q-value	0.634	0.589	0.423	0.403	0.618

, respectively. In Table 7, the q-values of altitude, slope, and aspect were all below 0.1, indicating that their influences on forest fragmentation were relatively weak, and the distance from rivers and annual rainfall were identified as the most significant natural factors triggering forest fragmentation, with their corresponding explanatory power values (q-values) of 0.406 and 0.401, respectively. The q-value of accumulated sunshine hours was 0.233, which indicated that it had a moderate influence among the natural environmental factors. In Table 8, the distance from roads, land use diversity, and the areal ratio of cropland were found to be the most influential anthropogenic variables on forest fragmentation, with corresponding q-values of 0.634, 0.618, and 0.589, respectively. The q-value of the areal ratio of construction land was 0.423, and ranked fourth among all the factors, indicating that the expansion of construction land had a strong influence on forest fragmentation. The q-value of distance from consumption centers was estimated to be 0.403, which was the lowest one among all the human factors, but still surpassed most natural factors (Table 7). The average q-values of human and natural factors were 0.533 and 0.19, respectively, indicating that the impact of human activities on forest fragmentation was notably greater than that of natural environmental factors in this particular region.

Forest fragmentation is the result of multifactorial interactions. The detection results of the interactive effects of different combinations of factors on forest fragmentation are shown in

Table 9. Detection result of the interactive effects of different combinations of factors on fragmentation.

X1	0.061
X2	0.086	0.007
X3	0.128	0.067	0.034
X4	0.244	0.187	0.208	0.401
X5	0.244	0.202	0.218	0.531	0.233
X6	0.084	0.033	0.131	0.467	0.327	0.406
X7	0.082	0.038	0.114	0.309	0.323	0.460	0.634
X8	0.214	0.193	0.256	0.470	0.369	0.612	0.840	0.589
X9	0.229	0.201	0.229	0.400	0.224	0.628	0.827	0.704	0.423
X10	0.089	0.065	0.098	0.283	0.246	0.473	0.696	0.468	0.538	0.403
X11	0.099	0.044	0.121	0.191	0.257	0.481	0.756	0.768	0.548	0.438	0.618
	X1	X2	X3	X4	X5	X6	X7	X8	X9	X10	X11

Note: X₁: altitude; X₂: slope; X₃: aspect; X₄: annual rainfall; X₅: accumulated sunshine hours X₆: distance from rivers; X₇: distance from roads; X₈: areal ratio of cropland; X₉: areal ratio of construction land; X₁₀: distance from consumption centers; and X₁₁: land use diversity.

. Regarding human factors, the interactive effects between the distance from roads and other factors showed the highest impact on fragmentation. The interactive effect of the distance from roads $\cap$ the areal ratio of cropland had the strongest explanatory power, with their q-value reaching 0.84, followed by the distance from roads $\cap$ the areal ratio of construction land (0.827), and the areal ratio of cropland $\cap$ land use diversity (0.768). In addition, the interactive effects of land use diversity $\cap$ other human factors, areal ratio of cropland $\cap$ other human factors, and areal ratio of construction land $\cap$ other human factors were also significant, with their q-values ranging from 0.4 to 0.8. For example, the q-value of land use diversity $\cap$ areal ratio of cropland was 0.768, and the q-value of land use diversity $\cap$ areal ratio of construction land was 0.548. Regarding natural factors, the interactive effects between the accumulated sunshine hours $\cap$ annual rainfall had the highest explanatory power of 0.531. Regarding the interactions between natural and human factors, the interactions between the distance from rivers $\cap$ the areal ratio of cropland and the areal ratio of construction land $\cap$ the areal ratio of cropland had the highest explanatory power values of 0.612 and 0.628, respectively. Meanwhile, the average annual rainfall, distance from rivers, and accumulated sunshine hours also showed high interaction power with other factors. The highest q-values were found for annual rainfall $\cap$ accumulated sunshine hours, areal ratio of cropland $\cap$ annual rainfall and distance from rivers $\cap$ annual rainfall, which were of 0.531, 0.470, and 0.467, respectively. The average q-values of the interactive effects between natural factors $\cap$ natural factors, natural factors $\cap$ human factors, and human factors $\cap$ human factors were 0.21, 0.26, and 0.66, respectively. These results indicated that the most dominant interactive effects affecting forest fragmentation in Xinyu City was the interplay between human activities.

4. Discussion

4.1. Assessments and Application of Land Cover and Fragmentation

The SVM model has been proved to be an effective and easily implemented classifier in many studies [56,73,74]. We not only directly applied the SVM model for classification but also refined the input features to improve the accuracy of the classification results, supporting the analysis of forest fragmentation (Table 3). Feature collection is an important step within the classification process; for example, an excessive feature set size has a negative impact on the classification accuracy of the SVM classifier [75,76]. After multiple experiments, we ultimately selected NDVI, NDBI, the brightness, greenness, and wetness derived from the K–T transformation, and an optimized combination of green, near-infrared, and shortwave infrared bands (bands 2, 4, and 7 for Landsat 5 TM and bands 3, 5, and 7 for Landsat 8 OLI) as the input features to obtain accurate land cover maps of Xinyu City. Based on the land cover classification results, it was evident that the loss of forest area is severe, reaching to 250 km² from 1999 to 2020. Cropland and water slightly decreased, while the construction land and unused land exhibited the most significant growth. Combining the results of this study with relevant research, it can be concluded that urban expansion leading to the loss of forest land and policies aimed at compensating for loss in cropland were the primary reasons for the reduction in forest area [47,77]. For instance, it was observed that as the amount of compensated cropland increased, the ecological cost also increased in the study conducted in Wuhan, which suggested that policies aimed at compensating for changes in cropland area, such as conversion of other land uses into cropland, can have ecological consequences and potentially contribute to forest fragmentation [78].

We preprocessed land cover images into forest/non-forest binary images, and then used the MSPA method to characterize forest fragmentation. Previous studies on forest fragmentation mainly focused on the quantitative expression of forest fragmentation, which has limited application compared with the results based on the MSPA model in this study [14,15,28]. For example, Dadashpoor et al. used landscape expansion indices, such as patch density (PD) and edge density (ED), to study landscape fragmentation, which can reveal the overall fragmentation situation, but cannot provide specific locations where different fragmentation situations occur [79]. However, MSPA has the capability to map spatial patterns at the pixel level, making it accurate in characterizing pattern features and sensitive to changes over time. Based on this, scholars have used MSPA to investigate changes in the connectivity of riparian habitats under different climate change scenarios [80]. Therefore, we utilized MSPA based on previous work to classify forest fragmentation components, enabling a more accurate revelation of the fluctuations in forest fragmentation. The results should be reliable and comparable relative to the existing similar studies.

Forest fragmentation is scale-dependent. Smaller landscapes are more sensitive to higher spatial frequency patterns, while larger landscapes are more sensitive to lower spatial frequency patterns. For example, fewer forest patches were defined as core patches with landscape-scale expansion [22,30,81]. We chose the eight-neighbor rule rather than the four-neighbor rule in the hypothesis that the eight-neighbor rule is closer to the real-world representation of neighbor influences [32]. The setting of edge width affects the recognition of core patches [81]. The previously small core patches would be recognized as islets as the edge width increases, while the previously narrower core area would be recognized as a bridge [82]. Kang recommended that the forest corridors should be set at minimum width of 35 m to meet the need of bird-supporting functions [83]. Considering the resolution of the image and the richness of species in the study area, we chose to set the edge width to 1 pixel (30 m) for forest fragmentation research, which is roughly in agreement with Kang’s argument. However, our study has certain limitations. For example, the MSPA approach categorized the landscape into forest and non-forest areas, potentially overlooking the influence of non-forest areas on forest patch connectivity [32]. In future studies, it would be worthwhile to simultaneously analyze the spatial distribution of non-forest areas to enhance our understanding of forest fragmentation dynamics. For instance, Jing et al. further classified non-forest regions into “holes” and “outer background” and proposed them as additional MSPA tools to understand the deforestation process in the Amazon [32], which provides more insights on forest landscape patterns and dynamics.

4.2. Factors Associated with Forest Fragmentation and Their Implications

Geo-detectors have been used in areas such as vegetation dynamics, land cover changes, and the effects of air pollution on health [40,70,84]. The utilization of a geo-detector benefits from its ability to consider the interactions between related factors [85], accommodate both numerical and qualitative data, and offer distinct advantages in assessing spatial variability. The results of the geo-detector showed that both natural and human factors accelerated the forest fragmentation in Xinyu City (Table 7 and Table 8). Meanwhile, the impact of human factors on forest fragmentation was greater than that of natural factors. Specifically, the distance from roads and land use diversity had the highest explanatory power q-values of 0.634 and 0.618, respectively, which indicated that the expansion of roads and town areas accelerated forest fragmentation. In terms of interaction detection, the interactions between human factors exhibited the highest impact on forest fragmentation, such as the q-value of distance to roads (X₇) $\cap$ the areal ratio of cropland (X₈) reached 0.84. This indicated that human activities were the dominant reason promoting forest fragmentation, which is consistent with related research results [86,87]. Meanwhile, in terms of the interactions between natural factors, the highest q-value of 0.531 was found between annual rainfall $\cap$ accumulated sunshine hours, which suggested that these two factors determined a good hydrothermal condition in this mid-latitude region to promote forest growth to achieve the purpose of a short rotational harvest of these plantations, thus accelerating forest fragmentation. The interactive effects of other natural factors were relatively small, especially for the topographic factor. On the one hand, this was attributed to the fact that forests in areas with high topographic variation remained relatively stable in the state; on the other hand, the study area was in the phase of accelerated economic and social development during the period 1999–2020, when human factors were principally responsible for forest fragmentation more clearly. Additionally, Xue et al. found that socioeconomic factors, such as road construction and rural expansion, have a greater impact on the fragmentation in oasis areas than natural factors [88]. Sulieman et al. found that the expansion of mechanized rain-fed agriculture, wood cutting, construction of infrastructure, and other human factors are the main causes of forest degradation and fragmentation in Sudan [89]. The study conducted in Lancang City found that the loss of habitat landscape connectivity was closely associated with urban expansion and distances to roads [47]. Our findings in terms of anthropogenic impacts on fragmentation in the current work are basically aligned with their specific observations. Given the high impact of human activities, it is necessary to plan human activities reasonably to mitigate their impact on forest fragmentation. For example, road construction in forests requires the clearance of vegetation, which can lead to deforestation. Research has found that 65% of forest cover change occurred within 2000 m of the nearest road [90]; proper planning of roads can help reduce their impact on forests. Meanwhile, the study on vegetation dynamics in Yunnan Plateau found that human factors had a weaker influence on vegetation dynamics, and soil type and elevation were the main influencing factors underlying vegetation dynamics [70], which differs from the specific observations regarding the driving forces in the current work.

4.3. Strategies for Mitigating Forest Fragmentation

Through the spatial plan of Xinyu City [91] and the related literature [92,93,94,95], we matched different fragmentation severity areas in 2020 to the local urban function positioning, and the corresponding results are shown in

Table 10. Functional positioning of different forest fragmentation severity zones in Xinyu City.

Fragmentation Level	Area Proportion (%)	Functional Positioning	Location
High	4.785	Ecological reserves	Northwest
Sub-high	20.574	Ecological reserves	Northwest and southern
Medium	33.493	Ecological restoration and safeguard zone	Western, central, and southern
Above-low	27.033	Industrial and agricultural development zone	Southwest and northeast
Low	14.115	Urban functional core zone	East-central

.

Specifically, the high and sub-high fragmentation areas are primarily situated in the northwest and south parts of the study area, which are mainly rural regions with dense forests, but relatively lack in forest protection measures and supervision. The spatial plan initiatives in these regions focus on soil and water conservation and ecological restoration, encompassing measures like reforestation, afforestation, and small watershed management. The medium fragmentation area is primarily located in the central, western, and southern parts of the study area, serving as a significant nature reserve and ecological barrier. Xinyu City has identified key bodies of water for protection, such as Xiannv Lake, and established core ecological barriers in areas like Mengshan Mountain in the northwest, Dagang Mountain and Jiulong Mountain in the south, and the Yangtiangang Mountain range in the central region. The above-low fragmentation area is mainly situated in the northeastern and southwest parts of Xinyu City, encompassing the primary distribution zones of permanent basic farmland in the city. Their main role is to support urban agricultural development. Additionally, these areas are dedicated to special agriculture and forestry industries, like oil tea, tangerine, and high-quality rice. The low fragmentation area is predominantly located in the east-central part of Xinyu City, with sporadic distribution in the central and western regions. This area serves as the core urban area of Xinyu City and the county town of Fenyi County, fulfilling essential functions for residents’ daily life, commercial development, and industrial production.

The focus of forestry development varies depending on societal needs and priorities, which differs from region to region. In order to carry out better forest management planning, we should consider the livelihood dependency on forests, forest income and benefits, and the capacity for forest management [96]. Based on the current forest fragmentation status, we have adjusted the forest area development strategy in Xinyu City to encompass three primary modes: ‘high-density development mode’, ‘intensive development mode’, and ‘conservation development mode’.

(1)

The high-density development model: This model is mainly applicable to medium fragmentation areas and focuses on the strict protection of forest resources, reduction in fragmentation, improvement of the overall quality of forests, and promotion of resource utilization. These specific measures include: implementing effective forest monitoring to prevent harm from invasive species, forest fires, and illegal encroachments; strengthening resource protection by categorizing and managing forests based on ecological redlines, urban development boundaries, and permanent basic farmland, among other limiting factors; developing ecotourism, sightseeing, and wellness-related ecological industries; and utilizing new technologies (such as remote sensing) for forest planning and management. This model aims to safeguard the original ecosystem and natural resources, achieve the sustainable development of forest resources, and strike a balance between ecological, economic, and social benefits.

(2)

The intensive development model: This model is primarily suited for regions characterized by high to sub-high levels of forest fragmentation. Its key focus is on safeguarding the integrity of existing forested areas, ensuring the connectivity of ecological corridors, and preserving the integrity of ecological patches. Its specific strategies include: promoting the return of cropland to forests on steep slopes in rural areas, afforestation, and the construction of protective forests; encouraging the conversion of agricultural land into forestry in hilly terrains and regions susceptible to soil erosion; and maintaining the existing forest ecosystems to ensure their continuity and health. This approach aims to preserve forest ecological environments, minimize conflicts between human activities and the environment, and strike a balance between environmental conservation and economic development. Meanwhile, the leakage effect should be taken into consideration when constructing certified or protected forests [97].

(3)

The conservation development model: This model primarily targets areas with low and above-low forest fragmentation, emphasizes the protection of forest resources, the expansion of urban green spaces, and increases public participation. Its specific measures include: increasing the availability of public green spaces and introducing aesthetically pleasing greenery to improve the urban landscape; protecting valuable old trees and implementing measures to manage forestry pests and diseases; and promoting public participation and raising awareness about nature conservation among the communities [98,99].

5. Conclusions

This paper aimed to characterize the spatio-temporal changes of forest fragmentation in Xinyu City, and to identify the underlying factors responsible for fragmentation, to offer recommendations for mitigating forest fragmentation. We employed Landsat images to generate land cover maps, which unveiled a substantial forest contraction of about 250 km² in Xinyu City from 1999 to 2020. The causes of forest fragmentation are complex. Through the MSPA and GDM models, the forest fragmentation conditions and causes were more accurately defined and quantified in the study area, allowing us to reveal the influence extents and degree of human and natural factors on forest fragmentation Following this, an appropriate forest resource regulation strategy was proposed based on the degree of fragmentation and local conditions. This approach not only broadens the perspective of forest fragmentation research but also enriches the practical theoretical framework for forest resource protection and development.