Delineating Urban Growth Boundaries with Ecosystem Service Evaluation

: China’s rapid urbanization over the past decades has been accompanied by ecological deterioration. This decline in the provision of vital ecosystem services now poses a signiﬁcant threat to urban area sustainability. Accordingly, the evaluation of ecosystem services has gained greater importance in ecological and sustainable development over the past decade. However, little information about ecosystem services is factored into urban planning and management decisions and limited studies to date have incorporated conservation prioritization when making decisions about urban growth boundaries. In this study, we proposed an initial framework to illustrate its application in Hangzhou. We modeled and mapped ﬁve ecosystem services (i.e., habitat quality as a proxy of biodiversity, carbon storage, water yield, sediment retention, nutrient retention) using the InVEST model and evaluated the overlaps among them. Zonation, a systematic conservation planning tool, was applied to explicitly spatialize conservation prioritization, and we proposed an analytical framework to deﬁne priority areas for ecosystem services conservation and delineated a rigid urban growth boundary. Our study integrated ecosystem service evaluations into the urban land-use decision-making process and addressed compromises in decisions regarding conservation prioritization.

( Figure 1). Today, it has approximately 9.86 million inhabitants and is expected to continue to grow rapidly in the near future. Hangzhou covers an area of approximately 16850 km , of which 4% is urban and rural settlements, 13% farmland, 63% forests, and 7% water bodies. The southwestern, western, and central areas of Hangzhou are hilly or mountainous, and covered predominantly by forests or forestry lands, within which is Qiandao Lake, the largest reservoir and important water supply in Zhejiang Province. The Qiandao Lake is the upper reach of Qiantang River that flows from the southwest to the northeast into the East Sea. The eastern and northeastern areas are a mosaic of urban built-up areas, agricultural land, and a multitude of rivers and wetlands. Important habitats, such as wetlands and forests, have disappeared in the suburbs because of urban expansion, and hilly areas in the southwest have suffered serious soil erosion attributable to unsuitable land use.

Data
The data used in ES assessments in the InVEST model are listed in Table 1. Land-use and landcover (LULC) maps are the most important input to all modules in InVEST. We obtained a raster grid LULC map (30 × 30 m grid cells) from the Hangzhou Planning Bureau that was the result of the 2016 annual land-use investigation ( Figure 2).

Data
The data used in ES assessments in the InVEST model are listed in Table 1. Land-use and land-cover (LULC) maps are the most important input to all modules in InVEST. We obtained a raster grid LULC map (30 × 30 m grid cells) from the Hangzhou Planning Bureau that was the result of the 2016 annual land-use investigation ( Figure 2).  [27]. * HQ: habitat quality; CS: carbon storage; WY: water yield; SDR: sediment retention; NR: nutrient retention. In general, InVEST can be run on spatial units of any resolution, although finer-resolution data are recommended if available. In this study, InVEST modules were run using the "30 × 30 m" resolution data. However, this resolution was not feasible computationally in Zonation, as there are more than 18 million grid cells. To balance computational feasibility, network compactness [28],  In general, InVEST can be run on spatial units of any resolution, although finer-resolution data are recommended if available. In this study, InVEST modules were run using the "30 × 30 m" resolution data. However, this resolution was not feasible computationally in Zonation, as there are more than 18 million grid cells. To balance computational feasibility, network compactness [28], and model accuracy, we chose to run Zonation at the resolution of the "100 × 100 m" grid cell. All InVEST outputs were aggregated to this cell size.
Modeling and mapping ES contribute to our understanding of ecosystems and, thus, is an important way in which to incorporate ES into decision-making. The ES assessment depends on the availability of spatially explicit information on the state of the ES and the trends [29]. Various models, such as InVEST, SolVES (Social Values for Ecosystem Services [30]), and ARIES (Artificial Intelligence for Ecosystem Services [30]), have been developed and applied widely. The InVEST model is a suite of tools that the Stanford and Minnesota Universities, World Wildlife Fund, and The Nature Conservancy developed. As Nelson et al. [31] stated, "InVEST is a suite of service models that use production functions to convert maps of land-use and land-cover (LULC), land management, and biophysical conditions into maps of service supply."

InVEST
In this study, InVEST was applied to model carbon storage, water yield, sediment delivery, nutrient delivery, and habitat quality (as a proxy of biodiversity). (1) The carbon storage module uses LULC maps together with stocks in four carbon pools (i.e., aboveground and belowground biomass, soil, and dead organic matter) to estimate the amount of carbon stored in a landscape currently. (2) The water yield in a given pixel was calculated as precipitation minus evapotranspiration (ET). The ET partition of the annual water balance was approximated using a Budyko curve, plant available water content, average annual precipitation, and a seasonality factor, Z, that represents seasonal rainfall's amount and distribution [32]. (3) The sediment retention module uses the Revised Universal Soil Loss Equation (RUSLE) [33] to predict the average annual rate of soil erosion in a particular area. The rate of soil erosion is a function of the area's LULC, soil type, rainfall intensity, topography, crop/vegetation, and management (C factor), and support practice (P factor). Regions with lower potential soil losses received higher scores for the sediment retention service. (4) The amounts of nitrogen (N) retained in each pixel were calculated by aggregating the amount of nutrient (N) load input to a pixel minus those being transported to the next pixel. Finally, (5) habitat quality is a function of four factors: each threat's relative effect, each habitat type's relative sensitivity to each threat, the distance between habitats and sources of threats, and the degree to which the land is protected legally [23]. In this study, six threats that represent anthropogenic drivers in human-dominated landscapes were considered (i.e., urban, rural residential, railway, primary roads, secondary roads, and cropland). Based on LULC maps and threats maps, the model generated a habitat quality map that showed habitat quality values between 0 and 1, in which 1 represents the highest habitat quality or the highest possible level of biodiversity.

Zonation
Zonation is the tool used most often in conservation prioritization. It can operate on species, ES, or any such biodiversity feature, and can be applied to landscapes with up to tens of millions of elements of conservation feature data [20]. Conservation prioritization with Zonation also can account for biodiversity, ES, and socioeconomic factors [7]. Zonation begins with the assumption that protecting everything would be the best conservation strategy. It then proceeds to rank sites iteratively, removing the spatial unit that leads to the smallest aggregate marginal loss in conservation features at each step while accounting for the total and distributions of features, weights given to features, and feature-specific connectivity. In this process, the least useful sites are ranked the lowest and the areas most valuable for conservation are ranked the highest. Zonation has three different cell removal rules depending on the way marginal loss is calculated: additive benefit function (ABF), core-area zonation (CAZ), and target-based zonation (TBZ) (more details about cell removal rules can be found in the Zonation user manual [34]. We used the ABF in Zonation, as it produces a higher performance on average over all conservation features (compared to the core-area removal rule) and no strict target for conservation is required (compared to the target-based planning cell removal rule) [34]. First, the program calculates each feature's loss of representation, as cell i is removed. The cell's δ i (marginal loss in conservation value) simply is the sum over feature-specific declines in value (∆V j ) following the loss of cell i, as shown in Equation (1): in which R j (S) is the representation of feature j in the remaining set of sites S, and {S − i} indicates the set of remaining cells minus cell i; V j , the value of feature j in the reserve network, is an increasing function of R j ; w j is the weight of feature j, and c i is the cost of adding cell i to the reserve network. The cell with the smallest δ i value is removed. Zonation produces a priority rank map and performance curve for each run. According to the iterative cell ranking's order, each cell in the priority rank map has a value between 0 and 1: cells with values close to 0 were removed first because of their low conservation value; while cells with a value close to 1 were retained because of their high conservation values [7]. The performance curve shows the proportion of features remaining in the landscape.
Connectivity is an important factor in conservation prioritization. Small and isolated patches are perceived to be as unfavorable for species' persistence and migration. Moreover, highly fragmented conservation sites are hard to manage and expensive [34]. Well-connected and moderately aggregated sites will be beneficial to conservation actions. In Zonation, some settings account for connectivity and aggregation. We chose the boundary length penalty (BLP), as it is the most common way of dealing with connectivity and aggregation in spatial prioritization. When using the BLP, a penalty is given for the boundary length of the reserve, and the marginal loss δ * i is calculated as below: in which δ i is the marginal loss in Equation (1); ∆(BL/A) is the change in the ratio of boundary length (BL) to area (A) of the reserve network following removal of cell i; β is a constant defining the strength of the boundary length penalty. The appropriate value for β needs to be found by experimentation.
In this study, we found that β = 0.1 could generate a favorable reserve network.

Correlation and Priority Area Overlap Analysis
We developed two methods to evaluate the spatial correspondence of ES. First, we calculated the Pearson's r among ES to measure their correlations. Second, we performed the calculation in Equation (3) to evaluate the pair-wise overlap among individual ES priority areas, in which o e is the expected overlap, o a is the actual overlap, n c is the number of cells that co-occur for each combination of services, n T is the total number of cells, and f is the fraction of cells in the priority area. If o a is larger than o e , these two services' priority areas have positive spatial associations; if o a is smaller than o e , they have negative spatial associations.

Spatial Distribution of Ecosystem Services
Using the data in Table 1 as inputs in InVEST, we quantified habitat quality, carbon storage, water yield, sediment retention, and nitrogen (N) retention in the InVEST model, as shown in Figure 3A-E. The nutrient retention module output was the result of the nutrient load and nutrient export in each pixel, and the quantity of N retention was calculated as N load minus N export. The absolute amount of N retention was directly related to the amount of N load in the landscape and high values for N retention were in farmland, which is the main source of N load. As in conservation prioritization, what really matters is the relative capacity rather than the absolute amount. Therefore, we calculated the N retention efficiency, as shown in Equation (4): In this way, we obtained the spatial distribution of N retention efficiency across the landscape, as shown in Figure 3F.
Ecosystem services have different spatial distributions and share similarities at the same time. The southwestern areas of Hangzhou are of high value for multiple services and the northeastern areas are of low value for many. Generally, minor human activities, natural forest cover, high precipitation, and rough topography characterize the mountainous regions. These factors account for high values in habitat quality, carbon storage, and water yield. High values for sediment retention were found largely in the valleys and gentle slopes in the southern, western, and central areas of Hangzhou. Similarly, intensive human activities, impervious surface, and low natural land cover characterize urban areas and, thus, they had low values for all ES. For N retention, the distribution of high-value In this way, we obtained the spatial distribution of N retention efficiency across the landscape, as shown in Figure 3F. Ecosystem services have different spatial distributions and share similarities at the same time. The southwestern areas of Hangzhou are of high value for multiple services and the northeastern areas are of low value for many. Generally, minor human activities, natural forest cover, high precipitation, and rough topography characterize the mountainous regions. These factors account for high values in habitat quality, carbon storage, and water yield. High values for sediment retention were found largely in the valleys and gentle slopes in the southern, western, and central areas of Hangzhou. Similarly, intensive human activities, impervious surface, and low natural land cover characterize urban areas and, thus, they had low values for all ES. For N retention, the distribution of high-value areas was consistent with that of farmland, yet high values for N retention efficiency appeared in smooth topography and natural vegetation cover.
We calculated the proportion of ES in each LULC type, as shown in Table 2. Forest is the dominant LULC type and covers 61% of Hangzhou. Forest demonstrated significant values for habitat quality (average: 0.91) and N retention efficiency (average: 0.86), and accounted for a great proportion of carbon storage, sediment retention, and water yield, as well as relatively high values for N retention. Agricultural land showed high values for N retention, while its retention efficiency was much lower (0.68). We calculated the proportion of ES in each LULC type, as shown in Table 2. Forest is the dominant LULC type and covers 61% of Hangzhou. Forest demonstrated significant values for habitat quality (average: 0.91) and N retention efficiency (average: 0.86), and accounted for a great proportion of carbon storage, sediment retention, and water yield, as well as relatively high values for N retention. Agricultural land showed high values for N retention, while its retention efficiency was much lower (0.68). To estimate ES correlations, we calculated correlations (Pearson's r) among ES values, as shown in Table 3. Spatial correlations among the five ES were complex, ranging from low negative (−0.47) to high positive (0.84). The N retention was negatively correlated with all other services. The highest positive correlation was between carbon storage and habitat quality (0.84), followed by that of carbon storage and water yield (0.71). The highest negative correlation was between N retention and habitat quality (−0.47), followed by N retention and carbon storage (−0.45). When negative correlations were excluded, habitat quality and carbon storage showed medium to high positive correlation with other services, reflecting the importance of natural vegetation cover (forest particularly) to carbon storage, water provision, soil retention, and nutrient filtration; N retention efficiency showed low to medium correlations with other services, with an overall correlation of 0.36. Existing protected areas (PAs) in Hangzhou (Figure 4), including forest parks and natural reserves, ecological forests, and ecological red lines, cover an area of 5720 km 2 -34% of the total area. They are scattered spatially and designated for different conservation objectives, such as species-specific reserves, drinking water supply (35% of PAs are in Chun'an county to protect the Qiandao Lake as a drinking water source), and historical cultural sites (such as the West Lake scenic area). To evaluate the individual representation of the ES within PAs, we calculated the proportion of the ES retained in it. The results are shown in Table 4. Habitat quality in PAs was high (0.85), and PAs provided more than half of the total sediment retention service and more than 40% of total water yield and carbon storage. The N retention had low representation in PAs (19%); nevertheless, N retention efficiency was high (0.82), which we believe was attributable to a lower proportion of farmland as the main source of nitrogen, as well as a high level of filtration by vegetation in PAs.

Identifying Priority Areas for Individual Ecosystem Services
The identification of priority areas for individual ES is helpful in targeting policy making for specific ES. We performed a prioritization in Zonation for five individual ES. The N retention efficiency was used as a representative of nutrient retention service in the following analysis. Priority rank maps for five ES are shown in Figure 5 with five color gradations. We took the top-ranked 60% as priority areas for each ES. Habitat quality and carbon storage shared priority areas in the south and mid-west; priority areas for water yield were concentrated in the southwest and in part of the west; priority areas for sediment retention appeared to be scattered and did not have larger patches, but the color gradation showed that the southern and western areas ranked higher, while priority areas for nutrient retention were largely in the central and western areas.  The ES' performance curves ( Figure 6) showed that, as cells were continuously removed from the landscape, the proportion of ES remaining declined at different rates. Compared to the curve of landscape loss plotted in the black, dashed line, the sediment and nutrient retention services shared the most similar trends and declined at a steady rate; the other three ES' decline first appeared slowly and rapidly thereafter, and the turning point occurred after 40% of cells were removed from the original landscape. The curves indicated that the top-ranked 60% of the landscape was capable to provide a large proportion of each individual ES. Ranked from high to low, the remaining proportions of ES in the landscape were carbon storage, habitat quality, water yield, nutrient retention, and sediment retention. The proportions of ES remaining were all greater than 80% when 20% of the landscape was removed, among which carbon storage and water yield maintained a proportion of approximately 90%. When 40% of the landscape was removed, more than 80% of carbon storage, 70% of water yield and habitat quality, and 65% of nutrient and sediment retention We performed a pair-wise network overlap analysis among individual ES priority areas, as shown in Table 5. The priority areas of habitat quality and carbon storage had the highest overlap (96.34%); the overlap between priority areas of habitat quality, carbon storage, sediment retention, and water yield were all approximately 80%. The overlap between nutrient retention and other ES was slightly lower. The overall overlap was 53.07%, greater than expected (36%). The results indicated that ES exhibited different degrees of synergy and it was possible to align them in a network. The ES' performance curves ( Figure 6) showed that, as cells were continuously removed from the landscape, the proportion of ES remaining declined at different rates. Compared to the curve of landscape loss plotted in the black, dashed line, the sediment and nutrient retention services shared the most similar trends and declined at a steady rate; the other three ES' decline first appeared slowly and rapidly thereafter, and the turning point occurred after 40% of cells were removed from the original landscape. The curves indicated that the top-ranked 60% of the landscape was capable to provide a large proportion of each individual ES. Ranked from high to low, the remaining proportions of ES in the landscape were carbon storage, habitat quality, water yield, nutrient retention, and sediment retention. The proportions of ES remaining were all greater than 80% when 20% of the landscape was removed, among which carbon storage and water yield maintained a proportion of approximately 90%. When 40% of the landscape was removed, more than 80% of carbon storage, 70% of water yield and habitat quality, and 65% of nutrient and sediment retention services were retained in the landscape.

Delineating Urban Growth Boundaries based on Conservation Prioritization
To evaluate the ability to conserve ES in existing PAs, we performed two runs in the Zonation analysis and prioritized PAs in the second run. In situation one, we identified priority areas based completely on the results of the ES assessment; in situation two, we identified priority areas for conservation actions outside of the existing PA network, by using a hierarchical mask that identified existing PAs' location. This guaranteed that the highest priorities were located in existing PA. Then, we analyzed and compared the extent to which ES were retained in conservation networks in the two situations. Figure 7 shows the priority rank map Zonation produced in two runs, with a color gradation indicating the value between zero and one assigned to each cell. The southwestern, western, and southern parts of Hangzhou ranked higher in both situations, which was consistent with the spatial distribution of high ES values. However, in situation one, areas in the west-central areas ranked higher, but the Qiandao Lake in Chun'An County ranked lower. The eastern and northeastern areas obtained very low ranks in both situations, largely 0~20%. However, in situation one, areas ranked between 0~20% were concentrated more in the east and north, while they were relatively scattered in situation two, with several small fragments distributed in the south and west. Figure 8 shows the fraction of the landscape remaining in the solution plotted against the mean and minimum proportions of remaining ES. It provides a quantitative measure of the way the existing PAs compared to the unconstrained Zonation solution. While they were similar, Zonation without PAs generated higher proportions of remaining ES when the fraction of landscape was in the range of 20%-70%.

Delineating Urban Growth Boundaries Based on Conservation Prioritization
To evaluate the ability to conserve ES in existing PAs, we performed two runs in the Zonation analysis and prioritized PAs in the second run. In situation one, we identified priority areas based completely on the results of the ES assessment; in situation two, we identified priority areas for conservation actions outside of the existing PA network, by using a hierarchical mask that identified existing PAs' location. This guaranteed that the highest priorities were located in existing PA. Then, we analyzed and compared the extent to which ES were retained in conservation networks in the two situations. Figure 7 shows the priority rank map Zonation produced in two runs, with a color gradation indicating the value between zero and one assigned to each cell. The southwestern, western, and southern parts of Hangzhou ranked higher in both situations, which was consistent with the spatial distribution of high ES values. However, in situation one, areas in the west-central areas ranked higher, but the Qiandao Lake in Chun'An County ranked lower. The eastern and northeastern areas obtained very low ranks in both situations, largely 0~20%. However, in situation one, areas ranked between 0~20% were concentrated more in the east and north, while they were relatively scattered in situation two, with several small fragments distributed in the south and west. Figure 8 shows the fraction of the landscape remaining in the solution plotted against the mean and minimum proportions of remaining ES. It provides a quantitative measure of the way the existing PAs compared to the unconstrained Zonation solution. While they were similar, Zonation without PAs generated higher proportions of remaining ES when the fraction of landscape was in the range of 20-70%. and west. Figure 8 shows the fraction of the landscape remaining in the solution plotted against the mean and minimum proportions of remaining ES. It provides a quantitative measure of the way the existing PAs compared to the unconstrained Zonation solution. While they were similar, Zonation without PAs generated higher proportions of remaining ES when the fraction of landscape was in the range of 20%-70%.  The proportion of a given ES can vary considerably as cells are removed from the landscape. This proportion is determined by an environmental service's distribution and the extent to which its provision area overlaps with other ES. An ES with a priority area nested within the distribution of other ES will have the highest proportion of service remaining when any given proportion of the landscape is selected. The lowest proportions are retained for an ES that is distributed widely and evenly, which results in an approximately linear decline in the landscape retained, as cells are removed [21]. Figure 9 shows the proportion of the landscape remaining in the solution plotted against the proportion of ES remaining in situation two. Compared with the performance curve in the individual ES prioritization, the curve of carbon storage and habitat quality showed a swifter decline. Moreover, the curves' turning point occurred when 30% of cells were removed from the original landscape, which was earlier than that in individual ES prioritization (40%), reflecting the existence of compromises among ES. The proportions of ES remaining all were over 80% when 20% of the landscape was removed, among which carbon storage, water yield, and habitat quality maintained a proportion of approximately 90%. When 40% of the landscape was removed, more than 70% of carbon storage, water yield, and habitat quality, and 65% of sediment retention services were retained in the landscape. The curves indicated that the top-ranked 60% of the landscape was able to provide more than 60% of multiple ES. The proportion of a given ES can vary considerably as cells are removed from the landscape. This proportion is determined by an environmental service's distribution and the extent to which its provision area overlaps with other ES. An ES with a priority area nested within the distribution of other ES will have the highest proportion of service remaining when any given proportion of the landscape is selected. The lowest proportions are retained for an ES that is distributed widely and evenly, which results in an approximately linear decline in the landscape retained, as cells are removed [21]. Figure 9 shows the proportion of the landscape remaining in the solution plotted against the proportion of ES remaining in situation two. Compared with the performance curve in the individual ES prioritization, the curve of carbon storage and habitat quality showed a swifter decline. Moreover, the curves' turning point occurred when 30% of cells were removed from the original landscape, which was earlier than that in individual ES prioritization (40%), reflecting the existence of compromises among ES. The proportions of ES remaining all were over 80% when 20% of the landscape was removed, among which carbon storage, water yield, and habitat quality maintained a proportion of approximately 90%. When 40% of the landscape was removed, more than 70% of carbon storage, water yield, and habitat quality, and 65% of sediment retention services were retained in the landscape. The curves indicated that the top-ranked 60% of the landscape was able to provide more than 60% of multiple ES. We chose to conduct conservation prioritization based on the results of situation two, Zonation with existing PAs, for two main reasons: first, the performance curve ( Figure 8) shows that when a large fraction (more than 60%) of the landscape was reserved, the minimum and mean proportions of multiple ES were very close; second, as existing PAs already have a considerable foundation for conservation, it would be easier to expand them than to designate a new reserve. Figure 10 shows the expansion of PAs. Newly expanded areas were distributed largely in the mid-west, which had high values for multiple ES and was adjacent to existing PAs. According to the ability to provide ES, land can be zoned as suitable either for conservation or for development ( Figure 11). We identified areas vital for providing multiple ES such as the topranked 60% of the landscape, potential ES provision areas such as those ranked between 20%-40%, and areas suitable for urban development such as those ranked between 0-20%. Rigid UGBs were delineated in areas suitable for development, as shown in Figure 12. As current urban areas cover only 5.6% of Hangzhou, there is considerable room for future urban development. Thus, urban development should be confined within UGB and take place in lower-ranked cells in sequence to mitigate its adverse effects on ES. We chose to conduct conservation prioritization based on the results of situation two, Zonation with existing PAs, for two main reasons: first, the performance curve ( Figure 8) shows that when a large fraction (more than 60%) of the landscape was reserved, the minimum and mean proportions of multiple ES were very close; second, as existing PAs already have a considerable foundation for conservation, it would be easier to expand them than to designate a new reserve. Figure 10 shows the expansion of PAs. Newly expanded areas were distributed largely in the mid-west, which had high values for multiple ES and was adjacent to existing PAs. We chose to conduct conservation prioritization based on the results of situation two, Zonation with existing PAs, for two main reasons: first, the performance curve ( Figure 8) shows that when a large fraction (more than 60%) of the landscape was reserved, the minimum and mean proportions of multiple ES were very close; second, as existing PAs already have a considerable foundation for conservation, it would be easier to expand them than to designate a new reserve. Figure 10 shows the expansion of PAs. Newly expanded areas were distributed largely in the mid-west, which had high values for multiple ES and was adjacent to existing PAs. According to the ability to provide ES, land can be zoned as suitable either for conservation or for development ( Figure 11). We identified areas vital for providing multiple ES such as the topranked 60% of the landscape, potential ES provision areas such as those ranked between 20%-40%, and areas suitable for urban development such as those ranked between 0-20%. Rigid UGBs were delineated in areas suitable for development, as shown in Figure 12. As current urban areas cover only 5.6% of Hangzhou, there is considerable room for future urban development. Thus, urban development should be confined within UGB and take place in lower-ranked cells in sequence to mitigate its adverse effects on ES. According to the ability to provide ES, land can be zoned as suitable either for conservation or for development ( Figure 11). We identified areas vital for providing multiple ES such as the top-ranked 60% of the landscape, potential ES provision areas such as those ranked between 20-40%, and areas suitable for urban development such as those ranked between 0-20%. Rigid UGBs were delineated in areas suitable for development, as shown in Figure 12. As current urban areas cover only 5.6% of Hangzhou, there is considerable room for future urban development. Thus, urban development should be confined within UGB and take place in lower-ranked cells in sequence to mitigate its adverse effects on ES.

Discussion
The identification of conservation priority areas and management strategies based on ES is invaluable for urban sustainability. Conservation prioritization is often designed to meet conservation targets with priority areas of minimum size or to optimize the conservation targets and maximize suitability. However, it is difficult to determine an appropriate target. As is often the case, targets are set too low to constrain development, and always lead to island effect protection. In this study, we did not set a quantitative target for conservation in advance. Instead, we ranked sites according to their ability to provide multiple ES. With this information, land managers and policymakers will know whether their decision about land use will result in the reduction of vital ES. The systematic conservation tools used here guarantee a suitable reserve network with compact patches and continuous boundaries that allow effective conservation action.
The synergies and compromises among multiple ES are complex and add to the uncertainties in conservation prioritization. A better understanding of ES' spatial patterns and compromises is necessary to make cost-benefit analyses and sustainable decisions in conservation [35]. Generally, trade-offs often exist between providing services (e.g., food or timber production) and other ES (supporting and regulating services). However, it also depends on the region and the type of

Discussion
The identification of conservation priority areas and management strategies based on ES is invaluable for urban sustainability. Conservation prioritization is often designed to meet conservation targets with priority areas of minimum size or to optimize the conservation targets and maximize suitability. However, it is difficult to determine an appropriate target. As is often the case, targets are set too low to constrain development, and always lead to island effect protection. In this study, we did not set a quantitative target for conservation in advance. Instead, we ranked sites according to their ability to provide multiple ES. With this information, land managers and policymakers will know whether their decision about land use will result in the reduction of vital ES. The systematic conservation tools used here guarantee a suitable reserve network with compact patches and continuous boundaries that allow effective conservation action.
The synergies and compromises among multiple ES are complex and add to the uncertainties in conservation prioritization. A better understanding of ES' spatial patterns and compromises is necessary to make cost-benefit analyses and sustainable decisions in conservation [35]. Generally, trade-offs often exist between providing services (e.g., food or timber production) and other ES (supporting and regulating services). However, it also depends on the region and the type of

Discussion
The identification of conservation priority areas and management strategies based on ES is invaluable for urban sustainability. Conservation prioritization is often designed to meet conservation targets with priority areas of minimum size or to optimize the conservation targets and maximize suitability. However, it is difficult to determine an appropriate target. As is often the case, targets are set too low to constrain development, and always lead to island effect protection. In this study, we did not set a quantitative target for conservation in advance. Instead, we ranked sites according to their ability to provide multiple ES. With this information, land managers and policymakers will know whether their decision about land use will result in the reduction of vital ES. The systematic conservation tools used here guarantee a suitable reserve network with compact patches and continuous boundaries that allow effective conservation action.
The synergies and compromises among multiple ES are complex and add to the uncertainties in conservation prioritization. A better understanding of ES' spatial patterns and compromises is necessary to make cost-benefit analyses and sustainable decisions in conservation [35]. Generally, trade-offs often exist between providing services (e.g., food or timber production) and other ES (supporting and regulating services). However, it also depends on the region and the type of services selected. In this study, habitat quality as a proxy of biodiversity is a supporting service, while carbon storage and sediment and nutrient retention are regulating services. Water yield was the only providing service in our assessment. The Pearson's correlations among ES showed that habitat quality, carbon storage, water yield, and sediment retention services were correlated positively, yet N retention was negatively correlated with all other ES. The positive correlation between water yield and other ES was unsurprising, as water yield is much more relevant to precipitation. According to the multi-year average precipitation data we obtained, the annual precipitation in the mountainous area in the southwest of Hangzhou is much higher than that in the urbanized area in the northeast. Negative correlations were eliminated when N retention efficiency was introduced to represent the relative capacity rather than the absolute amount of nutrient retention. According to Chan et al. [19], the correlation among services does not convey the full extent to which conservation activities would align, because correlations at low levels of provision are unlikely to be selected as priority areas, and correlations among services neglect grid cells' spatial distributions needed to meet suitability for conservation. Thus, the overlap between priority areas is more relevant for conservation, and the results showed that there was a certain degree of compromise among different ES. Conservation prioritization with Zonation considered the synergies, trade-offs, and suitability, and generated a solution in which more than 65% of multiple ES can be preserved in 60% of landscape remaining, indicating the potential to design a cost-effective network despite compromises between ES.
Based on the Zonation prioritization, we identified areas for vital ES provision, potential ES provision, and urban development and delineated rigid UGBs in Hangzhou. Moreover, the information contained in rank maps can be used to determine land-use type and the time sequence of future development. A comparison between the Hangzhou UGB 2016 in eight districts and the rank maps ( Figure 13) can identify areas with higher conservation values within the UGB and areas with lower values outside the UGB. Figure 13 shows that some patches in the northwest and north were inside the UGB, while lower-ranked areas in the eastern area were outside the UGB. The UGB can be adjusted accordingly to have the minimum adverse effect on ecosystem provision. This provides a possible way in which ES can be integrated into spatial planning decisions. services selected. In this study, habitat quality as a proxy of biodiversity is a supporting service, while carbon storage and sediment and nutrient retention are regulating services. Water yield was the only providing service in our assessment. The Pearson's correlations among ES showed that habitat quality, carbon storage, water yield, and sediment retention services were correlated positively, yet N retention was negatively correlated with all other ES. The positive correlation between water yield and other ES was unsurprising, as water yield is much more relevant to precipitation. According to the multi-year average precipitation data we obtained, the annual precipitation in the mountainous area in the southwest of Hangzhou is much higher than that in the urbanized area in the northeast. Negative correlations were eliminated when N retention efficiency was introduced to represent the relative capacity rather than the absolute amount of nutrient retention. According to Chan et al. [19], the correlation among services does not convey the full extent to which conservation activities would align, because correlations at low levels of provision are unlikely to be selected as priority areas, and correlations among services neglect grid cells' spatial distributions needed to meet suitability for conservation. Thus, the overlap between priority areas is more relevant for conservation, and the results showed that there was a certain degree of compromise among different ES. Conservation prioritization with Zonation considered the synergies, trade-offs, and suitability, and generated a solution in which more than 65% of multiple ES can be preserved in 60% of landscape remaining, indicating the potential to design a costeffective network despite compromises between ES. Based on the Zonation prioritization, we identified areas for vital ES provision, potential ES provision, and urban development and delineated rigid UGBs in Hangzhou. Moreover, the information contained in rank maps can be used to determine land-use type and the time sequence of future development. A comparison between the Hangzhou UGB 2016 in eight districts and the rank maps ( Figure 13) can identify areas with higher conservation values within the UGB and areas with lower values outside the UGB. Figure 13 shows that some patches in the northwest and north were inside the UGB, while lower-ranked areas in the eastern area were outside the UGB. The UGB can be adjusted accordingly to have the minimum adverse effect on ecosystem provision. This provides a possible way in which ES can be integrated into spatial planning decisions. Our study had some limitations. First, a LULC map was used as a proxy to quantify ES provision; however, the existing systems of land-cover classifications in China simplify the human influence on the landscape somewhat. Moreover, not all ES can be assessed adequately based on land-cover classes because they depend on landscape elements' particular qualities. Data on the Our study had some limitations. First, a LULC map was used as a proxy to quantify ES provision; however, the existing systems of land-cover classifications in China simplify the human influence on the landscape somewhat. Moreover, not all ES can be assessed adequately based on land-cover classes because they depend on landscape elements' particular qualities. Data on the actual composition, management practices and intensity of use, and knowledge of their effects on ES provision should be involved in further assessments and decision-making [36]. Second, ES in many landscape features were underestimated in our framework because of an incomplete selection of ES and data restrictions. For example, the riparian areas we mapped did not have significant values for carbon storage, water yield, and sediment and nutrient retention, but they are very important in recreational, cultural, and local climate regulation services provision, all of which have high values, particularly in urban regions. Third, the InVEST model can provide a rapid and relatively credible assessment under time and data constraints; however, it is based on general assumptions not on site-specific quantified relations. Many parameters in the model depend on the literature or expert knowledge and, thus, the results can suffer from subjectivity to a certain extent. Moreover, some of our inputs were obtained from secondary data; and some data had a coarse resolution (e.g., climate data), factors that introduce error in the ES assessments. They can be improved by using a more robust method, finer resolution data, empirical parameters in ES mapping, and a more direct biodiversity index in biodiversity mapping. Fourth, we neglected off-site effects in our framework, which might lead to an overestimation of carbon storage's effects, as it is a global service [35].

Conclusions
It is highly important to identify priority areas and formulate urban growth management policies to conserve ES threatened by rapid urbanization. In this study, we proposed a framework to integrate ES into urban spatial planning in Hangzhou. We modeled and mapped five ES (i.e., habitat quality as a proxy of biodiversity, carbon storage, water yield, sediment retention, and nutrient retention) with the InVEST model. Zonation, a systematic conservation planning tool, was applied to provide an explicit spatial representation for conservation prioritization. We proposed an analytical framework to define priority areas for ES conservation and delineated a rigid UGB. Our study addressed the integration of ecosystem service assessments into the decision-making process in urban land use and dealt with trade-offs among ES in conservation prioritization. According to priority rank maps, policy-makers or other stakeholders can easily identify areas suitable for further exploitation and those that require conservation. The identification of priority areas can provide a basis for establishing guidelines for decision-makers who face different policy options for urban expansion and natural capital management. With such tools, urban development can be guided properly and adjusted to ensure urban sustainability [11].
This study was a beneficial attempt to integrate ES into urban planning decision making and implementation. Despite the great progress made in ES conservation theories, techniques, and tools, most studies have simply identified and evaluated priority areas, and few have integrated conservation prioritization into urban planning practice and policy making. Traditionally, conservation planning has been inclined to identify areas abundant in biodiversity or valuable in providing certain ES (e.g., drinking water, recreation sites), and establish PAs as a barrier to threats. However, these PAs are far from sufficient to constitute a comprehensive conservation network and restrain urban sprawl. Consequently, landscape fragmentation and ecological space encroachment continue to worsen and will ultimately lead to islanding protection. In this study, we applied the theories and methods of conservation planning in an urban context. We carried out conservation prioritization for multiple ES and presented its potential to guide urban development and land-use policy formulation. In addition, our framework introduced scientific criteria into UGB delineation. For a long while, planners have overlaid important ecological features that are weighted subjectively in GIS and ranked sites according to their overlay values. This method of UGB delineation suffers from a high degree of subjectivity and arbitrariness in selecting and weighting ecological features and leads to a fragmented solution that fails to consider compactness, connectivity, and other landscape structural factors, which, in turn, reduces conservation action's effectiveness. The involvement of ES assessment and a spatially explicit analysis approach to UGB delineation demonstrated in this study can provide urban planners with more sophisticated information than simple rules of thumb [21]. Evaluating multiple ES in a widely used ecological model such as InVEST can provide a more scientific information base than weighting ecological features subjectively and lay a good foundation for the following policy determination. The Zonation spatial prioritization tool accounted for connectivity and aggregation, which are crucial for effective ecosystem management. Moreover, quantification and visualization of ES distributions and explicit representation of synergies and trade-offs can serve as a great communication tool for policy-makers, beneficiaries, and the public, and can provide clear guidance for the implementation of conservation actions. As the data we used in this case study are available publicly, it can be replicated to delineate rigid UGB in different areas.
Because of the inherent complexity of ecosystems, the interrelationships among ES remain a key point in ES research. Decision-makers need to acquire as much knowledge as possible about the relationships among ES in order that synergy, rather than trade-offs, can be obtained from ecosystem management decisions. The selected five ES appeared to be positively correlated in our analysis, which is favorable for the following spatial conservation prioritization. Although it is still challenging to obtain the best synergy given the different degrees of overlap among individual ES priority areas, Zonation turned out to be effective in addressing this problem and made prioritization more transparent. Nevertheless, it can be expected that involving other ES (such as food production, wood supply, cultural entertainment, etc.) would bring in different trade-offs in the decision making and make it more complicated. Moreover, we addressed the trade-offs between ES only in a spatial sense, and concentrated on ecosystem conditions and potential at present. Addressing trade-offs in a temporal sense would complement our knowledge base of ecosystem dynamics and tendencies and facilitate better development decision making regarding ecosystem management, which has not yet been covered in our study.
There are many ways in which further studies can enhance our framework. First, a more accurate and robust approach for multiple ES assessments is needed, which is challenging in urban regions characterized by a high degree of complexity and multi-functionality. Second, our study addressed only the biophysical assessments that help identify areas to ensure effective ES conservation [14]. As cities include complex political, economic, cultural, and ecological factors inherently, their development corresponds with various socioeconomic factors, not simply ecology. This highlights the need to combine assessments of the biophysical, economic, and social contexts that factor implementation opportunities and constraints into stakeholders' strategy development, implementation, and management. Third, a more rigorous approach to assessing ES would involve a complete consideration of actual supply and demand, which is possible only when detailed socio-economic data is accessible to us, such as the scale, scope, and beneficiaries of each ES. Actually, the demand of regional, national or even global ES can hardly be quantified at municipal level, such as fresh water provision and climate regulation in this study. We presented an example of integrating ES information base into traditional planning tools [37], which was an initial stage of ES approach application. Integrating ES into urban planning necessitates a better understanding of urban ecosystems and human well-being, a further development of ES quantification method, a complete consideration of demands, stakeholders, and beneficiaries, and more planning practices both at the local and regional scales.