## 2 Materials and Methods
### Methods
#### 2.3.3 Mean-Variable-Point Analysis Method  
In accordance with the theory of geomorphic development, an optimal analysis area for a specific geomorphic type exists, where the maximum elevation difference remains relatively stable. Furthermore, this area can, to a certain extent, maximize the expression of geomorphic integrity [28]. Hence, it is crucial to ascertain the optimal analysis window for terrain characterization parameters. Consequently, the mean-variable-point method [29] was utilized to determine the size of the optimal statistical window for extracting terrain feature parameters in the study area, using the degree of terrain relief as an example.  
First, using the Focused Statistics tool of the ArcGIS Neighborhood Analysis Toolbox with the rectangle analysis window type, the statistical cell sizes (2 \(\times\) 2, 3 \(\times\) 3,..., 48 \(\times\) 48) were sequentially calculated to correspond to the relief degree of the land surface. Then, the logarithmic fitting of the terrain undulation degree and the area corresponding to the window was carried out using Excel2021 software, and the fitting equation was obtained: \(y\) = 21.273 ln(x) \(-\) 160.49, as well as the coefficient of determination: \(R^{2}\) = 0.9501, which was a good fit.  
The calculation process of the mean-variable-point analysis method is as follows:  
1. The average terrain relief of \(N\) (2 \(\times\) 2, 3 \(\times\) 3,..., 48 \(\times\) 48) under incremental windows is obtained based on the terrain relief formula, and then the unit terrain relief \(T_{i}\) is calculated sequentially using the formula: \[T_{i}=\frac{t_{i}}{A_{i}}\] (1) where \(i\) is the rectangular analysis window; \(T_{i}\) is the topographic relief under the analysis window; \(t_{i}\) is the average relief under the analysis window; and \(A_{i}\) is the area under the analysis window (m\({}^{2}\)).
2. \(T_{i}\) is taken logarithmically (ln \(T\)) to obtain the sample series \(X\{x_{t},t=2,3,\ldots,48\}\).  
Figure 3: Fitting curve of river network density and catchment area threshold.  
3. The arithmetic mean \(\overline{X}\) and the sum of squared deviations \(S\) for \(X\) is computed: \[\overline{X}=\sum_{t=j}^{N}\frac{x_{i}}{N},S=\sum_{t=j}^{N}(x_{i}-\overline{x})^ {2}\] (2) where \(\overline{X}\) is the arithmetic mean of all samples; \(N\) is the total number of samples; \(S\) is the total sum of squared deviations.
4. Let \(i\) = 2, 3, \(\ldots\), \(N\), and for each \(i\) divide the sample into 2 segments: \(x_{1},x_{2},\ldots,x_{i-1}\), and \(\frac{x_{i}}{X_{i_{1}}},\frac{x_{i+1}}{X_{i_{2}}}\), and \(S_{i}\): \[S_{i}=\sum_{t=1}^{i-1}\bigl{(}x_{t}-\overline{X_{i_{1}}}\bigr{)}^{2}+\sum_{t=i }^{N}\bigl{(}x_{t}-\overline{X_{i_{2}}}\bigr{)}^{2}\] (3) where \(S_{i}\) is the difference between the sum of squared deviations of the two sample segments.
5. Calculate the expected value: \[E(S-S_{i})=E\Bigl{(}N^{-1}(i-1)(N-i+1)\bigl{(}\overline{X_{i_{1}}}-\overline{ X_{i_{2}}}\bigr{)}^{2}\Bigr{)}.\] (4)  
A scatterplot was made with serial number \(i\) as the horizontal coordinate and \(S-S_{i}\) as the vertical coordinate (Figure 4), where the vertex, which is the change point, is the point with the largest \(S-S_{i}\) value.  
As shown in Figure 4, the variable point is serial number 15. Therefore, the point corresponds to an optimal analysis window of 16 \(\times\) 16 and an optimal statistical area of 0.04 km\({}^{2}\).