• Geographical Detector Model

The geographical detector model (Geodetector), which is based on the theory of spatial stratified heterogeneity [67], was used to analyze the interaction between landscape metrics and various influencing factors in this study. First, we obtained the spatial data of independent variables through discrete classification for various influencing factors using the geometrical interval method [68–70], and then analyzed the influence of these variables on each landscape metric of the same spatial scales via Geodetector.

Specifically, the factor force q in Equation (5), ranging from 0 to 1, quantified the effect of different influencing factors on the spatial distribution of landscape metrics [32], and reflected the degree of spatial stratified heterogeneity of the metrics. The larger the q became, the more heterogeneous the landscape pattern.

$$\mathbf{q} = 1 - \frac{1}{\mathcal{N}\sigma^2} \sum\_{i=1}^{L} \mathbf{N}\_i \sigma\_i^2 \tag{5}$$

where N and σ 2 stand for the number of units and the variance of the dependent variable, respectively; i = 1 · · · L is the stratification of the dependent or independent variable; N<sup>i</sup> and σ 2 i stand for the number of units and the variance of the dependent variable in stratification layer i, respectively.

#### **3. Results**

#### *3.1. Spatiotemporal Variations of Land Use Types*

As shown in Figure 2, the overall landscape of the LXH is dominated by forests and croplands, while other seven land use types (including forest, shrub, orchard, grassland, water, floodplain and unused land) occupy a relatively small portion, accounting for less than a quarter of the total basin. Clearly, it is noticed that the proportion change in land use types were not large in general, but their temporal and spatial differences were obvious. Therein, the proportion changes of cropland and construction land were more prominent than any other land use types, particularly in the middle and lower watershed during 2000 and 2010. Temporally, the decreases of cropland and increases of construction land in this decade were more than 50% of that in total 35 years (See subfigures b in Figure 2). Spatially, it can be seen from the results of (c) in Figure 2 that the cropland decreased by 16.96% in the lower watershed and 2.70% in the middle watershed, respectively. The construction land increased by 17.78% in the lower watershed and 3.33% in the middle watershed, respectively. Their changes were all less than 1% in the upper watershed.

**Figure 2.** The types of land use in the LXH (1980–2015) (**a**) Spatial distribution of land use; (**b**) Percent coverage of the land use types; and (**c**) percent coverage of the land use types in the upper, middle and lower watershed of the river basin.

More importantly, land use changes, characterized by the transition from one type to another, were extremely prominent. From the land use conversion matrix between 1980 and 2015 in Table 2, it can be calculated that the total exchange area is 578.93 km<sup>2</sup> , or 24.71% of the total catchment area. Specifically, the conversions among cropland, forest, orchard, water, and construction land comprised 94.78% of the total exchange. Among them, construction land has increased by 158.67 km<sup>2</sup> or 127.47%. Cropland and forest have decreased by 146.33 km<sup>2</sup> or 20.00% and 39.24 km<sup>2</sup> or 3.08%, respectively. Other changes were all less than 20.00 km<sup>2</sup> . It can be observed that cropland and forests primarily contributed to land use exchanges and were the major land use types encroached on by urbanization. Of the 158.67 km<sup>2</sup> increases in construction land, 89.17% resulted from conversion of cropland and 10.35% resulted from conversion of forests.


**Table 2.** Land use types conversion matrix between 1980 and 2015 (km2).

Besides, there also exist great temporal and spatial differences in land use exchanges. Temporally, the Sankey diagram in Figure 3 visualizes exchanges of each land use type over different time periods. Therein, exchanges in almost all land use types during 2000–2010 are the most significant. Notably, the increase of construction land in this decade accounted for more than 50% of the total increase in the entire period, and 78.11% of the construction land increase came from conversion of cropland and 10.10% from conversion of forest. Spatially, Figure 4 shows that the conversions mainly occurred in the middle and lower reaches of the basin where a large amount of cropland was converted into construction land. Particularly, the lower reaches experienced the most drastic changes in the river basin regarding croplands and construction areas.

**Figure 3.** Comparison of exchanges of land use types in four time periods.

#### *3.2. Spatiotemporal Variations of Landscape Patterns*

The different landscape pattern indices were calculated at both landscape and land use class levels on the specific scale of 1000 × 1000 m. The former represents the overall spatial arrangement characteristics of each landscape patch, and the latter reflects the spatial arrangement characteristics of each landscape patch in different types. That is, the landscape pattern of each landscape patches at the class level can determine it at the landscape level, but not vice versa.

Figure 5 shows that changes of those indices are all obvious at the landscape level. For example, the landscape fragmentation index PD has an increasing trend and increased by 8.03%, while AI shows a decreasing trend and decreased by 2.25%. Meanwhile, the interference index AWMPFD and dominance index LPI decreased by 0.71% and 4.86%, respectively. The diversity index SHDI shows an increasing trend and increased by 5.25%. Temporally, these indices altogether indicated an increasing fragmentation and homogenization in landscapes, intensified human interference, and a weakening dominance of the

once-dominant landscape (forest and cropland) in the early years. Similar to the changes in land use types, the decade from 2000–2010 experienced the most significant changes in landscape pattern.

**Figure 4.** Spatial variations of land use type transition in different time intervals (1980–1990, 1990–2000, 2000–2010, 2010–2015 and 1980–2015).

**Figure 5.** Changes of landscape metrics (**a**) Shannon's diversity index (SHDI), (**b**) patch density (PD), (**c**) aggregation index (AI), (**d**) area-weighted mean patch fractal dimension (AWMPFD), and (**e**) largest patch index (LPI) at the landscape level from 1980 to 2015.

Moreover, similar to land use changes, landscape pattern also showed distinct spatial characteristics, reflecting the landscape heterogeneity in the LXH. As shown in Figure 6, the PD and SHDI are relatively small while AI, AWMPFD, and LPI are relatively large in the upstream area. However, this pattern reverses in the middle and lower watershed. This indicated that the landscape in the middle and lower watershed was more fragmented and heterogeneous than the landscape in the upper watershed, with a higher interference and lower dominance. Moreover, by comparing the values of these landscape metrics in 1980 and 2015 shown in subfigures c in Figure 6, we find that they are almost unchanged in most regions, except for part of the middle and lower watershed where cropland had been largely converted to construction land.

**Figure 6.** Spatial distribution of the changes in landscape metrics between (**a**) 1980 and (**b**) 2015, and (**c**) showed the results of subtraction of the indices in 1980 from 2015.

At the land use class-level, it can be found that the PD and AI value of cropland and forest changed more significantly than the other types during the total study period. The PD of cropland and forest increased by 42.26% and 41.59%, and the AI decreased by 1.68% and 0.21%, respectively. Meanwhile, the LPI values of cropland and forest decreased by 67.45% and 3.07%, respectively. The LPI of construction land in 2015 was 12 times that of the LPI in 1980, second only to the LPI of cropland. Among all land use types, the AWMPFD of cropland decreased most significantly, by 4.39% in total, while the net change of the AWMPFD of other types was less than 1.00%. That is, the degree of fragmentation, homogenization, and human alteration of cropland in this river basin was the most significant of all terrains. At the same time, the degree of landscape dominance of cropland was greatly reduced, while that of construction land was greatly improved. This indicated that the fragmentation and homogenization of cropland was mainly contributed by occupation of the construction land.

#### *3.3. Impacts of the Anthropogenic Factors on Temporal Landscape Changes*

Since the effect of natural factors on landscape changes was minor compared with anthropogenic factors in a short period of time, this paper selected eight anthropogenic factors involving demographic factors, socioeconomic factors, and urbanized activities to reflect their impacts on landscape changes temporally (refer to Figure 7 for detailed indices and their change trends). Meanwhile, due to the difficulty of quantitatively expressing the related policies, they were not included in the following quantitative analysis.

**Figure 7.** Temporal trends of influencing factors (**a**) Total population, (**b**) Proportion of nonagricultural population, (**c**) Gross domestic product, (**d**) Proportion of primary industry, (**e**) Proportion of secondary industry, (**f**) Proportion of tertiary industry, (**g**) Annual per capita income and (**h**) Total investment in real estate development during 1980 to 2015.

First, from the Figures 7 and 8, it can be found that the temporal changes of some influencing factors and land use types are non-linear; we chose the method of grey cor-

relation analysis in this context. According to the results from the relational analysis in Table 3, all the grey correlation coefficients are greater than 0.55, indicating that the changes of these eight anthropogenic factors are all significantly correlated with the changes of various land use types in the LXH temporally. More specifically, the four most correlated factors on each land use type were the TP, the PNAP, the PSI, and the PTI, which were demographic and socioeconomic factors and their correlations were all above 0.73. Among them, for cropland and each natural land (forest, shrub, grassland, water, floodplain, and unused land), the most correlated factors were PSI and TP; while for construction land, the most correlated factors were PNAP and PTI. From the Figure 7, it is clear that the TP, PNAP, and PTI all increase dramatically from 1980 to 2015 in the LXH. The PSI shows a trend of first increasing and then decreasing, but it also increases on the whole. Therefore, combined with the characteristics of the proportion changes in various land use types in Figure 8, it can be concluded that the increase of construction land was mainly correlated with the increase of non-agricultural population and the continuous development of the tertiary industry in the LXH. The decreasing of cropland and each natural land was mainly correlated with the increase of TP and the changes in secondary industry. This might be attributed to the growing population (especially the growth of urban population) and the transformation of industry (especially the growth of tertiary industry), which has accelerated the encroachment on cropland and natural lands to meet the demands for more construction land in the LXH [71,72]. Alternately, the other four influencing factors were also crucial to the changes of different land use types, but they had less of an effect compared with these main influencing factors.

**Figure 8.** Temporal trends of each land use type (**a**) Proportion changes of Cropland, (**b**) Proportion changes of Forest, (**c**) Proportion changes of Shrub, (**d**) Proportion changes of Grassland, (**e**) Proportion changes of Floodplain, (**f**) Proportion changes of Unused land, (**g**) Proportion changes of Orchard, (**h**) Proportion changes of Water and (**i**) Proportion changes of Construction land during 1980 to 2015.


**Table 3.** Grey correlation coefficients of influencing factors on land use changes.

#### *3.4. Impacts of Anthropogenic and Natural Factors on Spatial Landscape Changes*

Considering the fact that landscape pattern experienced the most significant changes from 2000 to 2010, we took this period as an example to study the influencing factors of spatial heterogeneity toward landscape pattern in this river basin. Table 4 manifests the significant spatial correlation between each landscape metric and each investigated influencing factor in 2000 and 2010 respectively. Among them, topographic elements DEM and SLOPE were all spatially negatively correlated with PD and SHDI, and positively correlated with AI, AWMPFD, and LPI. This indicated that in areas with low elevation and gentle slopes, the degree of landscape fragmentation, landscape interference, and landscape homogenization was stronger, and that the landscape dominance was weak. Table 4 also shows that all anthropogenic influencing factors (GDP, TP, and NLD) are positively correlated with PD and SHDI, and negatively correlated with AI, AWMPFD, and LPI spatially. This illustrates that the degree of landscape fragmentation, landscape interference, and landscape homogenization is relatively strong in more developed regions. In terms of the impact of land use changes brought by rapid urbanization on landscape patterns, the LUIN was positively correlated with PD and SHDI spatially and negatively correlated with AI, AWMPFD, and LPI. This reflected that areas with high LUIN were usually accompanied with a relatively stronger degree of landscape fragmentation, landscape interference, and landscape homogenization.

**Table 4.** Bivariate Moran's I correlation analysis between landscape metrics and influencing factors in the spatial dimension in 2000 and 2010.


Note: Permutation test was used to test in this study, and the P value of each group of variables was equal to 0.001, indicating that the spatial correlation was significant under 99.9% confidence.

Geographic detector analysis showed that the interpretation of DEM on the spatial distribution of each landscape metric was the largest among all influencing factors in 2000 or 2010, with the highest average q value being 0.28 or 0.23, followed by GDP, TP, and NLD, while Slope had the lowest q value (Figure 9). Results indicated that the spatial distribution of elevation was the key factor that induced the spatial heterogeneity of landscape pattern, and the spatial distribution of socioeconomic level, population density, and urbanized activities also played an important part. The analysis of interactions between the influencing factors and landscape metrics (Figure 10) showed that the interpretation of spatial distribution characteristics of landscape pattern by any two influencing factors was greater than that of any single influencing factor, indicating that the formation of spatial heterogeneity was the result of interactions between various influencing factors. Specifically, from the subfigures a in Figure 10, the interaction between DEM and LUIN is

the strongest among other factors. The interaction between LUIN and other four influencing factors on each landscape metric was almost stronger than that between any other two influencing factors in 2000. This indicated that the spatial differences of DEM and LUIN jointly resulted in the spatial heterogeneity of landscape pattern in 2000, stronger than the interactions between DEM and GDP, TP, NLD, Slope comparatively. However, compared with the results from 2000, since the interaction between DEM and TP, GDP, and NLD was strengthened, the interaction between DEM and LUIN was no longer the strongest among the other factors in 2010 (see subfigures b in Figure 10). These indicated that the spatial differences of the topographic elements and other influencing factors also jointly contributed to the spatial heterogeneity of landscape pattern in 2010.

**Figure 9.** The force q among each influencing factor on each landscape pattern metric in (**a**) 2000 and (**b**) 2010.

**Figure 10.** The force q among any two influencing factors on each landscape pattern metric interactively in (**a**) 2000 and (**b**) 2010.

#### **4. Discussion**

#### *4.1. Spatiotemporal Changes of Land Use and Landscape Pattern*

The results of this study showed that in the LXH, there exist large spatial and temporal differences in land use changes and landscape pattern changes. These changes appeared to be more prominent in the middle and lower watershed, and their changing rates were fastest during 2000 to 2010. Specifically, the land use change was featured by the increasing transition of cropland and forest to construction land, and the fragmentation and homogenization of landscape pattern was contributed to the encroachment of construction land on forest and cropland. That is, the decrease of cropland and forest was accompanied with the decreased degree of the cropland and forest landscape dominance and the increased degree of the cropland and forest landscape fragmentation and homogenization. These further proved the synchronization characteristics and interaction relationship between land use

changes and landscape pattern changes proposed by scholars [17,19]. Thus, it is possible for researchers to use the temporal change of a land use type to reflect the temporal change pattern of a certain landscape type in a period of time.

Besides, our findings of the land use and landscape pattern changes are consistent with the previous research in the whole Guangzhou city. For example, Zhang et al. [73] and Gong et al. [74] also found that the increase of construction land in new urban areas of Guangzhou city mainly came from cropland, forest, and other ecological land, especially after 2000. Gong et al. [41] also confirmed that the fragmentation and homogenization of cropland in Guangzhou was mainly contributed by the expansion of construction land. However, their research paid more attention to the urbanization expansion pattern of Guangzhou by comparing the changing differences of certain land use types and landscape patterns in different jurisdictions, instead of focusing on these spatiotemporal changes brought by urbanization [38,75]. Researches on the analysis of land use and landscape pattern changes in a watershed under the urbanization expansion pattern also exist, which provide theoretical basis and method reference for this study. But watersheds selected in their studies are relatively large, spanning multiple cities [9,21]. As a case study in this paper, the LXH was relatively small and in the range of Guangzhou city. Its middle and lower watershed is adjacent to the central urban area of Guangzhou, while the upstream area is far away from the central urban. This pushed the gradual widening of the difference between the northern and southern parts of the watershed influenced by urbanization. Our results also found that changes of land use and landscape pattern were different between the northern and southern parts.

Moreover, analysis results above reflected that the time period from 2000 to 2010 and the southern parts of the LXH with the most prominent changes should be taken seriously by relevant stakeholders. First, changes of land use and landscape pattern in southern parts of the LXH should be slow down and controlled, and the northern parts should be protected timely under the rapidly urbanizing trends. Then, the special time period from 2000 to 2010 needs to pay much attention to in related researches about the LXH. It means that this study gives not only a supplement to previous studies in these regions, but also is of great value for managers, planners, and scholars to make appropriate strategies.

#### *4.2. The Temporal and Spatial Influencing Factors*

In terms of the influencing factors of the changes of land use types and landscape patterns, previous studies mainly discussed the reasons for land use type conversion at different locations [32,76] and in different time periods [21,29], but few analyzed the factors responsible for the spatial heterogeneity of landscape patterns in river basins, nor did they comprehensively quantify the factors that contributed to the spatiotemporal change of land use types and landscape patterns. In this study, considering the fact that various land use type changes emphasized the transition of different landscape patches in different time periods, and that changes of landscape patterns reflected the difference of spatial configuration characteristics in different landscape patches, we analyzed the influencing factors on the temporal change of land use types and spatial heterogeneity of landscape pattern, respectively. We found that there was a greater difference in the spatiotemporal influencing factors of land use and landscape pattern changes in the LXH. Thus, it is very important to propose a targeted protection and development strategy, which can meet the current needs of the different regions in the LXH. Temporally, we found that the demographic factors, socioeconomic factors, and urbanized activities were important in shaping the temporal variations of land use types in this river basin, and changes of major land use types were more sensitive to the increasing of non-agricultural population and transformation of industry than any other factors. This was consistent with the findings in other studies on the influencing factors of land use types changes in other river basins [15,28]; they also found that the growth of urban population and changes of industries contributed to the increase of construction land [4,77]. Moreover, similar studies elsewhere underlined that government policy also played an important

role in the change of land use types during different periods [78–80]. Although there was no appropriate method for analyzing the impact of the related policies on land use changes quantitatively, we found that various land use types in the LXH have undergone significant changes in the last three decades after 1980, and these changes were particularly dramatic after 2000. Zhang et al. [73] found that the implementation of China's reform and development policy in 1978 was an important driving force for economic development and population migration of Guangzhou, pushing the continuous expansion of construction land to gradually occupy the cropland, forest, shrub, and grassland in the suburbs. Thus, the observed expansion of the construction land in the LXH after the early 1980s may be attributed to the implementation of this policy. In addition, the overall urban development master plan of Guangzhou in 2000 had put forward the strategy of expanding the urban area to the north and built Guangzhou into an international metropolis by 2010, which could further accelerate the expansion of construction land if practices continue. Correspondingly, the land use types in the LXH changed significantly after 2000 compared with the pre-2000 practices, accounting for more than 50% of the total variation in 35 years. Moreover, Baiyun, Huadu, and Conghua districts in the north had successively merged into the jurisdiction of Guangzhou in the year 2000, 2010, and 2015, respectively. The different speed of urbanization in different regions altered the variation characteristics of land use types. We also found that the lower watershed that contains the Baiyun and Huadu districts had the largest proportion of cropland conversion to construction land, which was five times higher than that of the upper and middle watershed. Therefore, apart from the demographic factors, socioeconomic factors, and urbanized activities, the relevant government policy, which is difficult to quantify, also significantly affected the variations of land use types.

Focusing on the impact of various influencing factors on landscape pattern changes in spatial dimension will be very useful in identifying and controlling the major driving forces, guiding the watershed protective management and sustainable planning. However, most of the current studies adopted the classification method to describe the spatial characteristics of landscape pattern and their influencing factors of different regions, and seldom analyzed the spatial relationships among variables quantitatively [9,33]. In the research of Ju et al. [70], the applicability of the geographic detector model in analyzing the driving force of construction land expansion was proved, providing a quantitative method for the analysis of the interaction among various spatial factors. But their research did not conduct a comparative analysis of the spatial driving relationships in different periods. Here, using the model of geographic detector by Wang et al. [69], this study compared the impact of different influencing factors on spatial landscape heterogeneity during 2000 and 2010, when the most dramatic land use changes happened. It can be found that the spatial distribution of LUIN and elevation were the two critical factors for the formation of landscape heterogeneity in 2000 compared with other factors, while the interaction between elevation and other human factors was strengthened in 2010; this illustrated that elevation was always a basic factor that directly determined the spatial distribution of landscape pattern. Liu et al. [61] also proposed that it was difficult for people to break through existing natural obstacles in the hilly regions of southern China, and this difficulty had largely restricted human activities. From Figure 11, it is clear that great differences in the terrain conditions exist in this study. Because areas with high elevation or greater slopes were difficult to develop and not suitable for urban construction, they were seldom disturbed by human activities [33], so that the degree of landscape fragmentation and homogenization in the upper reaches was low and the degree of landscape dominance was high. Therefore, elevation was a prerequisite for the impact of anthropogenic factors to occur. On the other hand, based on the difference of elevation in different regions of LXH, the human influencing factors, such as population and socioeconomic and urbanized activities, played an increasing role in the formation of the heterogeneous characteristics of the landscape pattern after 2000. This may be due to the fact that the spatial difference

of human influence factors increased significantly after 2000, as it shown in subfigure c in Figure 11.

**Figure 11.** Spatial distributions of the selected influencing factors selected in this study in (**a**) 2000 and (**b**) 2010, and (**c**) shows the results of the subtraction of each influencing factor between 2000 and 2010.

The above discussion illustrates that it is necessary to control the increasing trends of non-agricultural population and the continuous development of secondary and tertiary industry in the future, thus the demand for more construction land will be decreased. Meanwhile, relevant policies should try to meet these demands. We also should pay much attention to the southern parts of the LXH, and strengthen its adjustment ability to deal with the intensive population density, higher GDP, and greater urban construction. For example, the urban occupancy rate in these areas can be increased, artificial green land can be increased, the native forest and grassland must be strictly protected, etc. This means that establishing the spatiotemporal change trends and causes of land use and landscape pattern in a rapidly urbanizing watershed is very important for guiding the diagnosing of urbanization problems, clarifying the main protection areas and main control factors.

#### *4.3. The Limitations and Potential Outlooks*

This study produces a quantitative estimate of the spatiotemporal variations in land use types and landscape patterns and analyzes the dominant influencing factors leading to these changes in LXH quantitatively, which provides a systematic integration and deepening of previous studies. The main land use maps used in this study interpreted by the common method of visual interpretation, and their errors mainly came from the personal subjective judgment of the interpreters and the similarity of the tones and textures of the satellite image. Although these errors in the interpretation process were considered and improved through some reference maps (including topographic maps, vegetation maps, ground truth data at different sample points, and local resident interview data), there also existed uncertainties in the data measurement and description [36]. These errors will also affect the accuracy of the research results to some degree. Therefore, it is necessary to compare any two measurement methods to improve the analysis accuracy of land use maps and express them on different scales as much as possible in the future. Besides, due to the scarcity of historical landscape information, e.g., land use maps from 1980 to 2015, Topography data and Nighttime Light Image Data in 2015, we were unable to accurately establish the relationship between each land use type and different influencing factors, nor could we compare the driving factors of spatial heterogeneity of landscape patterns in each period. In addition, there are some factors that cannot be quantified, such as policy factors, which make the analysis of influencing factors still not comprehensive enough. In the future, it may be possible to construct a comprehensive model combining qualitative and quantitative analysis on all possible influencing factors of land use changes.

Resolution of available spatial data set of spatial influencing factors was also a limitation. We conducted the driving analysis between the influencing factors and landscape indices at a 1 km grid. Although on this scale the influencing factors also had good spatial correlation, the resolution limitation of influencing factors might affect the accuracy of the analysis results to some degree. Therefore, spatial data with appropriate accuracy at a higher resolution and longer periods could substantially improve the accuracy when analyzing the change of landscape patterns and their driving forces. In that sense, models of land use change would be a good alternative for such studies by simulating the land use in different years to increase the range and length of land use data [81], and hence guide the urbanization development of this region by analyzing the change of land use types and landscape patterns in the future. Moreover, the extensive establishment of the real-time monitoring data platform of different spatial influencing factors such as social economic activities, population density, and urbanization activities in the future with higher resolution will improve the accuracy of spatial analysis, thereby realizing its dynamic analysis.

Moreover, based on our comprehensive analysis of the spatiotemporal changes and causes of landscape pattern in the LXH and its ecological and hydrological effects in related researches [43], it is more urgent to establish specific strategies to guide the sustainable development of LXH in the future. The Hellwig classification and measurement method introduced by Hellwig in 1968 provides a decision-making method for formulating sustainable development strategies based on the evaluating of urban development [82]. This method was first applied in the sustainable decision-making process of urban green space biodiversity management in Lublin, eastern of Poland, thus the main ecological areas that should be protected can be established [83]. Then, other scholars used and extended this method at different scales in European Union to formulate sustainable development strategies based on the different goals. These application of the Hellwig method in different researches prove its effectiveness in evaluating the level of development of different regions in different fields at different scales, which provide a new direction for the future research of establishing the sustainable development strategy in the LXH based on the analysis of its changes about land use, landscape pattern, hydrological and ecological conditions under the rapid urbanization. After that, the specific areas of the LXH under the rapid urbanization process, in which its land use transition and landscape pattern fragmentation should be extremely controlled can be found.

#### **5. Conclusions**

In this study, we analyzed the spatiotemporal changes of land use types and landscape pattern of the LXH from 1980 to 2015 under the rapid urbanization of Guangzhou city, as well as quantified the major influencing factors temporally and spatially. The main conclusion can be concluded as one sentence that there exist great spatiotemporal differences in land use and landscape pattern changes and its causes in the LXH during the past 35 years. Specifically, it can be drawn as follows:

• The most obvious land use change was characterized as the large transition from cropland to construction land, bringing about the fragmentation of cropland that was encroached on by the construction land. The landscape pattern showed an increasing trend of landscape fragmentation, homogenization, and landscape interference, and a decreasing trend in landscape dominance. These changes mainly occurred in the lower watershed, particularly between 2000 and 2010. Therein, these changes were more than 50% in this decade compared with total 35 years.


The findings are of great significance for review and outlook of the ecological protection and sustainable development of the watershed around the rapidly urbanizing areas. It can not only allow decision-makers to clarify their main problems, but also guide them to clarify the key protection areas and control indicators. However, the analysis of the landscape patterns above was limited to the period from 1980 to 2015, and the comparison of influencing factors on spatial landscape configurations focused only on 2000 and 2010. Nonetheless, results in this study are insightful, although they could be more generalized with the analysis over a longer period.

**Author Contributions:** Conceptualization, Z.Z. and B.L.; methodology, Z.Z.; software, Z.Z.; validation, Z.Z., H.W., and M.H.; formal analysis, Z.Z. H.W., and M.H.; investigation, Z.Z. and B.L.; resources, B.L.; data curation, B.L.; writing—original draft preparation, Z.Z.; writing—review and editing, Z.Z., B.L., H.W., and M.H.; visualization, Z.Z.; supervision, B.L.; project administration, B.L.; funding acquisition, B.L. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the National Natural Science Foundation of China (Grant Nos. 51879289), the Guangdong Basic and Applied Basic Research Foundation (2019B1515120052) and the Innovation Fund of Guangzhou City water science and technology (GZSWKJ-2020-2).

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** The data presented in this study are available on request from the corresponding website.

**Acknowledgments:** We would like to thank the editor Madalina Buzatu and three anonymous reviewers for their constructive comments.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **Appendix A**


**Table A1.** Description of land use types in the LXH.

Table A1 gives a detailed description of the content about each land use type in this study, which can better display the classification standard of land use types.

**Table A2.** Analysis results of spatial autocorrelation of each landscape metric in different scales in 2000 and 2015 of LXH


Note: Permutation test was used to test in this study, and the P value of each landscape metrics in different scales was equal to 0, indicating that the spatial correlation was significant under 99.9% confidence. The Z value of them were all >1.96, reflecting that there exists extremely significant spatial autocorrelation among these landscape metrics in different spatial scales.

Table A2 reflects the degree of spatial autocorrelation about each landscape metric selected in this paper. It confirms that the spatial autocorrelation of these landscape metrics is extremely significant in different spatial scales ranging from 500 to 1200 m. Therefore, it can be proved that these spatial scales are all appropriate for analyzing the bivariate spatial correlation between each landscape metric and each influencing factor.

#### **Appendix B**

This Figure A1 is provided to screen the notable landscape metrics from 23 frequently used metrics. It demonstrates the correlation between any two kinds of landscape metrics among these 23 metrics above. Clearly, most of them were highly correlated with each other. Thus, when the absolute value of correlation coefficient between two indices is more than 0.9, only one is used; second, indices representing different aspects of landscape characteristics were selected to reduce the information redundancy among them. Finally, five representative indices were selected in this paper, including patch density (PD), aggregation index (AI), largest patch index (LPI), area-weighted mean patch fractal dimension (AWMPFD), and Shannon's diversity index (SHDI), which represent different aspects of landscape characteristics.

**Figure A1.** Results of the factor analysis among 23 common metrics.

Figure A2 presents the changing value of five selected landscape metrics at landscapelevel in 14 different granularities (including 30 m, 50 m, 100 m, 200 m, 300 m, 400 m, 500 m, 600 m, 700 m, 800 m, 900 m, 1000 m, 1100 m, 1200 m). It proves that 500~1200 m was the common characteristics interval of these landscape metrics. Thus, the spatial scale for analyzing the spatial autocorrelation of each landscape metric has been established.

**Figure A2.** Analysis results of the characteristic scale interval in 14 different granularities. **Figure A2.** Analysis results of the characteristic scale interval in 14 different granularities.

#### **References References**

	- driving forces: Implications for urban-rural coordination development in Wuhan, Central China. *Habitat Int.* **2015**, *47*, 113–125, doi:10.1016/j.habitatint.2015.01.012. 5. Meyfroidt, P.; Lambin, E.F.; Erb, K.H.; Hertel, T.W. Globalization of land use: Distant drivers of land change and geographic displacement of land use. *Curr. Opin. Environ. Sustain.* **2013**, *5*, 438–444. [CrossRef]
	- Hanoi Central Area, Vietnam under historical and extreme rainfall conditions. *Sci. Rep.* **2018**, *8*, 12623, doi:10.1038/s41598-018-30024-5. 8. Edge, C.B.; Fortin, M.J.; Jackson, D.A.; Lawrie, D.; Stanfield, L.; Shrestha, N. Habitat alteration and habitat fragmentation differentially affect beta diversity of stream fish communities. *Landsc. Ecol.* **2017**, *32*, 647–662. [CrossRef]

coastal watershed of southeastern China. *Glob. Ecol. Conserv.* **2020**, *23*, e01177, doi:10.1016/j.gecco.2020.e01177.

*Landsc. Urban Plan.* **2004**, *69*, 335–354, doi:10.1016/j.landurbplan.2003.10.033.

*viron. Manag.* **2018**, *61*, 954–967, doi:10.1007/s00267-018-1019-4.

9. Tang, J.; Li, Y.; Cui, S.; Xu, L.; Ding, S.; Nie, W. Linking land-use change, landscape patterns, and ecosystem services in a

10. Xiao, J.; Shen, Y.; Ge, J.; Tateishi, R.; Tang, C.; Liang, Y.; Huang, Z. Evaluating urban expansion and land use change in Shijiazhuang, China, by using GIS and remote sensing. *Landsc. Urban Plan.* **2006**, *75*, 69–80, doi:10.1016/j.landurbplan.2004.12.005. 11. Li, X.; Yeh, A.G. Analyzing spatial restructuring of land use patterns in a fast growing region using remote sensing and GIS.

12. Li, C.; Zhang, Y.; Kharel, G.; Zou, C.B. Impact of Climate Variability and Landscape Patterns on Water Budget and Nutrient Loads in a Peri-urban Watershed: A Coupled Analysis Using Process-based Hydrological Model and Landscape Indices. *En-*

