Spatial-Temporal Driving Factors of Urban Landscape Changes in the Karst Mountainous Regions of Southwest China: A Case Study in Central Urban Area of Guiyang City

: Rapid urban expansion has signiﬁcantly altered the regional landscape pattern, posing a serious threat to the sustainable development of natural and social ecosystems. By using landscape patterns indices and an area transfer matrix, this study analyzed the spatial-temporal changes of landscape patterns in the karst mountainous cities of southwest China from 2000 to 2020, by taking the central urban area of Guiyang City (CUAG) as the study area. This study explored the spatial and temporal driving factors of landscape pattern changes by using stepwise multiple linear regression and geographic detector methods. The results show: (1) CUAG’s landscape types altered changed drastically, with the area of forestland and construction land rapid increment and cultivated land decrement signiﬁcantly. (2) The patches of construction land and forestland tended to be aggregated, the degree of fragmentation was reduced, and the shape was complex; cultivated land fragmentation was intensiﬁed. The connectivity of the landscape was improved, while the level of landscape diversity declined, the trend of landscape homogenization was obvious. (3) Socioeconomic and geographical endowment drivers have determined landscape pattern changes. The ﬁndings of this study may be used to interpret other similar landscapes worldwide and may imply the protection of urban ecosystem and sustainable development.


Introduction
Landscape patterns, which are defined as spatial configurations of landscape patches, can assist in quantifying changes in the compositions and arrangements of land features with landscape metrics [1][2][3]. Understanding the internal mechanism of landscape patternprocess-scale-service, supporting the construct of landscape ecological security patterns, and encouraging the long-term development of regional social economies require research into landscape pattern alterations [4][5][6]. Over the last few decades, urbanization has been one of the most irreversible human activities driving worldwide landscape pattern change [7][8][9]. It is also a direct manifestation of the interaction between economic and social conditions and natural conditions [10,11]. The word "urban expansion" causes a rise in the area of urban construction land as a result of the overall influence of the urbanization process and changes in urban landscape type [12,13]. However, uncontrolled spatial growth during the urbanization process has resulted in a slew of issues, such as forest destruction, waste of land resources, the disappearance of native ecosystems, and biodiversity threats [14,15]. Landscape pattern heterogeneity is a prevalent geographical phenomenon in urban expansion that increases landscape fragmentation, homogenizes urban biodiversity, and creates a serious danger to the stability and ecological security of the urban landscape ecosystem [16][17][18][19]. Analysis of landscape pattern changes and their driving forces could aid in coordinating the relationship between urbanization and landscape ecosystem protection, as well as promote the sustainable development of both the urban social economy and the landscape ecosystem [20]. Especially when ecological conservation has become a global advocacy issue, the sustainability of the landscape is an important prerequisite for the development of the economy [21,22], and it is important to explore the driving force of the landscape pattern under the combined effect of landscape conservation and economic development for planning urban ecosystems, building ecological civilization and realizing ecological economic development [21].
Several scholars have studied the landscape pattern changes and driving forces of urbanization [23][24][25]. The relevant studies have been conducted to assess regional landscape patterns and driving forces in intercontinental [26], national [27,28], provincial [29], and civic [30,31] contexts. In general, using various landscape indices to measure the spatial-temporal heterogeneity of urban landscape patterns and to quantitatively reveal the landscape components and distribution characteristics is the conventional paradigm of landscape pattern analysis [32][33][34]. In addition, numerous methods and models for detecting the driving forces of the landscape pattern evolution from natural, social, economic vantages, as well as policy factors, have been applied, such as principal component analysis [35], correlation analysis [32], and linear regression analysis [31]. Although the aforementioned statistical analysis methods can investigate the relationship between the characteristics of urbanization landscape patterns and the influencing factors, there is a lack of ability to directly detect the heterogeneity of urban expansion landscape patterns. The geographic detector is a theoretical method and technique for detecting spatial heterogeneity and revealing the factors driving it [20]. It is feasible to calculate the impact and geographical variation of various elements on geospatial phenomena by using the models that have been built [20]. The geographic detector provides many advantages and has been applied widely in landscape pattern driving force exploration [36]. However, there are few studies on the spatial-temporal changes of landscape patterns in karst mountainous cities under the rapid urbanization; meanwhile, the driving mechanism from the perspectives of spatial and temporal scales to quantify evaluation of urban landscape pattern changes is still weak [37].
The karst area in southern China, centered on the Guizhou Plateau, is the site of the most typical, complex, and diverse landscape forms of karst development in the world, as well as the largest and most densely populated ecologically vulnerable area [38]. Owing to its special geomorphological features, heavy population pressure, and limited land resources, the region is extremely susceptible to human activities and environmental changes [39,40]. Throughout this region, the ecological environment consists of complicated geomorphological types with wide spatial distribution patterns, all of which have a direct impact on land use conditions and landscape pattern design. Since economic development, the karst mountainous region has witnessed dramatic urbanization, resulting in radical changes in urban landscape pattern and a series of ecological-resource environmental problems that have become an urgent proposition to be solved for the optimization of landscape spatial pattern, enhancement of ecosystem services, and sustainable socio-economic development in the karst mountainous region. The spatial-temporal changes of urban landscape patterns and their driving mechanisms in the context of coupled rapid urbanization and ecological civilization city construction are the primary scientific questions to decipher the relationship between karst mountain urbanization development and ecological environmental protection. Since the 21st century, Guiyang city has accelerated its urbanization process, with an urbanization expansion intensity of 18.14% from 2005 to 2015, entering a stage of rapid urbanization, while at the same time, the strategic goal of building an ecological civilization city was set in 2007, and an array of ecological environmental protection and restoration measures were subsequently implemented. Under the dual background of rapid urbanization and ecological civilization city construction, urban development and ecological protection were mutually coupled, which led to drastic changes in the urban landscape pattern of Guiyang [41]. However, what are the patterns of urban landscape Sustainability 2022, 14, 8274 3 of 18 pattern evolution? What are the dominant driving factors in different spatial and temporal scales? These scientific questions are still lacking in systematic research, and there is a shortage of scientific support for optimizing the spatial patterns of the urban landscape in Guiyang and even in karst mountains in the new era, enhancing urban landscape ecosystem services, promoting sustainable landscape development, and consolidating the achievements of ecological civilization city construction.
By using landscape indices, the stepwise multiple linear regression model, and the geographic detector method, studies explored the spatial and temporal driving forces of landscape pattern changes in the central urban area of Guiyang City (CUAG). The targets of this study were as follows: (1) to examine the spatial and temporal changes in various landscape types from 2000 to 2020 as a result of rapid urbanization; (2) to identify the evolution of landscape patterns in different periods; (3) to reveal the dominated driving forces for landscape pattern changes at the spatial and temporal scales, respectively.

Study Area
The central urban area of Guiyang City (CUAG), a typical mountainous city in southwest China (106 • 30 -106 • 55 E, 26 • 21 -26 • 48 N), is located in the middle section of Guizhou Province, with a total area of around 1178 km 2 ( Figure 1). The region has the most typical and complex landscape types of karst development in the world [38]. The elevation of the study area is ranged from 982 m to 1589 m with the average elevation of 1200 m, and the environment is a subtropical humid and mild climate with average annual temperatures of 15.3 • C and relative humidity of 77%. Guiyang City has grown rapidly in recent decades, with a total urban area of 8034 km 2 and a population of 4.32 million in 2010 [42]; the intensity of urbanization expansion reached 18.14% during 2005-2015. Guiyang is an illustration of a resource-dependent city with abundant natural resources. Large-scale resource-dependent and process industries are formed as a result of such resource quality, and Guiyang was designated as the national circular economy pilot in 2005 [42]. changes in the urban landscape pattern of Guiyang [41]. However, what are the patterns of urban landscape pattern evolution? What are the dominant driving factors in different spatial and temporal scales? These scientific questions are still lacking in systematic research, and there is a shortage of scientific support for optimizing the spatial patterns of the urban landscape in Guiyang and even in karst mountains in the new era, enhancing urban landscape ecosystem services, promoting sustainable landscape development, and consolidating the achievements of ecological civilization city construction. By using landscape indices, the stepwise multiple linear regression model, and the geographic detector method, studies explored the spatial and temporal driving forces of landscape pattern changes in the central urban area of Guiyang City (CUAG). The targets of this study were as follows: (1) to examine the spatial and temporal changes in various landscape types from 2000 to 2020 as a result of rapid urbanization; (2) to identify the evolution of landscape patterns in different periods; (3) to reveal the dominated driving forces for landscape pattern changes at the spatial and temporal scales, respectively.

Study Area
The central urban area of Guiyang City (CUAG), a typical mountainous city in southwest China (106°30′-106°55′ E, 26°21′-26°48′ N), is located in the middle section of Guizhou Province, with a total area of around 1178 km 2 ( Figure 1). The region has the most typical and complex landscape types of karst development in the world [38]. The elevation of the study area is ranged from 982 m to 1589 m with the average elevation of 1200 m, and the environment is a subtropical humid and mild climate with average annual temperatures of 15.3 °C and relative humidity of 77%. Guiyang City has grown rapidly in recent decades, with a total urban area of 8034 km 2 and a population of 4.32 million in 2010 [42]; the intensity of urbanization expansion reached 18.14% during 2005-2015. Guiyang is an illustration of a resource-dependent city with abundant natural resources. Large-scale resource-dependent and process industries are formed as a result of such resource quality, and Guiyang was designated as the national circular economy pilot in 2005 [42].

Data Acquisition and Processing
Three major types of data sources are required for this study. One was the remote sensing images of the CUGA, which included Landsat-5 TM images for 2000, 2005, and 2010, as well as Landsat-8 OLI images for 2015 and 2020, which were provided by the International Scientific and Technical Data website of the Chinese Academy of Sciences (http://www.gscloud.cn/ (accessed on 28 January 2022)) with the cloud content less than 5%. The CGCS2000_3_Degree_GK_CM_108E projection coordinate system was used for all images, with a spatial resolution of 30 m × 30 m. All remote sensing images were pre-processed with the ENVI 5.3 software platform (ITT Visual Information Solutions, Herndon, VA, USA), which included radiation calibration, atmospheric correction, geometric correction, image cropping, etc. Additionally, according to the China's land use classification standard (GB/T 21010-2017), the landscape types were classified into cultivated land, forestland, grassland, water body, construction land, and unused land by using the support vector machine supervision classification method. The classification accuracy of each landscape type and the overall classification accuracy of each period were more than 85%, with the kappa index greater than 0.85.
The second was to obtain the digital elevation model (DEM) data provided by the International Science and Technology Information Website of the Chinese Academy of Sciences (http://www.gscloud.cn/ (accessed on 28 January 2022)). The GDP and population density data belong to the kilometer grid data; vector data of the river, railway, and highway were available at the Data Center for Resource and Environmental Sciences of the Chinese Academy of Sciences (https://www.resdc.cn/ (accessed on 28 January 2022)). Based on ArcGIS 10.5 software (Environmental Systems Research Institute, Inc., Redlands, CA, USA), the surface analysis tool was used to obtain the slope (x 3 ) and elevation (x 4 ) of the study area, and the neighborhood analysis and raster calculator were used to extract the topographic relief (x 5 ). To facilitate the calculation of the spatial relationship between the landscape types of forestland, cultivated land, and construction land, as well as the distances to river (x 6 ), railway (x 7 ) and highway (x 8 ), respectively, the spatial neighborhood analysis method was used on ArcGIS 10.5 platform to obtain the distances of the three landscape types grid transformed into river, railway, and highway, and then the spatial relationship between each driver and the three landscape types was evaluated one by one with the spatial fusion connection method. Among them, the distance to river reflects the accessibility of the geographic endowment of the study area, and the distance to railway and highway reflects the accessibility factors of the study area. A total of 8 factors in 3 dimensions form the spatial drivers of CUAG, as shown in Table 1. The resolution of all above data was 30 m × 30 m, and reproject all coordinate systems to the CGCS2000_3_Degree_GK_CM_108E projection coordinate system.

Dimension Name Description
Economic level Gross Domestic Product, it is the result of the productive activity of all resident units in a country (or region) during a certain period of time.
Per capita GDP (x 2 ) The ratio of the GDP of a region to the resident population of the region during the accounting period.

Indices of gross domestic product (x 3 )
It is a relative number that reflects the trend and extent of changes in GDP over a certain period of time.

Total investment in fixed assets (x 4 )
It is the workload of construction and acquisition activities of fixed assets expressed in monetary terms.
Total retail sales of consumer goods (x 5 ) Refers to the amount of physical goods sold by enterprises (units) to individuals and social groups through transactions, not for production or business use, and the amount of revenue earned from the provision of food and beverage services.
Local financial revenue (x 6 ) Consists of the region's fiscal revenues. Local fiscal revenues include local budget revenues and extra-budgetary revenues.

Population size
Population density (x 7 ) Number of people per unit of land area.
Natural population growth rate (x 8 ) Ratio of natural increase in population in a given period to the average total population in the same period.

Dimension Name Description
Industrial structure Ratio of primary industrial structure (x 9 ) Value added and employment in the primary industrial as a share of GDP and total labor force.
Ratio of secondary industry structure (x 10 ) Value added and employment in the secondary industry as a share of GDP and total labor force.
Ratio of tertiary industry structure (x 11 ) Value added and employment in the tertiary industry as a share of GDP and total labor force.

Dimension Name Description
Socioeconomic factors GDP (y 1 ) Kilometer grid data to measure the level of economic development of a region.
Population density (y 2 ) Kilometer grid data, indicating the population density of the area.

Natural factors
Slope (y 3 ) The degree of steepness of the surface unit, the ratio of the vertical height of the slope to the distance in the horizontal direction.
Elevation (y 4 ) The vertical distance above sea level at a point on the ground.

Topographical relief (y 5 )
The difference between the elevation of the highest point and the elevation of the lowest point in a particular area.

Accessibility factors
Distance to river (y 6 ) Distance between urban construction land and river.
Distance to railway (y 7 ) Distance between urban construction land and railway.
Distance to highway (y 8 ) Distance between urban construction land and highway.
Finally, combined with the research experience of related scholars and based on the principle of scientific through data acquisition, the socio-economic factors in terms of economic level, population size, and industrial structure of the CUAG were selected as temporal driving factors from the China Urban Statistical Yearbook and Guiyang Statistical Yearbook, as shown in Table 1.

Landscape Metrics Selection
Landscape metrics are widely used to describe how processes and spatial patterns are related. Four landscape metrics at the class level (patch density (PD), largest patch index (LPI), number of patches (NP), area-weighted mean fractal dimension index (FRAC AM)) and four landscape metrics at the landscape level (aggregation index (AI), contagion (CONTAG), Shannon's diversity index (SHDI), and Shannon's evenness index (SHEI)) were chosen based on the relevant literature (Table 2) [35,43]. The landscape metrics were generated in Fragstats 4.2 program by using landscape type datasets from various periods.  The structure and properties of landscape type changes, as well as the direction of change for each landscape type, might be visualized with a dynamic change transfer matrix [44]. This study utilized the area transfer matrix approach to show the spatial and temporal changes of each type of landscape in CUAG throughout time. The formula of the transfer matrix is as follows [45]: where, P represents the area of each landscape type in CUAG, n represents the total number of landscape types, and i and j represent the landscape types at the beginning and end of the study period. Each row element of the transfer matrix represents the flow of landscape type information, and each column element represents the source information.

Stepwise Multiple Linear Regression
The stepwise multiple linear regression (SMLR) model could remove variables that caused multicollinearity and determine the independent variables that have a significant effect on the target object [46]. In light of previous studies, SMLR was developed to reveal the main temporal factors (Table 1) influencing changes in urban landscape patterns of CUAG [28,47]. First, to overcome the effects of variation in dimensionality and to guarantee that the results reflect the significance of the independent variables, the values of the 11 temporal drivers on three dimensions needed to be normalized in the range of [0, 1] before conducting the SMLR, a step that was performed in SPSS 22.0 software. SMLR was conducted with the three main landscape types (construction land, forestland, and cultivated land) area of different periods as the dependent variables, and the socioeconomic factors of the corresponding periods as the independent variables. The model was expressed as [48]: where y is the landscape type (construction land, forestland, and cultivated land) area of CUAG; β 0 is the constant term of the model, β p is the regression coefficient, x p is economic level, population size or industrial structure factor.

Path Analysis
The method of path analysis was first proposed by Sewall Wright in 1921 [49,50]. The path analysis was able to effectively indicate the direct effect of the variable of interest on the dependent variable and to estimate the indirect effect of the independent variable on the dependent variable [49]. For an interrelated system with a dependent variable y and independent variables x i (i = 1, 2,..., n) with a linear relationship between them, the regression equation as follows: By bringing the actual values into Equation (3) and understanding the system of equations using the principle of least squares, it is possible to find the pass coefficient P yxj . The path coefficient is a biased regression coefficient normalized by the variables and indicates the relative importance of each cause on the outcome. Equation (3) can be mathematically transformed to create the regular matrix equation: where, r xiyj is the simple correlation coefficient of x i and x j , r xjy is the simple correlation coefficient of x i and y. By solving Equation (4), we can find the path coefficient P yxi . P yxi is: where, b i is the partial regression coefficient of y on σ i , σ xi , and σ y are the standard deviations of x i and y, respectively. P yxi denotes the direct path coefficient of x i on y, and r xiyj and P yxi denote the indirect path coefficient of x i on y through x i . The residual path coefficient P ye is denoted as: But the traditional steps of the residual path analysis are tedious, and after the research conducted by related scholars, the method of implementing the throughput analysis in SPSS 25.0 software (International Business Machines Corporation, Norman, Chicago, IL, USA) is proposed [49,50], and the formula for calculating the residual path coefficient is shown as follows [50]: where, ρ e represents the remaining path coefficient, and r represents the correlation coefficient. The formula represents that the fitting model is better and the selected independent variables have a greater effect on the dependent variables while ρ e is less than 0.1; however, the formula represents that there are still other factors influencing the dependent variable, and the quality of the fitting model is poor while ρ e is greater than 0.2.

Geographic Detector
Geographic detectors are statistical tools for detecting spatial variability and revealing the mechanisms underlying it. Factor detection, risk detection, interaction detection, and ecological detection are the most common geographic detectors [51]. The dispersion variances between data in the entire study region and the strata of variables are compared in factor detection. When compared to standard statistical analysis, the interaction detector can better identify the interaction between two explanatory variables [52,53]. The geographical driving forces of landscape pattern changes in different times were evaluated in this study with factor identification and interaction detection models, based on prior studies, from socioeconomic, physical, and accessibility variables ( Table 1). The following are the factor and interaction detection formulas: where, h and L are stratifications of variable Y or factor X; the N h and N are the number of units in the h stratum and the whole area, respectively; σ 2 h and σ 2 are the variance of the Y values of the h stratification and the whole region, respectively; SSW is the sum of the in-stratification variance, SST the region total variance. The value range of q is [0, 1]. The stronger the explanatory power of the independent variable X to the attribute Y, and vice versa, the higher the q value.

Spatial-Temporal Variation of Landscape Types from 2000 to 2020
Cultivated land, forestland, and construction land were the most widespread landscape types, accounting for around 90% of the study area's total area. As a result of the effects of increasing urbanization and environmental restoration efforts, the landscape types in CUAG changed dramatically from 2000 to 2020. These findings show that the area of construction land and forestland significantly increased by 224.26 km 2 and 142.17 km 2 respectively, whereas the area of cultivated land continuously decreased by 327.73 km 2 during the study periods ( Figure 2). The area transfer matrix of landscape types in different periods shown the direction and degree of changes among different landscape types (Figure 3). This indicated that the transfer of landscape types in CUAG mainly concentrated on forestland, construction land, and cultivated land with the trend of the "transfer-out" for cultivated land and the "transfer-in" for forestland and construction land. However, the degree of transfer varied in different periods. From 2000 to 2015, the cultivated land was mainly transferred to construction land with an area of 296.13 km 2 , accounting for 86.4% of the total "transfer-out" area of cultivated land in the past 20 years. Since 2015, the rate of cultivated land transfer to construction land was slowed, whereas the area of cultivated land transfer to forestland increased significantly with the area of 196.60 km 2 . The area transfer matrix of landscape types in different periods shown the direction and degree of changes among different landscape types (Figure 3). This indicated that the transfer of landscape types in CUAG mainly concentrated on forestland, construction land, and cultivated land with the trend of the "transfer-out" for cultivated land and the "transfer-in" for forestland and construction land. However, the degree of transfer varied in different periods. From 2000 to 2015, the cultivated land was mainly transferred to construction land with an area of 296.13 km 2 , accounting for 86.4% of the total "transferout" area of cultivated land in the past 20 years. Since 2015, the rate of cultivated land transfer to construction land was slowed, whereas the area of cultivated land transfer to forestland increased significantly with the area of 196.60 km 2 .
land, and cultivated land with the trend of the "transfer-out" for cultivated land and the "transfer-in" for forestland and construction land. However, the degree of transfer varied in different periods. From 2000 to 2015, the cultivated land was mainly transferred to construction land with an area of 296.13 km 2 , accounting for 86.4% of the total "transfer-out" area of cultivated land in the past 20 years. Since 2015, the rate of cultivated land transfer to construction land was slowed, whereas the area of cultivated land transfer to forestland increased significantly with the area of 196.60 km 2 .   Figure 4 shows the changes in landscape metrics for three dominant landscape types at class level from 2000 to 2020. It was found that the LPI and FRAC_AM indices of construction land increased steadily, whereas the NP and PD indices of construction land decreased notably; however, the opposite patterns were displayed for cultivated land. This indicated that construction land patches tended to concentrated distribution with the patch fragmentation reduced, but the patch shape was more complex. On the contrary, due to the dramatic decrement of cultivated land during the study periods, the patch distribution of cultivated land was more fragmented and dispersed with its patch shape regularization. In addition, the changes of LPI, NP, and PD indices for forestland were similar to that for construction land; however, the variation of these indices was lower than that of construction land; the FRAC_AM index for forestland was fluctuated and remained relatively stable level from 2000 to 2020; this demonstrated that the degree of fragmentation and heterogeneity for forestland patches also decreased.

Landscape Pattern Changes from 2000 to 2020
At the landscape level, the landscape pattern in CUAG also fluctuated significantly ( Figure 5). The values of AI and CONTAG indices dropped from 86.03 and 47.65 in 2000 to 85.74 and 46.46 in 2005, respectively; and then, these two indices increased continuously from 2005 to 2020. These findings indicated that the landscape lumpiness, aggregation, and connectivity in CUAG were enhanced with the trend of a concentrated distribution during the study period ( Figure 5a). However, the values of SHDI and SHEI indices increased significantly from 2000 to 2010, which was the period of frequent changes of landscape types in CUAG, changes that led to an upward trend of landscape heterogeneity and diversity level. From 2010 to 2020, with the rapid urbanization and the implementation of ecological restoration measures, construction land and forestland gradually became the dominant landscape types in CUAG, and the values of SHDI and SHEI indices were declined with the trend of landscape diversity decrement and landscape homogenization enhancement (Figure 5b). Table 3 shows the fitting models of temporal driving factors for three dominated landscape types (cultivated land, forestland, and construction land) changes in CUAG. The model-fitting results were reasonable with the determination coefficient of these models greater than 0.94 (p < 0.05). These indicated that industrial structure was the dominant temporal driving factor of the three dominated landscape type changes in the study area. Ratio of tertiary industry structure (x 11 ) had a negative driving effect for cultivated land changes and had a positive driving effect on the changes in forestland and construction land. In addition, population density (x 7 ) and the ratio of primary industrial structure (x 9 ) had a negative driving effect on cultivated land and construction land, respectively. patch fragmentation reduced, but the patch shape was more complex. On the contrary due to the dramatic decrement of cultivated land during the study periods, the patch dis tribution of cultivated land was more fragmented and dispersed with its patch shape reg ularization. In addition, the changes of LPI, NP, and PD indices for forestland were simila to that for construction land; however, the variation of these indices was lower than tha of construction land; the FRAC_AM index for forestland was fluctuated and remained relatively stable level from 2000 to 2020; this demonstrated that the degree of fragmenta tion and heterogeneity for forestland patches also decreased.  (Figure 5a). However, the values of SHDI and SHEI indices in creased significantly from 2000 to 2010, which was the period of frequent changes of land scape types in CUAG, changes that led to an upward trend of landscape heterogeneity and diversity level. From 2010 to 2020, with the rapid urbanization and the implementa tion of ecological restoration measures, construction land and forestland gradually be came the dominant landscape types in CUAG, and the values of SHDI and SHEI indice were declined with the trend of landscape diversity decrement and landscape homogeni zation enhancement (Figure 5b).  Table 3 shows the fitting models of temporal driving factors for three dominated landscape types (cultivated land, forestland, and construction land) changes in CUAG. The model-fitting results were reasonable with the determination coefficient of these models greater than 0.94 (p < 0.05). These indicated that industrial structure was the dominant temporal driving factor of the three dominated landscape type changes in the study area. Ratio of tertiary industry structure (x11) had a negative driving effect for cultivated land changes and had a positive driving effect on the changes in forestland and construction land. In addition, population density (x7) and the ratio of primary industrial structure (x9) had a negative driving effect on cultivated land and construction land, respectively.   Figure 6 displayed the attribution analysis of different factors for the spatial evolution of three dominated landscape types in different periods by using the factor detection model. The dominant spatial driving factors of each landscape type in different periods were varied. Slope (y 3 ), topographical relief (y 5 ), elevation (y 4 ), and distance to railway (y 7 ) had significant effects on the spatial changes of construction land (p < 0.05), whereas GDP (y 1 ), slope (y 3 ), topographical relief (y 5 ), and distance to railway (y 6 ) had significant effects on the spatial changes of cultivated land and forest land (p < 0.05). On the whole, the explanatory power of physical factors (i.e., topographical relief (y 5 ), elevation (y 4 ), and slope (y 3 )) on spatial change of landscape pattern were stronger than that of socioeconomic factors and accessibility factors, with the q values of physical factors 1. represent GDP, population density, slope, altitude, topographical relief, distance to railway, distance to expressway, distance to water system in turn. The green rectangle indicates that the driving result is significant at the 95% confidence level; yellow rectangular indicates that the driving result is not significant.

Temporal Driving Forces
The interaction effects of various factors presented a synergistic relationship between nonlinear enhancement or two-factor enhancement (Figure 7). It was found that the interactive two factors could enhance the explanatory power of the factors affecting landscape pattern changes. The interaction of topographic relief (y5) with distance from railway (y7), elevation (y4), and slope (y3) had a strong explanatory power to the spatial change of the three dominant landscape types during the study period. This demonstrated that the coupling effect of natural condition and spatial accessibility were the dominant spatial driving factors of landscape pattern changes in the study area. In addition, the interaction of GDP (y1) or population density (y2) with natural factors (i.e., slope (y3), elevation (y4), and/or topographical relief (y5) had a relatively higher explanatory power in some periods. detection, (c) cultivated land factors detection: y 1~y8 represent GDP, population density, slope, altitude, topographical relief, distance to railway, distance to expressway, distance to water system in turn. The green rectangle indicates that the driving result is significant at the 95% confidence level; yellow rectangular indicates that the driving result is not significant.
The interaction effects of various factors presented a synergistic relationship between nonlinear enhancement or two-factor enhancement ( Figure 7). It was found that the interactive two factors could enhance the explanatory power of the factors affecting landscape pattern changes. The interaction of topographic relief (y 5 ) with distance from railway (y 7 ), elevation (y 4 ), and slope (y 3 ) had a strong explanatory power to the spatial change of the three dominant landscape types during the study period. This demonstrated that the coupling effect of natural condition and spatial accessibility were the dominant spatial driving factors of landscape pattern changes in the study area. In addition, the interaction of GDP (y 1 ) or population density (y 2 ) with natural factors (i.e., slope (y 3 ), elevation (y 4 ), and/or topographical relief (y 5 ) had a relatively higher explanatory power in some periods. Interaction detection between the factors leading to major landscape types changes from 2000 to 2020: y1~y8 represent GDP, population density, slope, altitude, topographical relief, distance to railway, distance to high speed, distance to water system in turn.

Landscape Pattern Dynamics of Karst Mountainous Cities
Research has previously shown that urbanization can lead to changes in landscape types and landscape patterns [54,55]. The ecological sensitivity of karst region in China is high, and the landscape pattern changes are particularly dramatic [56,57]. Urbanization expansion has been the dominant trend and direction of land use structure change in recent years as a result of economic growth, rising population, and social development, which has had a significant impact on the landscape pattern [58][59][60].
In this study, the landscape pattern in CUAG saw a series of changes from 2000 to 2020 as a result of growing urbanization, with a tendency of clearly increasing building land and dramatically decreasing cultivated land, according to this study. This finding is in line with findings from other studies [36,60]. This illustrates that the main kind of landscape pattern change during urbanization is the conversion of cultivated land to construction land. In addition, the encroachment of ecological land (such as forestland and grassland) during the expansion of construction land is also one of the characteristics of landscape change in the process of urbanization (Figure 3). Wu [9] and Weng et al. [13] found that, with urban expansion, landscape pattern changes showed that ecological land (i.e., Figure 7. Interaction detection between the factors leading to major landscape types changes from 2000 to 2020: y 1~y8 represent GDP, population density, slope, altitude, topographical relief, distance to railway, distance to high speed, distance to water system in turn.

Landscape Pattern Dynamics of Karst Mountainous Cities
Research has previously shown that urbanization can lead to changes in landscape types and landscape patterns [54,55]. The ecological sensitivity of karst region in China is high, and the landscape pattern changes are particularly dramatic [56,57]. Urbanization expansion has been the dominant trend and direction of land use structure change in recent years as a result of economic growth, rising population, and social development, which has had a significant impact on the landscape pattern [58][59][60].
In this study, the landscape pattern in CUAG saw a series of changes from 2000 to 2020 as a result of growing urbanization, with a tendency of clearly increasing building land and dramatically decreasing cultivated land, according to this study. This finding is in line with findings from other studies [36,60]. This illustrates that the main kind of landscape pattern change during urbanization is the conversion of cultivated land to construction land. In addition, the encroachment of ecological land (such as forestland and grassland) during the expansion of construction land is also one of the characteristics of landscape change in the process of urbanization (Figure 3). Wu [9] and Weng et al. [13] found that, with urban expansion, landscape pattern changes showed that ecological land (i.e., forestland, grassland, water body) was transformed into construction land. Furthermore, forestland was the dominant landscape type in CUAG, and the area of forestland increased steadily during the study periods (Figures 2 and 4). This finding is consistent with the results of the study by Wang et al. [56] in Guiyang City. This is mainly benefited by the unique landscape matrix of karst mountainous cities. As a typical karst mountainous city, CUAG is characterized by large mountains stretching in patches and small hills sporadically distributed, and the land for construction is developed near the mountains, forming a unique urban-mountain mosaic landscape [54,56]. As a result, rich natural and near-natural forests depend on the mountains. On the other hand, the national and provincial ecological protection policies are also an important reason for the increasing area of forested land. For example, the Chinese government launched ecological restoration projects such as the Natural Forest Protection Program in 1998 and the Reforestation Program in 1999 [57], and the Guiyang city government has planned and built a mountain-river-wetland ecological network system mainly in urban mountains since 2000; and by 2020, the forest coverage rate of Guiyang city will reach 55% and the per capita park green space in the built-up area will be 13.2 km 2 [58].
Landscape patterns in rapidly urbanizing areas have presented a remarkable, highly fragmented feature. The single, continuous natural patches have become complex, heterogeneous, and discontinuous mosaics [60][61][62]. In this study, along with the construction land in CUAG expanded extensively, its patch displayed a trend of concentrated distribution with the low fragmentation (Figures 2 and 4). This finding is similar to the results of studies in the Mianzhu [63] and Saudi Arabian cities [64]. However, the patch shape in CUAG became more complex (Figures 2 and 4). Unlike plains cities, construction land continuously shifted outward to form a central area. The number of patches and shape complexity were declining, indicating that construction land was gradually converging in number and regular in shape; on the one hand, this may be related to the expressways construction and relocation of industrial areas [64]; on the other hand, Guiyang City had gradually strengthened its urban spatial form control measures. Nevertheless, the expansion of construction land was limited by topographic factors such as elevation and slope, and the shape of construction land patches was more complex in mountainous areas [65,66]. However, the growth of construction land in the study area experienced two stages rapid expansion from 2000 to 2015 and slow expansion from 2015 to 2020. It might be influenced by urban planning and land use policy. During 2000-2015, the urbanization pattern of multi-group development in Guiyang City has driven the intense expansion of construction land. The urbanization rate rose from 25% to 52% [67]. The disorderly spreading of construction land posed a serious threat to sustainable urban development, intensified ecological fragmentation, and overloaded resource and environmental capacity. In order to promote compact urban patterns and control urban scale, local governments have implemented policies to strictly limit the growth of construction land. In particular, after 2015, the Guiyang city government proposed layout optimization land remediation measures, and these policies led to a slowdown in the growth rate of construction land in the study area from 2015 to 2020.
During 2010-2020, landscape tended to be monolithic, with a more concentrated distribution of the phenomenon, showing the phenomenon of landscape homogenization, which is consistent with the findings of previous studies [26,52,68]. With regard to the impact of land-use change on the structure of the landscape as construction land and woodlands develop into more complex landscapes and as their spatial organization become more disorganized, urban planning may be used to rationally optimize the ecology of the surrounding region. Second, we see that homogeneity is a crucial process that goes hand in hand with the spatial urbanization process. The landscape had gradually closed due to the growth of the forest, and is now less functionally designed and more simply laid out. Landscape uniformity was also being fueled at the same time by the abandonment of farming, urbanization, and infrastructural developments [26]. Therefore, how to fairly distribute landscape types and spatial distribution in the course of urbanization may be a trend for future land use landscape planning and environmental decision-making.

Driving Force of Landscape Pattern Dynamics
The research of spatial-temporal drivers of urban landscape change is the foundation for revealing the inner mechanism of landscape pattern-process-function, and also guarantees the optimization of landscape space and the construction of landscape ecological security pattern [31]. In this study, the ratio of tertiary industry structure and the ratio of primary industrial structure were the significantly exploratory temporal factors for landscape pattern changes in CUAG. These findings indicate that industrial structure was a key determinant of landscape pattern. These findings are consistent with previous related studies elsewhere that have shown the upgrading and adjustment of industrial structure as driver of urban expansion [48,69]. It is argued that economic structural transformation could accelerate urban expansion, since this process increased land demand associated with the growing private and foreign enterprises [69,70] and cities have competed intensively to attract economic development through improving infrastructure and setting up development zones [71,72]. Moreover, the demographic variables, total population, and population density are recognized to be positive independent factors of landscape pattern change [71,[73][74][75]. This result also agrees with our study that population density is an essential determinant of landscape pattern. It makes intuitive sense that more space will be occupied as more people inhabit [71,76]. An increase in population density will inevitably give rise to higher demand on residence, public facilities, and transportation, which would accelerate landscape pattern change [71,77].
Among the spatial variables, natural factors and accessibility factors are frequently found to be exploratory factors of landscape pattern change during urban expansion [26,77]. Areas with low-lying terrain and flat topography are conducive to the formation of construction land [78][79][80]. It is argued that human activities will encounter fewer obstacles in areas with such topological characteristics [75]. In this study, topographic relief was considered the special factor in the spatial evolution of urban landscape patterns. Areas with large topographic relief have less human interference, and forests have significant advantages under this terrain [77,81]. Our results support such arguments. In addition, as for the accessibility factor, when the expansion of rapid urbanization slows down and the urban landscape develops in an aggregated manner, with the reorganization of the landscape as the main way, socio-economic factors and urban spatial accessibility gradually become the main factors affecting the spatial evolution of the landscape pattern, which corroborates the study of Sebastian [77]. In addition, high spatial accessibility brings about high land use intensity and expands urban space, the connection between districts is closer, towns are developed in succession. The construction of roads, highways, and railways organically connects distant towns. This promotes urbanization and changes the urban landscape pattern [71].
In this study, small cities, slow development, and a wealth of natural resources put the CUAG's geographic location in terms of economics and transportation at a major disadvantage. The local ecological and cultural environment will suffer from overdevelopment and city-scale expansion. In order to create a city with ecological functions while combining its unique geographic endowment, it is vital to strike a balance between ecological environmental conservation and urbanization expansion.

Conclusions
We conclude the following: (1) the landscape structure in CUAG from 2000 to 2020 was drastically changed. The area of forestland and construction land continuously increased and became the dominant landscape types since 2015. The decrease of cultivated land is the main contribution of landscape structure change to a trend of rapid urbanization and ecologization. (2) The patch aggregation of construction land and forestland tended to improve, with the degree of fragmentation decreased and the patch shape complicated. Meanwhile, the patch fragmentation of cultivated land was intensified. During the study periods, urban landscape connectivity was improved; however, the trend of landscape homogenization was obvious with the level of landscape diversity declined. (3) Industrial structure was the dominant temporal influencing factor for major landscape types changes, topographic relief and the interaction factors between topographic relief and regional spatial accessibility or social economy were the main spatial influencing factors for landscape pattern changes. The findings of this study on the temporal and spatial changes of landscape patterns and their driving factors in typical karst mountainous cities could provide a reference for landscape pattern planning and the research of the relationship between landscape patterns and processes under the background of rapid urbanization and ecologization in mountainous regions.