Quanjiao County (117°49′～118°25′N, 31°51′～32°15′E) is a suburban county of Chuzhou City, located in the southern part of Chuzhou, Anhui Province, China (Figure 1), with a total land area of 1,568 km². Quanjiao is a typically hilly suburb county of eastern China, and it has all kinds of natural elements such as hill, forest, grassland, farmland, wetland, and river. The total population of Quanjiao County was 452,000 by the end of 2019, with the urbanization rate of 56.9% (SOQC, 2020). It is also a part of the Yangtze River Delta (YRD) region (Figure 1), which is one of the regions with the highest urbanization levels in China. Quanjiao County is 48 kilometers to Nanjing with an 83.2% urbanization rate, 98 kilometers to Hefei with a 76.3% urbanization rate, and 360 kilometers to Shanghai with an 88.1% urbanization rate. Surrounded by these extensive urban sprawls, Quanjiao is faced with the pressure of being intensively urbanized, with LUCC and ecological corridor loss.

The following data sources were used in this study: 1) seven land cover types, including cultivated land, forest, grassland, water area, built-up land, rural settlement, and unused land, extracted from 2010, 2015, and 2018 land use maps obtained from the National Geomatics Center of China (NGCC) with a 30-m resolution. The Third National Land Survey of the study area was used to calibrate the land use map in 2018; 2) precipitation and rainfall erosivity data were extracted from Chinese National Average Annual Rainfall and the Chinese National Rainfall Erosivity Map (Xun and Zhang, 2021), which were downloaded from the Resource and Environment Science and Data Center, China (https://www.resdc.cn/); 3) soil data were obtained from The Second National Soil Census Data, China. Soil characteristics, for example, soil texture and organic matter content were used to determine the standard content of sand, silt, and clay and calculate soil erosion factor by the EPIC model formula; 4) a 20-m resolution digital elevation model (DEM) was downloaded from the Geospatial Data Cloud site, Computer Network Information Center, Chinese Academy of Sciences (http://www.gscloud.cn), from which were extracted or computed elevation, slope, and slope length.

A comprehensive framework for optimizing KENPRA by integrating ESP and LUC was proposed (Figure 2), and four main technical steps were identified. Multisource data collection and data set building were conducted before the main methodological steps. The first step was to identify and establish ESPs with three sub-steps. The first sub-step was to identify an ecological source through an evaluation system, including habitat quality, ecosystem services, and landscape connectivity. The second sub-step was to construct a resistance surface by assigning resistance coefficients to different land use types and geomorphological factors. Ecological corridors were generated with circuit theory by using the results of ecological source and resistance surface. The third sub-step was to develop ecological security patterns by determining main ecological sources and ecological corridors. Moreover, a multicriteria evaluation system was built to evaluate the intensity of LUC from different land use functions, for example, ecological land conflict, cultivated land conflict, and construction land conflict. Each land use function was assessed based on suitability and driving forces of land use conversion. The first sub-step of assessing habitat quality, ecosystem services, and landscape connectivity laid the foundation for evaluating ecological land suitability for this sub-step. The third step was to optimize KENPRA from the three perspectives of ecological security: importance, criticality, and timeliness. Key ecological nodes (KENs) and priority restoration areas (PRAs) were identified through the assessment results of ecological pinch, ecological barrier, and LUC index.

Land use conflict arises from the scarcity of resources with obvious spatial attribution (Zou et al., 2019). A hierarchical multicriteria evaluation system (Jiang et al., 2020) was developed for land use conflict (Table 2) from land function suitability and potential driving forces (Liu et al., 2014; Iojă et al., 2014). Cultivated LUC of study areas is mainly influenced by natural factors (Liu et al., 2014; Chen et al., 2021), for example, slope, surface soil texture, and the potential driving forces of land use conversions and mainly include social and location factors (Table 3). Developed LUC was assessed by suitability factors, including natural, social, and economic factors, and driving forces of land use conversion (Supplementary Appendix Table SB2). Ecological land suitability was evaluated by the results of ecosystem service provision and habitat quality. The LUC index can be calculated using the following Eq. 7 (Jiang et al., 2020): P i = ∑ W 1 × W 2 × W 3 × f i j , (7) where Pi refers to the LUC of grid I; W₁, W₂, and W₃ are the weights of indicators at different levels; and f _ij is the score of index j in grid i. We used the natural breaks method to divide the conflict into three scales: strong, medium, and weak. LUCs were divided into 27 categories (Jiang et al., 2020) and assigned corresponding to land use conflict intensities (Table 2). The LUC values ranged from 0 to 1, where 1 stands for the highest conflict and 0 indicates no conflict.

KENPRA identification is based on the importance, criticality, and time for ecological nodes and restoration areas to ESP (Fu et al., 2020; Wang et al., 2021). Ecological pinchpoint, ecological barrier, and LUC (as an index) were selected to indicate these three aspects (Zhang D. et al., 2020; Jiang et al., 2020). Currently, ecological pinchpoints are highly intensive areas, indicating that there are few or no alternative paths for species and other ecological processes, which are irreplaceable and important in landscape connectivity (Peng et al., 2018b). Ecological barriers, which have a high possibility of impeding species movement between patches of habitats, are critical areas for implementing restoration projects that could increase the likelihood of ecological corridor connectivity by restoring these nodes. If there are no timely protections or restoration projects involved, LUC areas will have a higher probability of having intense land use conflict and being urbanized.

Ecological pinchpoints and ecological barrier nodes and areas were identified by using the Pinchpoint Mapper Toolbox of ArcGIS and Circuitscape tool (http://www.circuitscape.org/linkagemapper). The ecological pinchpoints were generated using current intensity data (Peng et al., 2018b). Ecological barrier nodes were created by using the Barrier Mapper module of the Circuitscape tool, which found the value of important connectivity nodes by calculating the connectivity recovery value after the removal of barrier points (Mcrae et al., 2012). LUC node identification was based on the assessment of LUC. Each layer of ecological pinchpoint, ecological barrier, and LUC index was divided into four classes: 1) extremely high, 2) high, 3) medium, and 4) low. Priority restoration areas were extracted from each layer of extremely high areas and areas which were not less than 1 km². Key ecological nodes were optimized and extracted from the top 10% of the high-value areas with not less than 1 km² (Fang et al., 2020) and were transferred from polygon to point in ArcGIS 10.6.

To better understand the effects of KENPRA on future land use land cover (LULC) and conservation, a CA-Markov model was used to predict future LULC under different KENPRA-based scenarios. The CA-Markov model combines Cellular Automata and Markov Chain, which focuses on quantitatively assessing and predicting LULC (Ghosh and Chifos, 2017). In our study, a CA-Markov model was used to predict and assess the land use scenarios in 2030 based on the historical land use data. Two scenarios were developed based on ecological security source area identification under different schemes. In Scenario A, the protection of the ecological nodes and key areas of restoration, together with extremely high current density and LUC areas, should be restricted. In Scenario B, KENPRA identified and protected only extremely high current density areas.

The CA-Markov model predicted the LULC of future years based on the transition probability matrix from at least two past time periods. It was calibrated by using suitability maps with spatially explicit factors and constraint variables for each LULC type. We used the transition probability matrix created by LULC of 2015 and 2018 and suitability maps to predict the LULC in 2030 under Scenario A and Scenario B. LULC maps of 2015 and 2018 were prepared as raster layers with a cell size of 10 × 10 m. We converted LULC raster layers to rst form in IDRISI Selva, and the transition probability matrix was initially created by using a Markov model from LULC of 2015 and 2018 as the base years. A multicriteria evaluation module was applied to calibrate the model by creating suitability maps of all LULC types integrated with an initial transition probability matrix in CA-Markov to produce a future growth potential surface (Ghosh et al., 2021). Suitability maps were considered suitable factors and constraints for each LULC type. Constraint factors serve to limit the LULC under alternatives. Land policy and land condition slope, for example, permanent prime farmland have been considered to be constraints factors (Han et al., 2009; Ghosh and Chifos, 2017). All LULC and suitability layers were converted into 10 m*10 m grid cells according to their original resolution size. The spatial data preparation and raster analyses were implemented in ArcGIS v. 10.6 (ESRI, 2018), and the Markov and CA-Markov model were implemented in IDRISI Selva (Eastman, 2012).We used a transition probability matrix created in 2015 and 2018 LULC and suitability maps to predict the LULC until 2030 under Scenario A and Scenario B. The Kappa coefficient was used for verifying the prediction accuracy of the CA-Markov model (Pontius an d Millones, 2011), which was calculated by the CROSSTAB tool in IDRISI software. Models with a Kappa coefficient greater than 0.80 were considered achieving a high consistency (Feinstein and Cicchetti, 1990), which was found to be 0.861 and 0.836 under Scenario A and Scenario B, respectively, suggesting a good performance. K a p p a   c o e f f i c i e n t = ( T × C ) − G T 2 − G (8) Where, T is the test pixels, C is the correctly classified pixels observations, G is the sum of multiplied total value.

The results of habitat quality, ecosystem services, and landscape connectivity in 2018 showed spatial heterogeneity across the study area (Figure 3), while the extremely high areas of these three indictors were concentrated in the northwest comprising hilly areas. Specifically, areas with extremely high values of habitat quality and ecosystem services were mostly distributed in the northwest, accounting for 34.1 and 8.1% of the total area, respectively. The same as the result of landscape connectivity that extremely high and high areas were distributed among continuous forests in the northwest, which had a national forest park and a headwater reserve with high species movement connectivity capacity. Finally, 15 patches of ecological sources were selected in the study area (Figure 3), with a total area of 21,497.75 ha, which were mainly distributed in the northern hilly areas and accounted for 84.8% of the total area of ecological sources. To verify the identification result of ecological sources, we compared the selected ecological sources with the existing national and provincial nature reserves, the results of which showed good spatial consistency.

The ecological resistance surface modified by dominant resistance and recessive resistance is shown in Figure 4A. Overall, the built-up areas had the highest value of ecological resistance, and their value decreased gradually from downtown centers to other areas. The intensive cumulative current corridors mean high current density and value of connectivity as shown in Figure 4B. High and extremely high values of cumulative current constituted only 6.8 and 1.0% of the total length of the cumulative current surface, respectively, and were found only in the corridors of cultivated land in the southeast and corridors connecting cultivated land and forest in the middle of Quanjiao County. Finally, 31 ecological corridors were identified with a total length of 214.28 km (Figure 4C), with a maximum, minimum, and average length of 25.00, 0.73, and 8.93 km, respectively.

The assessment of LUC in 2018 showed that 82.56% of the study area had potential conflict, including extremely intense, intense, medium, and weak levels. Only 17.44% of the study area had no land use conflict (Figure 5A). Areas with extremely intense and intense LUC were mainly distributed in southeast areas (7,179.03 m² and 40,825.09 ha, respectively), which were dominated by cultivated lands. Apparently, areas of no land use conflict were distributed in the northern hilly areas, which were inconsistent with the spatial patterns of high habitat quality and multiecosystem services provision, indicating that these areas are highly suitable and critical for ecological protection.

LUC in ecological corridors is shown in Figure 5B. Areas with high values of LUC were mainly distributed in the southeast of Quanjiao County, which were close to the downtown and accounted for 5.1% of the total length of ecological corridors. These areas were mainly in ecological corridors EC-6, EC-8, and EC-18. Of special note, corridor EC-6 achieved approximately 50% of the total extremely intense areas of corridors.

A total of 34 key ecological nodes were found (Figure 6), including 17 key LUC nodes, 11 pinchpoints, and six barrier nodes. Key ecological nodes of LUC points were distributed mainly in the southeast and surrounded by downtown across intense LUC cultivated land. Surprisingly, we also found three nodes in the northern forest, adjacent to the north boundary, which may be due to the conflict between conservation and rural building development. By contrast, key ecological nodes of pinchpoints and barriers were mainly distributed in the north and middle of Quanjiao County, which played important roles in connecting ecological sources between cultivated land and forest.

In addition, a total area of 4,357.2 ha were identified as priority corridors (Figure 6), including 2,485.2 ha of intensive land conflict corridors and 1,872 ha of high probability of barrier corridors, comprising 28.2 and 21.3% of the total corridor areas, respectively. Dense priority restoration areas were located in southeast Quanjiao county, resulting in potentially intensive land use conflicts. Several priority restoration areas in the middle of Quanjiao were found due to its importance of connectivity for species migration and ecological process flow between southern cultivated land and northern hilly areas. The strategies for land management for each priority restoration area were proposed (Supplementary Appendix Table SB3).

The change in LULC showed that the built-up areas increased slightly from 0.6% in 2010 to 0.7% in 2015 and increased rapidly to 2.36% in 2018 (Supplementary Appendix Table SB4). In contrast, the overall vegetation covers of the study area decreased gradually by 2.6% from 2010 to 2018, indicating that the built-up expansion was at the cost of cultivated land, forest, and grassland. The changes in LULC showed a similar trend, with different changes in different LULC type under these two KENPRA-based simulation scenarios from 2018 to 2030 (Figure 7; Supplementary Appendix Table SB4). Specifically, built-up areas are projected to have a larger expansion under Scenario B (310.23 ha) than Scenario A (857.53 ha) during the period of 2018–2030. Moreover, 1.46 and 2.48 ha of ecological corridors are likely to be occupied by built-up land, while −6.19 and 17.7 ha of ecological corridors will expand to rural settlement under Scenario A and Scenario B, respectively. The area of forest showed a significant increase, accounting for 24.1%% and 22.31%% of land. In contrast, the area of cultivated land decreased to 63.81 and 64.46% under Scenario A and Scenario B, respectively. Moreover, the overall growth of forest, grassland, and water areas in ecological corridors will have a greater growth under Scenario A (3,479.05 ha) than under Scenario B (3,259.76 ha) during the 2018–2030 period. In brief, the simulation result revealed that ecological security patterns will achieve better conservation efficiency with less ESP destruction under Scenario A than Scenario B in the period between 2018 and 2030.

Our study demonstrated a way to develop ESPs through ecological sources and ecological corridors to facilitate KENPRA identification at Quanjiao County. The results suggested that ESP linking between ecological processes and landscape patterns could help us better understand ecological processes of movement and fluxes and improve ecological integrity by conserving and restoring core patches and ecological corridors, especially in suburban areas where there is great pressure of rapid LUCC and land expansion (Cook, 2002; Peng et al., 2018b)

We identified ecological sources from multiple aspects, including habitat quality, ecosystem services, and landscape connectivity, the results of which all showed good spatial consistency with existing high ecological value areas, for example, nature reserves. The result is consistent with Fu’s study (Fu et al., 2020), which found that assessing habitat maintenance, ecosystem service importance, and landscape connectivity showed a good performance in delimiting ecological patches in Loess Plateau, China. Compared with previous paradigms that addressed ecological source identification on representative species or target group habitat (Hepcana et al., 2009; Luque et al., 2012), our approach embraced multiple important ecological processes by assessing critical ecosystem services and internal landscape connectivity. Ecosystem services, which were translated from the structure and function of ecosystems into the provision of important services (Daily et al., 2009), could better predict critical ecosystem process and function by quantitatively assessing certain ecosystem services (Costanza and Folke, 1997; Daily, 1997). For example, assessment of water yield service tells us the capacity of water provision and flux process (Cai et al., 2017), which is important for wildlife, for example, fish and amphibian, ecosystems, and human well-being. The results also could be explained from source–sink theory. “Sources” are known as patches that support local species population growth and contribute positively to the ecological process (Chen et al., 2008). These critical ecosystem service areas play an importance role as “sources” to maintain and facilitate the ecological process. In addition, integrating internal landscape connectivity into ecological source identification to strengthen the importance of patches could be better understood from this theory on how patch connectivity promotes species migration across heterogenous landscapes (Hansen, 2011).

In our study, the principle of MCR was used in building an ecological resistance surface, whereas circuit theory was used to generate ecological corridors to overcome the shortcomings of the MCR model which often fails to clarify species spread randomly and identify a specific range of corridors (Peng et al., 2018b; Dai et al., 2021). It can be understood from circuit theory that electric currents share the random walk property, which is similar to the ecological flow of movements of animals and ecological process flux landscape patches (Mcrae et al., 2012). Ecological corridors could be generated in heterogenous landscapes by calculating the “resistance” or “current.” A similar result was found in Peng’s (Peng et al., 2018b) study, in which circuit theory was used to successfully identify ecological corridors and key ecological nodes in Yunnan Province. Our study proved that this method could also be applied effectively at the county scale.

KENPRA optimization is a critical step that moves land planning into restoration practice, which provides the scientific-based evidence for decision makers to implement restoration projects and environmental policy. Our study presented a comprehensive framework of KENPRA optimization by integrating ESP and LUC to adapt to rapid urbanization and future LUCC. This study showed that LUC was a promising indicator for quantitative assessment of potential LUCC. Incorporating LUC to KENPRA optimization could improve the identification of priority restoration areas. Jiang’s study (Jiang et al., 2020) demonstrated that areas with intense conflict had a high possibility of LUCC and tended to be converted into the type of higher competitiveness. Similarly, we found in our study that cultivated land areas, which were high intensive LUC areas, had high potential to be converted into built-up areas. This finding is consistent with the perspective of geography (Coccoli et al., 2018; Zou et al., 2019). LUC is a geographical phenomenon which originated from the scarcity of spatial land resources and the spillover of spatial functions. It is also an antagonistic phenomenon in the process of spatial resource allocation in the process of human–land relationships and conflict over spatial resources (Zou et al., 2020).

The results of the CA-Markov model indicated that by integrating LUC, KENPRA-based scenarios showed better response to urban expansion in 2030, when there will be less ecological corridor destruction and better conservation efficiency for ESP. Different from Zhang’s study in the Yangtze River Delta region, China (Zhang D. et al., 2020), in which they incorporated land use/cover predictions into ecological security source areas, our study integrated potential land use/cover change by integrating LUC into the overall ESP process to improve KENPRA identification. By doing so, it is more effective for decision makers to build ESP construction and KENPRA identification in urban or suburban areas where they have high pressure of rapid urbanization. Another advantage of integrating LUC is that the framework can be improved by taking natural conditions, socioeconomic activities, and environmental policy into account (Iojă et al., 2014). Compared with traditional ways of ESP construction which focused on natural conditions (Yu, 1996), this framework could be more resilient to future development change. Simultaneously, the forward-looking KENPRA identification will be helpful for managers in developing land planning agendas and in implementing restoration projects.

Our study explored how to bridge the gap between practical needs and research limitations by proposing a new approach of integrating ESP and LUC into KENPRA optimization, which will facilitate the ongoing and forthcoming China’s Ecosystem Restoration Projects at the basic administration unit of the county. It also could be applied at other regions and countries across different scales. In addition, our proposed framework of incorporating ESP and LUC, which combines essentially land spatial planning of land function suitability assessment and landscape pattern analysis, could be more approachable and efficient when applied in Chinese National Land Spatial Planning and Spatial Ecological Restoration Planning. It also could improve the understanding of the coupled human–nature system of suburbs by improving linkages among LUCC, ecosystem services, and land planning.

It is worth noting that the stimulation result indicates a significant increase of forest with cultivated land decline in the period of 2018–2030 under ENKAR-based scenarios. The stimulated LULC change could be explained by LULC constrain factor of ecological redline policy, which is carried out by Chinese government to strictly protect natural areas with high value of ecosystem services and ecological sensitive areas. The farmland in mountain regions could be used for afforestation with ecological redline policy and restoration projects. Therefore, undergoing and upcoming environmental project and land policy should fully consider in KENPRA for responding future LULC. In addition, based on our results from KENPRA identification and a field site visit, site-specific strategies for land management in the study area were proposed. Building and maintaining ecological corridors were suggested to be incorporated into urban planning. Developed land expansion was limited in PRA-1 to PRA-3 by delineating development control lines to limit built-up land expansion. Strict conservation policies, including the ecological redline policy, and prohibited activities of deforestation and mining were encouraged in core patches and ecological corridors in PRA-4 to PRA-6. Essential ecological restoration was suggested, including reclamation in PRA-5, wetland restoration in PRA-2, and building ecological corridor bridges in PRA-14 and PRA-15 to help wildlife migration when ecological corridors are interrupted by artificial roads. These refined strategies could better protect the ecological security pattern and balance ecosystem protection and economic development.

Our result is consistent with that of previous studies that suggests that considering land use transition in spatial planning is conducive to facilitating suburban spatial governance (Ge et al., 2020) and achieve SDGs (Vörösmarty et al., 2018). ESP was originally proposed in the urban context to facilitate landscape design and urban planning for conserving habitat and safeguarding important ecological networks (Yu, 1996; Waldheim, 2006). Though approximately 15.2% of the world’s terrestrial area has been established as protected areas to protect global wildlife habitats and ecological networks (Yang et al., 2019; Peng et al., 2020a), urbanization is a global trend with momentum. Therefore, our process to guide ESP establishment and restoration implementation in urban and suburban areas could effectively conserve and restore terrestrial and inland freshwater ecosystems and their services, which contribute to sustainable land use and the achievement of SDGs in semi-natural and artificial areas.

Our study was, however, limited to one single context studying a Chinese county. Future work could usefully apply this framework to other regions at different scales. With the increasing demand for “Landscape Planning for Ecological Restoration” and with ever more restoration programs being implemented, more such studies are needed. The proposed approach could be used to improve land planning and management by better adapting to rapid LUCC and could also be helpful in reducing landscape fragmentation and conserve wildlife habitat and biodiversity in semi-natural or artificial contexts (Peng et al., 2018a; Zhang D. et al., 2020). Furthermore, the spatial composition of ecological corridors need to be explored in the future studies; ecological corridors with different widths, for example, may have different influences on ESP maintenance (van Schalkwyk et al., 2020), which should be considered in optimizing KENPRA. Last, our research is limited by data availability. Future research would benefit from further spatiotemporal data enhancement of the variables analyzed here and the incorporation of variables absent in this study.