2.2.2. Moving Window Method

The moving window method in Fragstats 4.2 software was used to generate the landscape index raster map, and 1500 m was selected as the moving window radius after several calculations and comparisons [42].

#### 2.2.3. Determination of Evaluation Unit

Considering the scope of the Nanming River watershed and sampling workload, referring to existing studies [43], the watershed was divided into 2 km × 2 km evaluation units as ecological safety evaluation plots, with a total of 640 sampling areas. Based on this, the landscape ecological safety index of each plot was calculated separately (Figure 2).

**Figure 2.** Sampling areas for ecological security assessment of landscape pattern in the study area.

2.2.4. Construction of the Landscape Ecological Security Index

The ecological safety index (*ESIk*) of the land use landscape was calculated based on the landscape disturbance index and vulnerability index [44]. The equations are as follows:

$$ESI\_k = \sum\_{i=1}^{n} \frac{A\_{ki}}{A\_k} \times (1 - 10 \times LDI\_I \times LVI\_i) \tag{1}$$

$$LDI\_i = aC\_i + bF\_i + cF\_i \tag{2}$$

where *ESIk* in (1) is the landscape ecological safety index of the *k*-th evaluation unit, and a larger *ESIk* indicates a higher degree of ecological safety in the landscape and vice versa; *LDIi* is the index of landscape disturbance; and *LVIi* is the fragility index, based on the results of a previous study [45]. Each land use landscape fragility is specifically set to 5 levels: the value of constructed land is 1; a forest is 2; grassland is 3; cultivated land is 4, and water is 5. In this case, *n* is the number of landscape types; *k* is the number of evaluation units; *Aki* is the area of the class *i* landscape of the *k*-th evaluated cell; and *Ak* is the total area of the *k*-th evaluation cell. In (2), *Ci* is landscape fragmentation, and *Hi* indicates the diversity index. *Fi* denotes the number of watershed landscape sub-dimensions; *a*, *b*, and *c* are the weights of *Ci*, *Hi*, and *Fi*, respectively, which are assigned to 0.5, 0.3, and 0.2, following the existing research results and the conditions in the study area [46].

#### 2.2.5. Spatial Autocorrelation Analysis

We used spatial autocorrelation analysis to detect the spatial agglomeration of regional geographical phenomena [47]. Global Moran's I index was used to measure the overall spatial agglomeration characteristics of landscape ecological safety. The local spatial extent between regions was measured using the local Moran's I index [48,49]. The equations are as follows:

$$\text{GlobalNorm's I} = \frac{n \sum\_{i=1}^{n} \sum\_{j=1}^{m} \mathcal{W}\_{ij} \left(\boldsymbol{\chi}\_{i} - \overline{\boldsymbol{\varpi}}\right) \left(\boldsymbol{\chi}\_{j} - \overline{\boldsymbol{\varpi}}\right)}{\sum\_{i=1}^{n} \sum\_{j=1}^{m} \mathcal{W}\_{ij} \sum\_{i=1}^{n} \left(\boldsymbol{\chi}\_{i} - \overline{\boldsymbol{\varpi}}\right)^{2}} \tag{3}$$

$$\text{LocalNorm's I}\_{i} = \frac{(n-1)(\boldsymbol{x}\_{i}-\boldsymbol{\overline{x}})\sum\_{j=1}^{n}\mathcal{W}\_{ij}(\boldsymbol{x}\_{j}-\boldsymbol{\overline{x}})}{\sum\_{j=1}^{n}\mathcal{W}\_{ij}\left(\boldsymbol{x}\_{j}-\boldsymbol{\overline{x}}\right)^{2}}\tag{4}$$

where *n* is the number of grids; *x* is the average vulnerability in the study area; *xi* and *xj* are the attribute values of the *i*-th, *j*-th raster, respectively (*i* = *j*), where *i* = 1, 2, 3, ... , *n*; *j* = 1, 2, 3, . . . , m; *Wij* is the weight value; *Wij* = 1 when *i* and *j* are adjacent, and *Wij* = 0 when they are far apart.

### *2.3. Data Sources*

Landsat TM satellite 2000 images, 2010 images, and Landsat 8 satellite images of 2020 were used as the base data, which were obtained from the Geospatial Data Cloud platform (http://www.gscloud.cn/ accessed on 10 September 2022) with a spatial resolution of 30 m. First, the ENVI 5.3 software is used to pre-process the remote sensing images, such as geometric correction [50], atmospheric correction [51], and image enhancement [52], to complete the preparation and processing of fundamental geographic data. Second, a combination of supervised classification and human–computer interaction [53] was used to interpret and decipher the land use data into the waters, cultivated land, grassland, forest land, and constructed land according to the actual state of the watershed, which was not accounted for in the land use classification due to the tiny amount of unused land in the Nanming River watershed. Finally, the confusion matrix was utilized to rectify and validate the correctness of the results after land use classification [54], and the overall classification accuracy surpassed 85%, indicating that it could meet the analysis goals of this study. Furthermore, we obtained land use status maps for the study region in 2000, 2010, and 2020.

#### **3. Analysis and Results**

#### *3.1. Research Framework*

The changes in the characteristics of land use landscape were analyzed from 2000 to 2020 using GIS and RS technology. We also analyzed changes in the land use landscape pattern index from landscape patch scale and landscape scale. ArcGIS and landscape pattern index methods were used to construct the landscape ecological security index and to analyze the spatial and temporal variation mode of landscape pattern evolution and its ecological security (Figure 3).

**Figure 3.** Research framework.
