Multi-scale geographically weighted regression modeling of urban and rural construction land fragmentation—a case study of the Yangtze River Delta region

Landscape fragmentation, which can be natural or man-made, leads to decreasing species diversity and environmental degradation. The analysis of landscape fragmentation helps people understand the relationship between land and ecosystems in order to achieve sustainable development. Using the data from 306 counties in the Yangtze River Delta region, we calculated the landscape pattern index and constructed an evaluation index for urban and rural construction land fragmentation. We used a multi-scale geographically weighted regression model to analyze the determinants of land fragmentation and their spatial heterogeneity. We found that land fragmentation has significant spatial autocorrelation. The spatial distribution of the fragmentation is “ladder-like” and higher in the southwest and lower in the northeast. The determinants’ spatial scales from largest to smallest are average yearly precipitation, population (PS), temperature (TEM), gross domestic product (GDP) per capita (GPC), GDP per unit area (GPUA), proportion of secondary production (PSP), slope (SLO), and elevation (ELE). The mean contributions of the regression coefficient from highest to lowest were PER, ELE, TEM, location, GPC, SLO, GPUA, PSP, and PS. The evaluation index could help policymakers make urbanization plans and implement effective urban and rural spatial structures.


I. INTRODUCTION
Human activities related to economic and social development have caused great disturbances to urban and rural spaces; large amounts of arable and forest land are being transformed into urban and rural construction land, and urban sprawl is intensifying, leading to landscape fragmentation [1]. At present, the fragmentation and inefficiency of urban and rural land represent a significant challenge to comprehensive land management and the compact development of urban and rural construction spatial structures in the Yangtze River Delta. Urban and rural construction land is involved in economic and social activities. Thus, analyzing land fragmentation provides scientific support for comprehensive improvement and optimal land use layouts and for exploration of the economic efficiency of compact spatial structures in urban and rural areas. Fragmentation analysis is also essential for integrating regional development strategy. Therefore, studying the regional degree of fragmentation in urban and rural construction land, alongside related influencing factors in the process of urbanization, is of great significance to sustainable development and land management research.
A landscape's pattern mainly refers to its spatial arrangement, including the type, number, spatial distribution, and configuration of landscape units [2]. Landscape fragmentation is a multi-perspective concept, and related research can be classified in three main fields. The first constitutes studies into evaluation index systems for landscape pattern fragmentation. Landscape indicators have been widely used to quantify changes in landscape patterns; they can be divided into three levels: patch, category, and landscape [3]. Kamusoko and Aniya [4] categorized landscape fragmentation based on the study scale and six landscape indicators: maximum patch index, number of patches, patch density, average patch size, landscape shape index, and scatter and juxtaposition indices. They concluded that human activities are driven by agricultural expansion. To deepen their indicator system, Dewan et al. [5] added three more indicators: area-weighted mean fractal dimension, contagion, and Shannon's diversity index. Wang et al. [6] considered nine robust indicators for fragmentation metrics. In addition, Su et al. [7] identified three key criteria for selecting landscape metrics: 1) the metrics should represent different aspects of landscape characteristics, 2) the metrics should not be highly redundant, and 3) the metrics should be documented in relevant studies. Andrea et al. [8] designed and applied a landscape fragmentation composite metric to explore the influence of various factors on landscape fragmentation.
The second field comprises studies of spatial and temporal changes in landscape fragmentation patterns. These studies generally identify the distribution patterns and spatial effects of landscape fragmentation patterns, intending to find policy responses. Turner [9] conducted a landscape pattern analysis based on the spatial analysis program of raster cells. To identify policy solutions, Weng et al. [10] studied the spatial and temporal variations of landscape fragmentation patterns, by identifying the distribution patterns and spatial effects of landscape fragmentation patterns. Huang et al. [11] quantified the spatial and temporal characteristics of landscape patterns in the coastal bay area of southeastern China from 1988 to 2007. Canedoli et al. [12] used the effective mesh size method to assess landscape fragmentation in urban Milan. Zuzolo et al. [13] analyzed the extent of spatio-temporal fragmentation in the metropolitan city of Naples from 1990 to 2018; they used a landscape expansion index to determine the impact on the ecological mosaic.
The third field concerns the factors influencing landscape pattern fragmentation. These factors have been analyzed through global regression and geographically weighted regression (GWR) models. Many scholars have used GWR models to analyze the spatially variable relationships of landscape fragmentation with natural, economic, and social factors [14][15][16]. Some scholars have used a combination of ordinary least squares (OLS) and GWR models to analyze the changing relationships between land use and urban landscape patterns [1,17]. Miao et al. used a geo-detector approach to study changes in landscape fragmentation related to road buffers and road sections. That approach may satisfactorily explain the spatial and temporal characteristics of landscape fragmentation along the Qinghai-Tibet Highway [18].
The results of previous research have shown that GWR models outperform OLS models and can also detect spatial non-smoothness. However, the GWR models used in previous studies have two major limitations. First, when using GWR models, researchers have tended to ignore the multiple testing problem, resulting in questionable local parameters [19]. Second, a GWR model is a special type of multi-scale GWR (MGWR) model in which the parameters of all independent variables vary at the same spatial scale. Thus, GWR models can ignore the spatial scales of different influencing elements, leading to un-robust computational results. This is particularly important because the complex social, economic, and demographic factors associated with landscape fragmentation likely vary at different scales. The scaling nature of landscape spatial patterns dictates that there is a scale correlation between changes in landscape patterns and their influencing elements [20]. Landscape patterns have strong spatial properties, with spatial heterogeneity being the most common spatial property in the relationships between urbanization factors and landscape patterns [21,22]. From a spatial perspective, heterogeneous relationships between urbanization and landscape patterns vary at local, regional, and global scales [14,23,24]. Previous studies have highlighted the importance of multi-scale effects in the study of landscape patterns [25,26] alongside the relevance of the scaledependent features of landscape studies and multi-scale information regarding landscape patterns [27][28][29].
Fotheringham et al. [30] proposed the MGWR model based on a generalized additive model (GAM) [31]. This model uses the respective optimal bandwidth for each independent variable for regression, thus essentially solving the problem of there being different scales and bandwidths for different variables. While GWR models specify different optimal bandwidths for each variable (which usually reflects the average of the optimal bandwidths for all independent variables), the MGWR model builds on this by allowing the bandwidths of each independent variable to be different. By relaxing the assumption that the scales of the spatial processes are the same for all local coefficients, the MGWR model allows for a more powerful and explanatory model. Compared to the classical GWR model, the MWGR model has three main improvements. First, it allows for different levels of spatial smoothing for each variable and addresses the shortcomings of GWR models. Second, these covariate-specific bandwidths can be used as indicators of the spatial scale at which each spatial process acts. Third, the multi-bandwidth approach yields more realistic results. Thus, the model is more useful for spatial processes [32]. MGWR models are more intuitive and easier to interpret than spatially varying coefficient models. By treating the models as GAMs, Yu et al. [33,34] addressed the problem of statistical inference in MGWR models without local parameter estimates and obtained their hat matrices.
In addition, MGWR has some applications in geography, remote sensing, ecology, health, and other fields. Fotheringham et al. [30,35] applied MGWR to empirical datasets of the Irish Famine and Chinese air pollution, respectively. The results of the studies show that the MGWR model outperforms the GWR model and provides valuable information on the scale at which different processes operate. Oshan et al. [36,37] show that when the bandwidths estimated via MGWR strongly differ in magnitude from the bandwidth estimated by GWR, the parameter estimate surfaces and the associated tests of significance may also diverge. MGWR provides a richer but more parsimonious quantitative representation of the determinants of obesity rates. Additionally, Yang et al. [38] introduced MGWR into land surface temperature downscaling to establish a multiscale and nonstationary relationship between LST and biophysical indices and develop a hybrid method (MGWRK) coupled with area-to-point kriging. The results show that in LST downscaling, MGWR is capable of meticulously delineating nonstationary relationships between LST and multiple biophysical indices at different scales, which is expected to improve the temperature fidelity and spatial information enrichment in urban regions with diversiform landscapes and intricate land-atmospheric interactions.   [39] conducted a simulation study with empirical examples on obesity rates in Phoenix. The authors show that a quantitative measure of bandwidth uncertainty can be obtained from the confidence interval of the bandwidth utilizing red pool weights. Understanding bandwidth uncertainty provides important insights into the operational scales of different processes, especially when comparing bandwidths specific to covariates.   [40], using an MGWR model, investigates the relationship between landscape patterns and socio-economic factors in the urbanization process of Shenzhen from 2000 to 2015. The study showed that MGWR is a powerful extension of GWR, which not only reveals spatial heterogeneity but also measures the scale of operation of covariates. The above empirical results show that the MGWR model performs effectively than GWR.
In summary, the spatial and temporal distribution of landscape fragmentation, and the factors that influence it, are important issues in land studies. The MGWR model regresses each independent variable using its optimal bandwidth, thus basically solving the problem of different variables having different scales and bandwidths [30]. Therefore, this study applied seven indicators of the landscape pattern index of urban and rural construction land to 306 districts and counties in the Yangtze River Delta. The Integrated Fragmentation Index (IFI) of urban and rural construction land was calculated using the entropy weighting method, and spatial autocorrelation analysis was conducted using that index. Further, the factors affecting the IFI using the MGWR method were explored, revealing the spatial scales of different influencing factors.
To achieve these objectives, three lines of inquiry were explored: • How do landscape patterns change spatially in urban and rural areas? How can an evaluation system for landscape pattern fragmentation be constructed at a large regional scale?
• What are the spatial scales of action for the different factors influencing the MGWR analysis of the data point variables in the study area? • Why are there multi-scale effects of urban and rural construction land fragmentation? Through which mechanisms do these variables affect landscape fragmentation?

1) OVERVIEW OF THE STUDY AREA
The Yangtze River Delta is located at the lower reaches of the Yangtze River, bordering the Yellow Sea and the East China Sea. It lies between 27°0′ and 35°30′ N and between 114°54′ and 123°0′ E (Fig. 1), with a total area of 358,000 km 2 . The northern part of the Yangtze River Delta is low and flat, while the southern part is hilly and undulating. The elevation ranges from 0 to 255 m. It has a subtropical monsoonal climate with four distinct seasons. It includes 41 cities across the Shanghai, Jiangsu, Zhejiang, and Anhui provinces and was populated by 227 million people as of the end of 2019.

2) DATA SOURCES
The data sources for the study are divided into two main parts: (i) land use data from the "30 m global land cover data" from GlobeLand30 in 2020 (http://www.globallandcover.com/), with an accuracy of 30 m, and (ii) natural environment data and social and economic data. The natural environment and socio-economic data were obtained by reviewing existing literature. Following this process, the natural factors indicators selected were elevation (ELE), slope (SLO), average yearly precipitation (AYP), and temperature (TEM). The ELE data were obtained from the Geospatial Data Cloud (http://www.gscloud.cn). The SLO was derived from ELE for the Yangtze River Delta region. The AYP and TEM were obtained from WorldClim2.1 climate data for the period 1970-2020 (https://worldclim.org/data/index.html); the social data were obtained from the statistical yearbooks of 41 municipalities.

1) EVALUATION MODEL OF URBAN AND RURAL CONSTRUCTION LAND FRAGMENTATION a) SELECTION OF URBAN AND RURAL CONSTRUCTION LAND FRAGMENTATION EVALUATION INDICATORS
The landscape pattern fragmentation index was calculated using FRAGSTATS 4.2 software, which provides a brief description of these landscape metrics. There is no perfect metric for habitat fragmentation metrics, but many fragmentation metrics are useful under specific conditions and for specific biological problems [6]. After consideration of relevant references [2,[6][7][8]41,42] and the actual profile of the study area, seven landscape pattern metrics were combined, reflecting the area, density, shape, edge, and spatial relationships of land types. Principal component analysis was then used to reduce redundancy between these metrics and develop the IFI.

b) DATA STANDARDIZATION AND WEIGHTING OF INDICATORS
The data of the above seven indicators were calculated, and the data of each indicator were standardized. The weight of each indicator was calculated according to the definition of information entropy. The judgment matrix of each evaluation index was constructed to obtain the normalized judgment matrix. Information entropy can be defined as where i is the i th evaluation object, i.e., i area, and is the th evaluation index; The calculated weights of the indicators are shown in Table  II. By combining the indices of the seven evaluation indicators, the IFI of each district and county was derived and visualized (larger values indicate a higher degree of fragmentation).

2) SPATIAL AUTOCORRELATION ANALYSIS
The global Moran's I index can determine whether regions with similar IFI values are spatially clustered. Global spatial autocorrelation analysis of IFI was conducted to describe the overall distribution of the IFI in the Yangtze River Delta Edge density (ED), mean patch area (MPS), area-weighted mean shape index (AMSI), mean fractal index (MNFRAC), splitting index (SPLIT), landscape shape index (LSI), and aggregation index (AI).
The total length of the boundary per unit area, revealing the extent to which the spot is fragmented. Mean Patch Area 1 (MPS) = Reflects the degree of fragmentation; the smaller the value, the higher the degree of fragmentation.
The higher the value of this index, the higher the degree of fragmentation.
Mean Fractal Index (MNFRAC) Reflects shape complexity over the entire range of spatial scales.
Splitting Index (SPLIT) A measure of patch separation, with higher values indicating more fragmented patches.
Landscape Shape Index (LSI)

/√
A larger value for a measure of patch edges indicates a higher degree of fragmentation.
A smaller index of agglomeration indicates greater fragmentation.
region, and a global Moran's I index statement was obtained. It was calculated using the following formula: where N denotes the number of study subjects; denotes the spatial weight matrix; denotes the value of study object i (i.e., the IFI of a district cell); and � denotes the average of all objects (i.e., the average IFI of all districts and counties).
Local spatial correlation analysis (Anselin Local Moran's I) was performed on the IFI to generate the local indicator of spatial autocorrelation. It was calculated as follows: where N denotes the number of study objects; denotes the spatial weight matrix; is the value of the variable at position i (i.e., the IFI of a district cell); and ̅ denotes the average of all objects (i.e., the average IFI of all districts). The type of the spatial weight matrix is the adjacency matrix.

3) MGWR MODELS
GWR models can effectively account for spatial heterogeneity by allowing local parameters rather than global parameters [43,44]. This can be expressed as where ( , 1 ) are the spatial coordinates of the i th sample, is the number of the independent variables, represents the independent variables, ( , ) is the estimated coefficient of the th sample for the th variable, 0 ( , 1 ) is the intercept term, i is the error term, and is the dependent variable representing the landscape metrics.
The MGWR model was calculated as follows: where represents the bandwidth used for the regression coefficient of the th variable.
Each regression coefficient of the MGWR model, , is based on a local regression and has a specific bandwidth. This is the biggest difference between the MGWR and classical GWR models; in the latter, all variables of have the same bandwidth. The MGWR model's kernel functions and broadband selection continue to use the kernel function and broadband selection criteria of the classical GWR model [30,39]. Here, the most commonly used quadratic kernel function-corrected Akaike Information Criterion (AICc)was used. The variables in the model are normalized so that each variable has a mean of 0 and a standard deviation of 1. Thus, the bandwidth of the MGWR model is independent of the size and variability of each explanatory variable; this allows individual bandwidths to be explained and compared [30,45].
While the classical GWR model uses a weighted least squares estimation method, the estimation of the MGWR model is quite different; it can be seen as a GAM [31]: For GAMs, the back-fitting algorithm can be used to fit the individual smoothing terms. The back-fitting algorithm first requires the initialization of all smoothing terms. This means that the coefficients in the MGWR model need to be estimated upfront [30,33,34,46]. There are generally three options for initialization: (i) classical GWR estimation, (ii) semiparametric GWR estimation, and (iii) least squares estimation. Here, the classical GWR estimate was chosen as the initial estimate. After the initialization settings were determined, the difference between the true value and the predicted value obtained from the initialization estimate (i.e., the initialization residual) , was calculated as follows: This residual, ̂, plus the first additive term, ̂1 , were regressed against the first independent variable, 1 , using a classical GWR regression to find the optimal bandwidth, 1 . A new column of parameter estimates for ̂1 and ̂ was then calculated to replace the previous estimates. Then, the residuals plus a second additive term, ̂2 , were regressed on the second variable, 2 . Then, the parameter estimates for ̂2 and ̂ were updated for the second variable. This process was then repeated until the last independent variable (the k th independent variable, ) was reached. The above integral constitutes a step that can be repeated until the estimates converge upon a convergence criterion [46]. Here, the classical residual sum-of-squares variation ratio (RSS) was used as the convergence criterion: where represents the sum of the squared residuals from the previous step and represents the sum of squared residuals from this step. Here, MGWR2.2 software (https://sgsup.asu.edu/sparc/mgwr) was used to calibrate all of the MGWR models. Fig. 2 shows that the distribution of urban and rural construction land in the Yangtze River Delta region exhibited obvious clustering characteristics. Additionally, there was an obvious correlation between the amount of urban construction land and the city level. In the eastern part of the Yangtze River Delta region, urban and rural construction land along the Shanghai-Suzhou-Nanjing economic corridor showed obvious clustering and contiguous effects. Furthermore, there were strong contiguous characteristics along the Shanghai-Hangzhou economic corridor. Due to their location in the plain area, urban and rural construction areas in Anhui and northern Jiangsu were scattered in a star-like pattern, with high patch density; Zhejiang and southwestern Anhui exhibited a lower density. Due to southeastern Zhejiang's coastal location, urban and rural construction land was distributed in a strip-like pattern.

A. IFI DISTRIBUTION CHARACTERISTICS
The IFI of each district and county was calculated, and the degree of fragmentation of urban and rural construction land in the Yangtze River Delta was divided into 10 classes; the division results are shown in Fig. 3. The degree of fragmentation of urban and rural construction land in the Yangtze River Delta was distributed in a "ladder-like" pattern, gradually increasing from the southeast coast to the southwest.

B. ANALYSIS OF SPATIAL AUTOCORRELATION RESULTS
Spatial autocorrelation analysis of the IFI in the Yangtze River Delta region yielded the results of Moran ' s I analysis statements, and the Anselin Local Moran's I analysis map shown in Fig. 4. The global Moran's I index of 0.37, Z-value of 10.85, and P-value of 0 were all significant at the 5% level, indicating that there was a significant spatial clustering effect of the urban and rural construction land with the IFI values in the districts and counties (i.e., there was a spatial positive autocorrelation). As can be seen from Fig. 4, 21 districts and counties had high-high aggregation distributions that passed the test; these areas were distributed in the mountainous areas of southern Anhui and southern Zhejiang. The industries in these areas are mainly based on agriculture and tourism, and their complex topographic conditions make the distribution of urban and rural construction land relatively fragmented. All of these areas exhibited high IFI values. Sixty-six districts and counties exhibited low-low aggregation distributions; they were distributed in large and medium-sized cities in the Yangtze River Delta plain. These areas exhibited developed economies, strong social vitalities, high urbanization levels, and small IFI values. There were eight districts and counties with low-high distributions; they were mainly located in articulation zones between low mountain ranges and low hills (i.e., topographies mainly consisting of hills and plains). These areas were less influenced by topography and had smaller IFI values compared to their neighboring districts and counties.

1) MODEL COMPARISON
The OLS, GWR, and MGWR analyses were performed for IFI and its influences. As seen in Table III, the goodness-of-fit R 2 of the MGWR model was higher than those of the classical OLS and GWR models; its AICc values were smaller than those of the classical OLS and GWR models. This demonstrates that the MGWR model achieved a better fit than the GWR model. The MGWR model also had fewer effective parameters and a lower sum of squared residuals. This outcome indicates that the regression results that were obtained with fewer parameters were closer to the actual values. The classical GWR model ignores the diversity of each variable scale, which introduces large amounts of noise and bias in the regression coefficients. This leads to instability regarding the regression coefficients. The MGWR model also proved to be less prone to multicollinearity problems than classical OLS and GWR models. Therefore, based on the above results, the MGWR model was considered superior, because it explores the factors influencing landscape fragmentation.

2) SCALE ANALYSIS
The MGWR model can reflect the magnitude of the spatial scales of the action of different variables, and it can model the spatial variations of different processes at different spatial scales. As all variables are normalized to have a mean equal to 0 and a standard deviation equal to 1, the bandwidth calculated by the MGWR model can be used as a direct indicator of the spatial scale of variation in the individual relationship between the IFI and each covariate [30]. The classical GWR model's bandwidth is 94; the MGWR model's bandwidth calculated here is shown in Table IV. The ELE action scale was 43, accounting for 14.05% of the total sample. The SLO action scale was 49, accounting for 16.01% of the total sample; these were much lower than the action scales of other variables. This indicates that ELE and SLO affected IFI on a relatively local scale and that the values of these two variables were the same within the action scale. If this scale were exceeded, then the coefficients would change dramatically. This result also indicates that the IFI was very sensitive to the natural topography. Furthermore, the scale effect of gross domestic product (GDP) per capita (GPC) was also very small (87), which indicates that the variation of urban and rural construction land fragmentation varied more spatially with GPC. The TEM action scale was 285, which was close to the global scale, and the coefficients were stable in space. The intercept term, AYP, and PS all had a global scale of 306, indicating that there was no spatial heterogeneity for each of these variables.

3) COEFFICIENT SPACE PATTERN ANALYSIS
The local parameter estimates generated by MGWR reflect the possible spatial heterogeneity in the process of influencing IFI. To explore the spatial variations of the effects of each indicator, the MGWR model results were analyzed for each variable indicator; statistical descriptions of the MGWR coefficients are shown in Table V. For the regression coefficients of the intercept, ELE, SLO, AYP, TEM, proportion of secondary production (PSP), GPC, and GDP per unit area (GPUA) were more significant in the MGWR results, while the regression coefficients of PS were not significant.
The effect of location on IFI as reflected by the intercept was negative, and a clear circular structure was visible. As shown in Fig. 5a, the high values were mainly concentrated in the eastern plains of the Yangtze River Delta. The intercept values ranged between -0.40 and 0, with a mean value of -0.394 and a standard deviation of 0.006. As the intercept is a global variable, the overall impact of different zones on the IFI in each district and county did not vary greatly. In terms of the absolute value of the coefficient, the intensity of its effect was greater among all variables. In the MGWR model, under the control of a given set of spatial factors, the intercepts in the local models take into account the remaining spatial variation. After normalizing the variables, independently of each potential spatial variation effect of the other explanatory variables, it is possible to make local inferences regarding the remaining spatial variation. The intercept represents the   effects of different locations on the IFI when the other independent variables have been determined. This study was conducted with controlled natural, social, and demographic factors, so that the intercept reflects the effects of other locational factors such as the built environment and urbanization rate on the IFI. The parameter estimates of the local intercepts were meaningful because they indicated a lower (significantly negative) IFI. The lower IFI level in the eastern Yangtze River Delta may have been caused by the larger urban size and more compact urban built-up land in the plains, which led to a lower IFI.
The ELE had a significant negative effect, as shown in Fig.  5b. The ELE coefficient ranged between -1.33 and 0.59, with a mean value of -0.616 and a standard deviation of 0.409. In general, its coefficient value varied widely from region to region within the Yangtze River Delta. In terms of the absolute value of the coefficient, the intensity of its influence was greater and an important factor affecting the IFI. Lower values were mainly concentrated in the northern plains and metropolitan areas, such as Shanghai, while higher values were concentrated in the southwestern Yangtze River Plain. These differences can be explained in two ways. First, in the northern plains, rural settlements were situated in relatively flat areas with more rural PS. Thus, the number of rural settlements was higher, leading to a higher degree of landscape fragmentation. Furthermore, metropolitan areas (i.e., Shanghai) were located in the lower ELE areas, leading to a higher degree of landscape fragmentation due to the higher degree of urbanization and the larger scale of construction land. Second, due to the need for flood control in the southwestern Yangtze River Plain region, settlements were situated in areas with relatively high ELE. This placement led to a higher degree of landscape fragmentation in this region. In other areas, the relationship between ELE and IFI was not significant, probably because the elevation values were smaller or larger. Furthermore, there were fewer rural settlements in the region, and there were few large elevation changes.
SLO positively affected IFI, as shown in Fig. 5c. The SLO coefficient values range between -0.588 and 0.964, with a mean of 0.145 and a standard deviation of 0.187. Regarding the absolute value of the coefficient, the intensity of its effect was smaller among all variables. In general, its coefficient values were concentrated in smaller regions, and the magnitude of variation varied widely across the different regions. Previous studies have suggested that SLO positively affects landscape fragmentation in rural settlements [15,47]. The results of this study supported this conclusion, but only in the smaller areas in the southwestern part of the study area. IFI was influenced by SLO; the larger the SLO, the more dispersed the distribution of settlements and the greater the IFI. Within the study area, most of the built environments were located in areas with low SLO, whereas most of the green spaces were located in high SLO areas. In general, the higher the human activities, the higher the landscape fragmentation. Therefore, the landscape fragmentation pattern usually coincided with the SLO pattern. This occurred because rural areas are dominated by agricultural production activities, arable land is generally distributed in plains and mountain valleys, and (for the convenience of production) settlements are usually located near arable land. The larger the SLO, the more difficult the area is to cultivate. At the same time, in areas of larger SLO, the higher cost of engineering construction and the problems of soil erosion and landslides are not conducive to the formation and development of rural settlements, which tend to be smaller and irregularly shaped. AYP positively affected IFI; its regression coefficient distribution showed an obvious "step-like" structure from north to south, as shown in Fig. 5d. The value of the AYP coefficients ranged from 0 to 0.82, with a mean value of 0.755 and a standard deviation of 0.033. Regarding the absolute value of the coefficient, its influence was the strongest of all of the variables, and it was the most important factor affecting IFI. Overall, there was little difference in the impact of AYP across the different regions. The high values were mainly concentrated in the mountainous areas of southern Zhejiang, where rainfall was abundant and rural settlements were more scattered. Previous studies have shown that people prefer to live by the water and will choose to settle in places with abundant rainfall and humid climates. This explains the strong positive correlation between AYP and landscape fragmentation in the Yangtze River Delta. The weakening effect of precipitation on the northern plains occurred due to the decrease in AYP.
The TEM had a negative effect throughout the study area, and its regression coefficient distribution showed a clear "stepped" distribution from northeast to southwest, as shown in Fig. 5e. The coefficients ranged from -0.62 to 0, with a mean value of -0.456 and a standard deviation of 0.053. The coefficient values varied little from region to region, with high values concentrated in the southwest mountains and low values concentrated in the north China plain. This indicates that TEM affected IFI on a global scale; in other words, TEM had a similar effect on IFI in all regions of all study areas. In terms of the absolute values of the coefficients, the intensity of the effect was greater among all variables. TEM had a significant effect on the location of settlements; people preferred to establish settlements in areas with warm climates (with warm winters and cool summers). In the study area, the warm and humid climate attracted people to settle. At the same time, TEM was different in urban versus mountainous areas due to the existence of the "urban heat island effect" in the former.
PSP negatively affected IFI, as shown in Fig. 5f. The coefficients were in the range of -0.389 to -0.003, with a mean value of -0.088 and a standard deviation of 0.115. The significant impact area was small and was concentrated in the southwestern mountainous hilly area. This may have been because the secondary industry in the area mainly comprised agricultural and sideline products and because the processing plants were closer to residential areas. When the value-added ratio of the secondary industry was low, the industries were mostly labor-intensive, relatively capital-poor, and weak in expanding their reproduction. Thus, it contributed little to the expansion of construction land and even inhibited the formation of new construction land after losses and hardships. Thus, this effect was instead manifested in the reorganization of the landscape [48]. Urban planning and management policies also had an obvious impact on the layout of industrial areas, which were relatively compact; this led to a lower degree of landscape fragmentation.
GPC proved to be an important indicator of the intensity of human economic and social activities at the district and county scale. It negatively affected the IFI, as shown in Fig. 5g. The coefficient ranged between -0.92 and 0, with a mean value of -0.196 and a standard deviation of 0.265. This influencing factor varied widely from region to region. In terms of the absolute value of the coefficient, the intensity of its influence was small. The GPC bandwidth was 87 and its significant influence was concentrated in the southwestern mountainous region. Both of these pieces of evidence suggest that GPC influenced the IFI at a very local scale. A previous study compared 120 cities globally and concluded that cities with higher incomes were more dispersed [49]. Here, the contrary finding was obtained, but only in the western part of the study area. This may have been due to the different perspectives of the sampling scale. There was no significant association between GPC and the IFI in the eastern part of the study area. These differences can be understood in two ways. First, GPC was much lower in the western part of the study area than in the other regions. According to the marginal utility principle, the southwestern region would have been more sensitive to the impact of human economic activities. In the western part of the study area, there was relatively little economic development; rural settlements and small cities were dominant, and the structural pattern of urban land was relatively loose. As people's economic income increases, the need for urban expansion intensifies, thus reducing landscape fragmentation. Second, the central and eastern regions of the study area were more economically developed. There was no significant relationship between GPC and IFI in these economically developed areas, probably because of an increase in the income levels of residents compacted urban and rural settlements under the higher urbanization rate. GPC played more of a role in the reorganization of construction land, which promoted the clustering of patches and decreased fragmentation [50]. People in developed areas tend to be more aware of planning, and the government is more capable of controlling urban development and urban forms in these areas.
GPUA is an indicator of the intensity of human economic activity per unit area of land. GPUA positively affected the IFI, as shown in Fig. 5h. Its coefficients ranged from 0.013 to -0.271, with a mean value of 0.088 and a standard deviation of 0.074. The intensity of its impact was small in terms of the absolute values of the coefficients. The significant impact range was small and concentrated in the southwestern mountainous hilly area; the differences in the impact of GPUA were not significant between the different regions. Previous studies have shown that human activities have a complex role in landscape fragmentation [51]. These findings suggest that in underdeveloped areas, human economic activities are one of the causes of landscape fragmentation.
The socio-economic effects on landscape patterns were more complex than the geographical factors. An increase in PS decreased landscape fragmentation for built-up land and enhanced landscape clustering. This is similar to the findings of Che et al. (2020) [50], but here the regression coefficients were not significant. Considering the complexity of socioeconomic activities, there may have been complex interactions between the factors selected in this study. In addition, the strengths and mechanism of their effects still need to be further explored.

IV. DISCUSSION
The main objectives of this study were (1) to construct an evaluation system for the landscape fragmentation of urban and rural construction land; (2) to examine the intensity, nonstationary, and spatial dependence between the patterns of landscape fragmentation and natural, social, and economic factors in urban and rural construction land; (3) to explore the scales of spatial effects of influencing factors; and (4) to analyze the existence of multi-scale effects and influencing mechanisms on landscape fragmentation. Based on land use and other data, this study applied the entropy weight method and the MGWR model, alongside statistical inferences, to empirically study the factors that influence landscape fragmentation and analyze the spatial scales at which different influencing factors act. Natural and socio-economic factors were found to have profound and significant effects on landscape fragmentation. This study emphasizes the importance of considering the scales of action of influencing factors when studying landscape patterns. The results of the study can improve urban managers' understanding of landscape pattern fragmentation and can act as a reference for the management and formulation of policies according to local conditions.
The MGWR model can help to reveal the factors that influence landscape pattern fragmentation. Multicollinearity and endogeneity can be problematic for GWR models due to their complexity and multifactoriality [37]. MGWR can overcome these limitations and produce a better fit. some bandwidths in MGWR are approximately 306, a value close to the best possible maximum distance between any two samples. This finding suggests that the relationship between IFI and Intercept, AYP, and PS is likely to be global. Other bandwidths were restricted to between 43 and 285, reflecting low, medium, and high levels of heterogeneity in the relationships, suggesting a high level of spatial heterogeneity in the relationships between IFI and ELE, SLO, and GPC. Compared to previous studies [14][15][16][17][18], MGWR can detect the spatial scale of action of different influences well through the variation in bandwidth.
Compared to existing studies, the innovations of this paper are: 1. This paper introduces MGWR into the study of the impact mechanism of landscape fragmentation, further enriching the application cases of MGWR models. By comparing OLS regression, GWR, and MGWG models to examine the spatial relationships between landscape fragmentation and related influencing factors, the results demonstrate that MGWR produces a more realistic and useful model of spatial processes. The results validate both the findings of previous studies and the credibility of the MGWR results, and this technique is worthy of further use in future research. Compared to previous studies, this paper further enriches the research framework by selecting multi-dimensional indicators such as natural, social, and economic indicators to construct the research framework system. 2. The construction of an evaluation index system for landscape pattern fragmentation is of great significance for local land management. The scientific construction of an evaluation system to identify fragmented areas is an important prerequisite for land remediation and the formulation of land management policies. In this paper, by introducing information entropy, the weights of each indicator are calculated so that IFI can be calculated, making the evaluation results more objective. MGWR helps to identify the global, regional, and local determinants of potential factors affecting landscape fragmentation and improves the ability of the model to explain the local landscape fragmentation situation, thus facilitating more specific policy formulation to mitigate the negative impacts of landscape fragmentation. The impact on urban and rural planning is understood in terms of changes in locally estimated coefficients in the spatial dimension. In contrast to previous studies, MGWR can quantify and reflect the level of urbanization and ecological conservation, taking into account both operational scale and spatial variations. The results obtained from MGWR not only identify the various impacts of urbanization but also reflect the relationship between global and local social and economic factors and urbanization and the different landscape indicators. MGWR can generate a specific bandwidth for each driver, from which the effects of urbanization variables can be assessed as they operate over space. The results of the bandwidths can help managers to develop master planning policies or provide targeted and practical local guidance. Based on the spatial coefficients of variation obtained from the MGWR runs, it is possible to determine the different effects of different urbanization factors on landscape patterns. The local impact of urbanization factors on landscape patterns reflects the clear differences between the areas east and west of the Yangtze River Delta region, indicating that the overall spatial structure of towns in the Yangtze River Delta region develops from west to east. The intercepts calculated by the MGWR model have important implications for subsequent research. Unlike OLS models, the MGWR model's intercepts are statistically nonzero, and, after controlling for the variables in the model, the spatial heterogeneity in the parameter estimates identifies hotspots for high or low IFI values. These spatial patterns may include the effects of both geography and geographic patterns associated with omitted variables. For example, factors-built environment, distance from large cities, and car ownership-may play important roles in urban sprawl [49]. Intercept distances may help to identify other spatial determinants. Previous evidence indicates that factors such as distance from large cities and urbanization rates may influence the degree of landscape fragmentation [15,51]. Thus, it will be necessary to consider these determinants in future multi-scale analyses.
However, there are several limitations of this study. The first is the limitation of using conventional metrics to study landscape patterns. Many class-level pattern metrics are highly correlated with habitat richness, which makes their use as indicators of habitat fragmentation problematic [6]. Furthermore, describing landscape patterns only through derived metrics misses methods that can describe the overall degree and spatial distribution of fragmentation on any taxonomic land cover map. Second, due to data availability and modeling method limitations, this study only analyzed the study object from cross-sectional data. However, changes in landscape patterns are dynamic processes, and there is a lack of continuous data and methodological models for validation. Additionally, the influences of variables regarding regional humanities and planning policies on landscape patterns were not considered. Finally, the results of this study were influenced by several aspects such as the spatial resolution and sampling methods. District and county scale sample sizes were utilized here, so there was no fine-scale calculation of landscape fragmentation; larger-scale analysis may mask errors in model calculation, leading to inaccurate results. The sample size was also small, and the size of the sample capacity may have affected the accuracy of the conclusions. At the same time, this study was a case study, and its findings were limited to the Yangtze River Delta region. Case validations in other regions are needed to generalize these findings.
Therefore, future research should mainly focus on the following aspects: first, innovations in landscape pattern evaluation indices and evaluation systems. The characteristics of the habitat determine the species ecological niche [31], and the specificity of the ecological niche of the species under study should be considered when selecting the landscape pattern metric. Attempts should also be made to apply a holistic approach to describing landscape fragmentation using the general concepts of contagion, complexity, and spatial entropy. Spatial syntax, complex network modeling, and circuit theory could also be applied to the exploration of landscape pattern distribution characteristics. Looking ahead, multi-scale investigations across patch mosaics, gradients, and graphical network models should move beyond simply describing metric variations with scales to using scale relationships to identify pattern-process relationships [52].
The second area of focus should be exploring and applying a multiscale geographically and temporally weighted regression (MGTWR) model. The MGWR model should be extended to the panel scale to incorporate a spatio-temporal aspect. The MGTWR model could provide a new reference perspective regarding spatial heterogeneity for issues, such as spatiotemporal evolution mechanisms and the factors influencing landscape pattern changes. The third focus area should be spatial modeling analyses of landscape pattern changes based on spatial interaction models (such as Poisson's gravity, negative binomial gravity, spatial filter gravity, spatial metric interaction, Hurdle gravity models, etc.). These models could explain the spatial and temporal evolution mechanisms of regional landscape patterns or could predict future interregional landscape pattern evolution paths and patterns. Fourth, landscape pattern research should be conducted from the perspective of multiple scales. To accurately describe the landscape pattern, the observed scale must accurately correlate with the intrinsic scale. Therefore, in the future, the scaledependent relationship of the modified effective grid size should be further investigated to improve the pattern analysis and regression models in landscape ecology. Finally, there were some biases and unexplained phenomena in the results due to the limitations of socio-economic factors such as urbanization and landscape metrics. Therefore, future studies should aim to collect more data, including the socio-economic factors of urbanization, as well as more comprehensive landuse time-series data and more detailed classification information. This will allow for a more comprehensive analysis of landscape patterns and urbanization processes in the Yangtze River Delta.

V. CONCLUSIONS
In this study, the MGWR model was used to more precisely analyze the spatial heterogeneity of the degree of influence of each indicator on the fragmentation of urban and rural construction land in the Yangtze River Delta.
After empirical analysis, the following conclusions were obtained: 3. The fragmentation of urban and rural construction land in the Yangtze River Delta was spatially autocorrelated, with a global Moran's I index of 0.37. The degree of fragmentation was higher in the southwest and lower in the northeast, and the spatial distribution was "ladder-like." From a spatial perspective, high-high clusters were found in the mountainous areas in the southwest, while low-low clusters were concentrated near the metropolitan areas of large cities. The IFI was very sensitive to ELE and SLO factors at the regional and county scales, and its influence scale was the smallest. The intercept term, AYP, and PS were all global scale variables, and their effects are relatively smooth in space. The effect of PS on the IFI was not significant. The other significant variables all had some spatial heterogeneity, and from largest to smallest, their spatial scales were TEM, GPUA, PSP, and then GPC. 4. The IFI is sensitive to ELE and SLO factors at the regional and country-wide scales of the Yangtze River Delta, and its scale of influence is the smallest. The intercept, AYP, and PS are global scale variables, and their effects are relatively smooth in space. The effect of PS on IFI is not significant. The other significant variables all have some spatial heterogeneity, and their spatial scales are TEM, GPUA, PSP, and GPC, from largest to smallest. 5. SLO, GPUA, and AYP all positively influenced IFI, while ELE, PS, GPC, and TEM negatively influenced it. AYP was the most dominant influencing factor, followed by ELE. TEM, the intercept, and GPC had stronger influences on the IFI, while SLO, GPUA, PSP, and PS had weaker influences. The calculation results of the MGWR model were better than those of the previous GWR and OLS models and were also applicable to the influencing factors of urban and rural construction land. The results of this study can provide guidance for urban planning and the development of urbanization policies. In the future, multi-dimensional extended analysis of the dynamic causes of landscape fragmentation should be conducted, considering land-use intensity analysis and the composition and structure of landscape components. More ecological processes should also be added to the MGWR model to provide a more comprehensive research reference for sustainable development and biodiversity conservation, among other uses.