The UANSTM is located in the center of Eurasia (Figure 1). It expands from 83°24′ E to 91°54′ E in an east–west direction and 41°11′ N to 46°12′ N from south to north. The agglomeration is bordered by the Gurbantünggüt Desert in the north, the Tianshan Mountains in the south, Jinghe County of Bortala Mongolian Autonomous Prefecture in the west, and the Balikun Kazakh Autonomous County in Hami City in the east. The average elevation is approximately 1,000 m. The total area of the urban agglomeration is 215,400 km², accounting for 13% of Xinjiang’s total area, and the municipal administrative units include Urumqi, Wujiaqu, Changji Hui People Autonomous Prefecture, Turpan, Shihezi, Karamay, Kuitun, Wusu, and Shawan City,.

The UANSTM has the most developed economy and transportation system, the densest population, and the most concentrated industry in Xinjiang. It plays an essential role in constructing the core area of the Silk Economic Belt in the new era, is crucial to the national security, and serves as a key facilitator in achieving new urbanization in Xinjiang, a task set in the 14th Five-Year Plan and the 2035 Long-Term Goal. The climate is temperate continental climate, with large diurnal temperature variations and a large annual variation. The ecological environment in this region is fragile, and under the combined impacts of human activities and climate, it faces shortages of resources and land. As a result of salinization, glacier retreat, and the ecological fragility of this region, environmental degradation is serious. The fragile ecological environment largely restricts the sustainable development of the social economy of the UANSTM.

This study used geographical, climatological, and soil data (Table 1).

1) Geographical data included land use data from the Data Center for Resources and Environmental Sciences with a resolution of 30 m (https://www.resdc.cn/), digital elevation model (DEM) data from the Geospatial Data Cloud with a resolution of 30 m, from which slope and relief amplitude were extracted (http://www.gscloud.cn/), normalized difference vegetation index (NDVI) data from the China National Space Administration (http://www.cnsa.gov.cn/), and river basin data from the National Geomatics Center of China (http://www.ngcc.cn/ngcc/).

The UANSTM is a complex agglomeration that includes urban and natural landscapes. In recent years, due to population growth and the continuous expansion of construction land, various ecological and environmental problems, such as declining air quality and water pollution, have surfaced. Therefore, for ecological risk assessment, five factors were selected as the ecological risk factors in this region. According to the ecological risk characterization paradigm “probability and loss multiplication” thinking, the corresponding ecological risk probability is defined as the probability of the occurrence of factors that would result in the reduction of the ecosystem service values in a certain period.

We used Arc GIS 10.7 to classify the aforementioned risk probability indicators into five levels: 1, 2, 3, 4, and 5. The higher the level, the higher the probability of ecological risks; the combined maximum method was used to obtain a composite evaluation of the probability. That is, if any indicator in the grid has a higher level, the probability of risks occurring in that grid is high. The discrete coefficient method was used to construct a composite index of the aforementioned three ecosystem services, and the following formula was used to characterize ecological loss (Zhang et al., 2022). Finally, based on the comprehensive assessment results of probability and loss, a two-dimensional matrix of ecological risk assessment was constructed. Calculation was performed as follows: E S i = σ i j x i j = 1 x i j 1 N ∑ j = 1 N x i j − x ¯ i j 2 E S c o m = E S E S 1 , E S 2 , E S 3 , … E S n , (12) where E S i is the ecosystem service indicator of the i -th grid; x i j is the normalized value for the j -th ecosystem service, which is a normalized average; E S c o m is the composite indicator for ecosystem services; and Δ E S i indicates the degeneration of ecosystem services.

The bivariate spatial autocorrelation model is a type of exploratory spatial data analysis (ESDA). Generally, ESDA can reveal the spatial pattern and mechanism of a geographical object by describing its spatial distribution, visualizing it, and identifying its spatial layout, agglomeration characteristics, and anomalies (Bi et al., 2007; Acevedo Bohorquez and Velasquez Ceballos, 2008). The ESDA includes global and local spatial autocorrelation. While the first is used to study the overall spatial relationship, the latter is applied to investigate spatial variation (Kang et al., 2016). The ESDA was introduced through ArcGIS and the GeoDa spatial econometrics analysis software, and we constructed a bivariate exploration mechanism between the ecological risk index and ecosystem service value.

Global spatial autocorrelation is usually expressed in terms of Moran’s I exponent. The closer the absolute value of Moran’s I to 1, the stronger the spatial autocorrelation of the study unit. The local indicators of spatial association (LISA) map was used to identify the spatial association patterns of high and low values in ecological risk areas, using the following equations (Wu et al., 2021a): I = n ∑ i = 1 n ∑ j = 1 n w i j x i − x ¯ x j − x ¯ ∑ i = 1 n ∑ j = 1 n w i j ∑ i = 1 n x i − x ¯ 2 , (13) I i = n x i − x ¯ ∑ j = 1 n w i j x j − x ¯ ∑ i = 1 n x i − x ¯ 2 , (14) where I is the ecological risk indicator, i is the ecosystem service value, and j is the bivariant autocorrelation coefficient; x_i and x_j are the ecological risk indicator and ecosystem service values of a specific unit, respectively, n is the number of units in the study, and ω_ij is the neighboring weight.

Figure 3 and Figure 4 show the spatial distribution of the ecological risk probability index. The slope is steep in the northwestern and mid-eastern parts and moderate in the south and north. Steep-slope areas are mainly distributed in the high-elevation areas of the UANSTM, such as Wusu City, Shawan County, Manas County, Hutubi County, Changji City, the southern part of Urumqi County, the junction of Fukang City, Jimsar County, Qitai County, and Mulei Kazakh Autonomous County in the south, Dabancheng District, Gaochang District, and the northern part of Shanshan County. The northern and southern areas typically have lower slope values. High VFC values can mainly be found in the northwest and center parts as woodlands and grasslands are largely distributed in this area. High DPC values are mainly concentrated in the northwest and central parts of the UANSTM, whereas low values can be found in the north and southwest because the unused land in these parts of the region impedes ecological processes. Considerable wind erosion is mainly observed in the northeast and south, and wind erosion is particularly strong in Mulei Kazakh Autonomous County because there are more unused sites and less grasslands. Water erosion is high in the northwestern and central parts and low in the north and southeast; Shanshan County is most prone to water erosion, mainly because of its steep terrain and abundant precipitation. The frozen-thaw erosion areas in the UANSTM are relatively small and mostly consistent with the distribution range of water erosion areas, mainly because of the high elevation, low temperatures, and high radiation.

Habitat quality showed a spatial decline pattern from the northwest–southeast line to the northern and southern sides. This is because the mild climate, abundant water resources, and rich ecosystems in the northwest-southeast area provide favorable conditions for numerous species. From 2010 to 2020, the spatial distribution pattern of habitat quality in the UANSTM changed significantly; habitat quality decreased significantly in the northwest and improved significantly in the central and eastern areas, and the overall habitat quality increased by 0.13%.

The high water production in the northwest and the low production in the north and south are mainly due to the high VFC and abundant water resources in the northwestern areas. From 2010 to 2020, the water yield of the UANSTM decreased considerably in the southwest and increased slightly in the central and eastern regions. Specifically, the high-value areas were mainly concentrated in the southern mountainous areas such as Wujiu City, Shawan County, and Urumqi, as well as the mountainous woodland areas with generally high WY values. The low-value area was located in the southeast of the urban agglomeration (especially in Shanshan County), and the total water production decreased by 0.34%. Overall, the dry climate, limited rainfall, and high proportion of undeveloped and unused land resulted in a low WY in the UANSTM.

Soil conservation also decreased from the northwest–southeast line to the north and south, mainly because the predominant land use type in this area is grassland, which can intercept rainfall; consequently, more rain can infiltrate the soil, maintaining soil nutrients. From 2010 to 2020, the spatial distribution pattern of soil conservation in the UANSTM changed significantly, among which the soil conservation values in the southwestern, central, and eastern regions increased significantly, and the overall soil conservation increased by 0.12%. High soil conservation values were mainly observed in Urumqi County, Wusu, Shawan, and Changji City, whereas low values were found in the southeastern part of the UANSTM, mainly in Shanshan County, Sen County, and Toksun County.

Overall, the ecosystem services in the UANSTM showed a decreasing trend from the northwest-middle east line to the south and north sides, and the distribution characteristics were consistent with those of the ecological risk probability index (Figure 5A and B). Because the land use type in the northern and southern parts of the study area is mostly unused land, with sparse precipitation, dense evaporation, and large-scale desertification, the total values of ecosystem services in the southern and northern parts of the region were low. In addition, the central and western regions are rich in vegetation, with a wide spatial distribution of grassland and woodland. For instance, the southern part of Changji Hui People Autonomous Prefecture and the northwest of Turpan Prefecture are densely forested, and the total values of ecosystem services are high. From 2010 to 2020, the total ecosystem services in the UANSTM decreased by 0.43%. Spatially, ecosystem services showed drastic changes in the central region but varied less in the north and south (Figure 5C). The areas with the largest declines were Karamay City, central Changji Hui People Autonomous Prefecture, and central Urumqi. The largest increases were observed in the southern part of Changji Hui People Autonomous Prefecture and the southwestern part of Urumqi City.

The distribution characteristics of ecological risks in the UANSTM are consistent with those of ecosystem services. Based on the natural breakpoint classification method, six risk probability indicators are reclassified to form five levels. The higher the grade, the higher the probability of ecological risk. A comprehensive evaluation result of ecological risk probability was obtained by using the maximum synthesis method (Figure 6). The results of probability integration showed that the ecological risk of the UANSTM increased from the northwest to the middle east, and the distribution characteristics were consistent with those of the ecological risk probability index. The low- and medium-risk areas accounted for a relatively small proportion, whereas the high-risk areas accounted for a relatively large proportion.

Low and medium risks were concentrated in the northwest, mainly including the northern part of Changji Hui Autonomous Prefecture, the southern part of Karamay City, Yili Kazakh Autonomous Prefecture, and the northern part of Urumqi. In these areas, woodland and grassland are the main land use types, the ecosystem has good integrity and stability, the total ecosystem services are numerous, and the losses are small; consequently, the risk occurrence probability is relatively low. High-risk areas were located in the northern and southern parts of the study area. The VFC in this area was low as the predominant land use type was unused land, with poor landscape connectivity, high wind erosion intensity, few ecosystem services, and high losses, with a high ecological risk. Furthermore, slopes, freeze-thaw erosion, and water erosion in some parts of the central and northwestern parts were considerable, which also resulted in high ecological risks.

Moran’s I is frequently used to investigate the overall distribution and spatial aggregation of a region, but it cannot indicate the spatial correlation. Therefore, the LISA analysis was used to explore the correlation degree of ecological risks in the study area and determine whether there was spatial aggregation. The LISA cluster diagram of ecological risks and ecosystem service loss was used to show local spatial autocorrelation, and the two factors were positively correlated (Figure 7). Four significant autocorrelations emerged in the study area, namely, high-high (HH), low-low (LL), low-high (LH), and high-low (HL). The HH areas were mainly concentrated in the north and south, LL areas in the northwestern, central, and eastern parts, and HL and LH areas were scattered throughout the study region. Areas with higher ecological risks correspond to areas with greater losses of ecosystem services at the spatial scale.

In this study, the ecological risk was negatively correlated with the total number of ecosystem services and positively correlated with the loss of ecosystem services. Since the 1970s, research on ecological risk assessment has been developing gradually (Munns, 2006). The main purpose of ecological risk assessment is to provide decision support for ecosystem protection and management (Wong et al., 2015). However, traditional ecological risk assessment methods are easily affected by the complexity and variability of ecosystems, and it is difficult to take into account the multi-objects that need to be protected in ecosystems (Kang et al., 2016). When only ecological entities are protected, without considering human wellbeing, it is difficult to implement the specific decision-making level. With the continuous improvement of human living standards, ecological risk assessment is gradually becoming an important measure of ecosystem services (Rapport, 2007; Peng et al., 2015). Some authors (Faber and van Wensem, 2012) believe that ecological risk assessment will eventually move toward the relationship between ecosystem processes and services, and the assessment of ecosystem services will serve as the end point of ecological risk assessment, connecting ecological processes with ecological risk sources. An ecological risk assessment framework oriented to ecosystem services can reduce the error between the loss of ecosystem services and the management of complex ecosystems (Galic et al., 2012; Xu et al., 2014).

However, there are three challenges in integrating ecosystem services into regional ecological risk assessment. First, the existing assessment methods tend to focus on a single ecosystem service rather than considering multiple ecosystem services as a whole system (Qian et al., 2022). Second, the analytical framework needs to be further improved to implement comprehensive assessment, focusing not only on the impact on ecosystem services but also on the status of ecosystems (Xu et al., 2014; Wong et al., 2015; Xu et al., 2016). Some researchers suggested that the goal of protection should not only be the provision of diversified ecosystem services but also the maintenance of a high level of ecosystem health (Rapport et al., 1999). A healthy ecosystem is considered as the development goal of environmental management, emphasizing the integrity of the ecosystem and providing the basis for ecosystem assessment (Costanza, 2012). However, until now, few studies have established an ecological risk assessment framework oriented by the combination of ecosystem services and ecosystem health. Nevertheless, compensating for the deficiency of traditional ecological risk assessment methods and improving the efficiency of risk management, the introduction of ecosystem services into ecological risk assessment for development and improvement has become the frontier and hotspot of current ecological risk assessment research (Kang et al., 2016).

This study, using the UANSTM as the study area, evaluated ecological risks by selecting factors based on four aspects: terrain sensitivity, ecological resilience, landscape vulnerability, and ecological sensitivity. It used the InVEST model to quantify the three key ecosystem services, namely, water production, habitat quality, and soil conservation to evaluate the ES pattern. Subsequently, the correlations between ecological risks and ecosystem services were analyzed.

Regarding the methodology, in this paper, four ecological risk probability indicators were calculated, which is consistent with the ecological risk assessment method of Huang et al, (2022) but differs from that used by Qian et al, (2022), which was based on the landscape ecological loss index (LELI) model. Zhang et al. (2022) used two indicators, landscape disturbance and landscape fragility, and applied the landscape pattern index model to accurately evaluate the LER caused by different landscape patterns. The improved ERA model adopted in this paper facilitates the spatial visualization of multi-source ecological risks, provides a scientific basis for strengthening ecosystem management, and adds a response analysis of ecological risks and ecosystem services, making ecological risk assessment more ecologically indicative. In terms of the results, this paper used the data from 2010 to 2020 to evaluate and analyze the ecological risk and ecosystem service pattern of the UANSTM, concluding that ecological risk and ecosystem service loss are positively correlated spatially. However, Zhang et al, (2022) took the Xi’an metropolitan area as the research area and used remote sensing data from 2010 to 2020 to explore the relationship between landscape ecological risks and ecosystem services (such as carbon sequestration, soil conservation, water yield, crop yield, and habitat support). The authors concluded that there was no obvious linear relationship between landscape ecological risks and ecosystem services. The two studies produced different results, mostly because of two aspects. First, the Xi’an metropolitan area is relatively small, whereas the UANSTM covers a large area, and there are certain differences in the data on ecological risks and ecosystem services. Second, the methods used were also different; the present study calculated the ecological risk indicators, whereas Zhang et al. used the landscape pattern index to calculate the ecological risk. The landscape pattern index method only focuses on the fragmentation of static landscape patterns and the vulnerability of land use types, but it does not consider ecosystem processes and environmental conditions. Therefore, the evaluation results are not comprehensive. Aziz (2023) mapped the potential and realized ecosystem services intensity from 140 of Pakistan’s terrestrial protected areas (TPAs) using the Co$ting Nature model, with the potential services being positively and the TPA size being negatively correlated with the realized ecosystem services. The research methodology using in this paper is similar to that applied here. The present study provides a suitable method to investigate the relationship between ecological risks and ecosystem services and represents a scientific reference for subsequent related research.

Our study has some limitations. First, the spatial resolution of the data used needs to be further improved to achieve a more accurate assessment. The meteorological and soil data had a resolution of 1,000 m, which is relatively low. The research area was a relatively small local area, and the low resolution may cause microscopic deviations when making relevant calculations. Second, this paper does not consider the trade-offs and coordination relationships among ecosystem services, which may have a certain impact on the composite index of ecosystem services, which in turn affects ecological loss and risk assessment results. Third, we selected factors based on the four aspects, namely, ecological resilience, landscape vulnerability, terrain sensitivity, and ecological sensitivity, to evaluate the ecological risk, and the selection of factors was subject to certain subjectivity. With this in mind, the establishment of the index system can still be optimized.

According to our understanding and analysis of previous studies and the limitations of this study, the following aspects can be further improved: 1) when conducting small-area research, it is necessary to further improve the accuracy of the data, supplement relevant scientific data, and optimize the index system according to the data availability; 2) in follow-up research, long-term dynamics of key ecosystem services should be explored to make ecological risk assessment predictive; and 3) since the interaction between landscape pattern and ecological processes is highly uncertain, the ecological risk assessment results are not the only constant values, and the identification of areas of uncertainty in risk level should be strengthened to reduce the impact of subjective perception on the assessment results.