Optimizing Green-Gray Infrastructure for Non-Point Source Pollution Control under Future Uncertainties

Non-Point Source Pollution (NPS) caused by polluted and untreated stormwater runoff discharging into water bodies has become a serious threat to the ecological environment. Green infrastructure and gray infrastructure are considered to be the main stormwater management measures, and the issue of their cost-effectiveness is a widespread concern for decision makers. Multi-objective optimization is one of the most reliable and commonly used approaches in solving cost-effectiveness issues. However, many studies optimized green and gray infrastructure under an invariant condition, and the additional benefits of green infrastructure were neglected. In this study, a simulation-optimization framework was developed by integrated Stormwater Management Model (SWMM) and Non-dominated Sorting Genetic Algorithm (NSGA-II) to optimize green and gray infrastructure for NPS control under future scenarios, and a realistic area of Sponge City in Nanchang, China, was used as a typical case. Different levels of additional benefits of green infrastructure were estimated in the optimizing process. The results demonstrated that green-gray infrastructure can produce a co-benefit if the green infrastructure have appropriate Value of Additional Benefits (VAB), otherwise, gray infrastructure will be a more cost-effectiveness measure. Moreover, gray infrastructure is more sensitive than green infrastructure and green-gray infrastructure under future scenarios. The findings of the study could help decision makers to develop suitable planning for NPS control based on investment cost and water quality objectives.


Introduction
Rapid urbanization has resulted in increased impervious surface areas and increased surface pollutant load caused by increasing human production activities [1,2]. Rainfall washes off the surface which contains large pollutant load, resulting in large amounts of pollutants being carried in runoff. The untreated runoff that is discharged directly into receiving water bodies will cause serious Non-Point Source Pollution (NPS) [3]. The extreme rainfall resulted from climate change can further aggravate the challenges faced by urban drainage systems [4]. Climate change and rapid urbanization are the two most influential uncertainties in urban water system management [2,5,6].
Green infrastructure and gray infrastructure are considered as the main measures to control NPS [7,8]. Traditional gray infrastructure, such as storage tanks, can intercept polluted runoff into receiving water bodies at the terminal, which can effectively alleviate the runoff pollution through centralized measures [9,10]. For green infrastructure, when the runoff flows over its surface; meanwhile during realizing the runoff volume control, it can absorb, exfiltrate and deposit contaminants from runoff, so it can further treat polluted runoff [11]. Jia et al. found that the mixed green infrastructures were effective in removing pollutants from runoff, especially in removing NH3-N (ammonia nitrogen) (73%), Total Nitrogen (TN) (74%), Total Phosphorus (TP) (95%) and Total Suspended Solids (TSS) (35%) [12]. Maniquiz-Redillas and Kim evaluated six green infrastructure systems for the removal of heavy metals (Zn, Pb, Cu, Ni, Cr, Cd, and Fe) from polluted runoff and found that all green infrastructure systems were effective in removing particulate-bound metals (Fe, ZN, Pb) while dissolved metals (Cr, Ni, Cu, Cd) were more difficult to remove with green infrastructure systems [13]. Wilson et al. compared the water quality performance of the green infrastructure site and the conventional site from a whole perspective and found green infrastructure could significantly improve the water quality performance of the areas [14].
The implementation of green infrastructure has been regarded as a more sustainable stormwater management, such as Sustainable Urban Drainage Systems (SUDs), Sponge City, and Low Impact Development practices (LID), which can significantly control peak flow, runoff volume, and runoff pollution load at the source through decentralized measures [8,[15][16][17]. Although having additional aesthetic and ecological benefits compared to gray infrastructure, many studies have also proven that green infrastructure is more ineffective in controlling high-intensity rainfall of short duration, which makes it become an uneconomic measure correspondingly [18][19][20][21]. Since green and gray infrastructure have their own advantages, more and more studies have begun to focus on studying the combined green-gray infrastructure from a compound perspective [22][23][24][25].
The primary and widespread concern in urban drainage systems planning is how to realize the most cost-effectiveness, and Multi-Objective Optimization (MOO) can be used to obtain Pareto solution sets, which can be more intuitive and easier for decision makers to choose optimal planning scheme [26][27][28]. Therefore, a few studies have proposed frameworks that coupled MOO algorithms with hydrological models for cost-effectiveness optimization of drainage systems. Men et al. proposed a mathematical optimization method for LID practices layout and Sponge City planning by using coupled Preference-Inspired Co-Evolutionary Algorithm using goal vectors (PICEA-g) and SWMM [29]. Eckart et al. developed an integrated framework by coupling Borg Multi-objective Evolutionary Algorithm (Borg MOEA) with SWMM to simulate the optimal LID strategies for stormwater control [30]. Qiu et al. used NSGA-II built into The System for Urban Stormwater Treatment and Analysis Integration (SUSTAIN) model to optimize and evaluate various combined LID measures for NPS control [31]. However, most of the previous studies only considered optimizing green infrastructure rather than gray infrastructure or green-gray infrastructure, and a few other studies considered optimizing gray infrastructure, but only took pipe networks as the goal, and pay more attention to water quantity rather than water quality. Optimizing pipe networks has little practical value, for retrofitting pipe networks is extremely hard to realize in the built areas [32]. In addition, green infrastructure can provide significant indirect additional benefits while meeting direct goals [26], however, the results of many studies that integrated green and gray infrastructure to optimize was incomplete because these additional benefits were neglected [25,27,33]. Moreover, the optimization method of most studies were based on the unchanged condition and merely considered initial construction cost as a cost function [34][35][36]. However, it is necessary to recognize the impact of future uncertainties on optimization results and to consider the operation and maintenance (O&M) cost since both green and gray infrastructure have a long lifespan.
Consequently, the objectives of this study are to: (1) couple SWMM and NSGA-II to propose a simulation-optimization framework; (2) estimate the impact of various additional benefit levels of green infrastructure on optimization results of green-gray infrastructure; (3) optimize green, gray and green-gray infrastructures under future uncertainties.
A pilot site of Sponge City in Nanchang, China, was chosen as an example to achieve the above goals. Future climate change and urbanization scenarios were predicted by using General Circulation Models (GCMs), Change Factor Methodology (CFM), and district development plan, respectively.

Site Informations
The study area is located in the south of Nanchang City, with an area of 152.2 ha. This site has a subtropical monsoon climate with an annual average rainfall of 1700 mm, and the rainy season is from May to July during each year, accounting for more than 50% of the total annual precipitation. In this site, polluted and untreated stormwater runoff has been discharged directly into poorly flowing man-made open channels through municipal separate storm sewer systems, and NPS poses a serious threat to the water bodies. Land use types in 2018 were extracted from GIS by referencing satellite images, and land use types in 2030 were simulated based on the Master Plan of Nanchang Central District (2011-2030) formulated by authorities. The main land use types in 2018 were industrial land and green space since the site is a developing economic development zone. In 2030, a number of green spaces will be turned into industrial land ( Figure 1). Table 1 lists the specific information of land use types in 2018 and 2030.  To set up the SWMM model is the pre-processing step for optimization. The site was divided into 49 sub-catchments, 59 conduits, 54 nodes and 3 outfalls according to the extracted satellite images and the direction of the pipe networks (as shown in Figure 2). Horton infiltration was chosen for the hydrology module of SWMM, and the hydraulic module of SWMM was dynamic wave. TSS is one of the most representative pollutants in evaluating the water quality of runoff [37,38], and Assessment Standard for Sponge City Effect also requires to assess the water quality control effect by TSS removal rate [39], therefore, the pollutant in the water quality module of SWMM was set to be TSS, and the land uses in the water quality module were divided into roof, road/other impervious area and green space based on the surface features.

Design Rainfall
The design rainfall event with the duration of 2 h and the return period of 2 years based on the Chicago approach was adopted in this study. The latest storm intensity equation in Nanchang is shown as below: where i is rainfall intensity (mm/min); P is return period (year); t is rainfall duration (min).

Goodness-of-Fit Test
SWMM has many parameters that are difficult to be determined directly, such as Max/Min infiltration rate, Manning's coefficient for impervious/pervious area and conduits, depth of depression storage on impervious/pervious area, etc., in the hydrology module, buildup and washoff functions in the water quality module. Therefore, the parameters that are hard to be determined directly were first set referring to the SWMM Manual and past studies [25], and then, those were revised manually according to the past realistic rainfall events. Nash-Sutcliffe efficiency coefficient (NSE) is a commonly used method to test the goodness-of-fit between observed values and simulated results of hydrological models and was chosen to evaluate the accuracy of the model in this study. The NSE should be greater than 0.5, and the closer to 1 the NSE is, the more credible the model is.

Green Infrastructure and Gray Infrastructure
According to the land use types and surface condition in the study area, green roof, permeable pavement and vegetative swale were selected as green infrastructures to control the runoff from the source. Green roof, permeable pavement and vegetative swale were only placed on the roof, road and green space, respectively. Storage tanks were selected as the gray infrastructures and were set at the front of the three outfalls to intercept polluted runoff during the rainfall event at the terminal until the storage tank was full, and then the polluted runoff would discharge into the channel directly. The green-gray infrastructures scheme is shown in Figure 3.

Decision Variables and Optimization Objectives
Changing the scale of infrastructures is the most significant factor to affect cost and performance compared to changing other parameters [18]. Thus, in order to make the optimization feasible and save optimization time, the decision variables for green infrastructure were construction area (x 1 , x 2 , x 3 ) (see Figure 3) for each LID practice on a given land use type of each sub-catchment, and the other parameters were set as constant. The decision variables for gray infrastructure were construction area (x 4 ) (see Figure 3) of storage tanks, and the height of storage tanks was kept as constant either.
We chose Life Cycle Cost (LCC) as an indicator to evaluate the economic benefits of green and gray infrastructure in this study. LCC includes construction cost and future O&M cost, and the future O&M cost is usually calculated by Discounted Cashflow Model (DCF), which can discount the expected future cost to the present cost and has been used in many previous studies [33,40,41]. Additionally, green infrastructure will produce additional benefits in its life cycle compared with gray infrastructure, such as ecological, aesthetic and environmental benefits [8], and this is why it has become the more recommended stormwater management at this stage. DCF was also used to calculate these future Value of Additional Benefits (VAB) of green infrastructure in this study. However, it is beyond the scope of this study to accurately assess the VAB. Therefore, this study used the α coefficient to represent the discounted values of the expected additional benefits in the future, and it can be considered as a positive cost to discount the negative LCC [36]. When the α coefficient is higher, there are more additional benefits that are taken into account.
How to obtain a higher NPS control rate in a lower cost as much as possible is a primary task in this study. Therefore, the LCC of the green and gray infrastructure and the TSS reduction rate were set to be minimizing and maximizing objectives, respectively.
The objective functions are shown as below: where C green and C gray are the construction costs of green infrastructures and gray infrastructures, respectively; V green and V gray are the present discounted value of the O&M cost of green infrastructures and gray infrastructures, respectively; A i j is the area of jth land use in ith sub-catchment; P jk is the ratio of kth green infrastructure to jth land use; N k is the unit construction cost of kth green infrastructure; X g is the construction area of gth storage tank; H is the given height of storage tanks; G is the unit construction cost of storage tanks; r is the ratio of Annual O&M cost to the construction cost; v is the discount rate; T green and T gray are the life cycle of green and gray infrastructures; n is the number of subcatchments, n = 49; m is the types of land use, m = 3; o is the types of LID practices, o = 3; α is the coefficient representing the VAB of green infrastructure; TSS and TSS extracted from SWMM by using Python are the TSS load discharged from the three outfalls before and after constructed green and gray infrastructures, respectively. The cost information is shown in Table 2. The cost information of green and gray infrastructures was based on Technical Guide for Sponge City Construction and realistic projects in cities of the same economic level as Nanchang.

Future Scenario
Since both green and gray infrastructure have a long lifespan, it should not be neglected that urbanization and climate change in the future will have effects on its performance. Therefore, we designed four scenarios (C0U0, C30U0, C0U30, C30U30) in this study. The baseline scenario (C0U0) is land use in 2018, and the urbanization scenario (C0U30) is land use in 2030 according to the regional development plan. By comparing the differences between C0U0 and C0U30, we can conclude that green space will decrease from 35.38% (in 2018) to 19.67% (in 2030), while industrial land and road will increase from 32.18% (in 2018) to 45.64% (in 2018) and from 17.26% (in 2030) to 19.51% (in 2030), respectively (See Figure 1b). Therefore, the impervious ratio and land use type of the subcatchments in the SWMM model were changed correspondingly to simulate the urbanization scenario. Climate change scenarios (C0U0 and C30U0) were obtained by integrated Global Circulation Models (GCMs) and Change Factor Methodology (CFM), which has been validated and applied in the past studies [25,42]. MarkSim is a weather generator that emulates downscaling results from 17 GCMs for the four Representative Concentration Pathway (RCP) scenarios, including RCP 2.6, RCP 4.5, RCP 6.0, and RCP 8.5 [43]. It can create a random series of local daily precipitation in the future through given specific longitude and latitude coordinates. In this study, all of the GCMs were applied to obtain the ensemble average, and the RCP 6.0 scenario was selected because it is highly consistent with the current tendency and is the most likely to happen [44]. The process that uses CFM to estimate the local scaled future rainfall event is explained as follows [45]: The first step is to use the GCMs weather generator to obtain daily precipitations in 2015 (baseline) and in 2030 (future); Step 2 is to calculate the change factor (Equation (6)); Step 3 is to calculate the local scaled future scenario by applying change factors (Equation (7)); where CF is change factor; GCMs b and GCMs f are the values of annual average precipitation from GCMs baseline and GCMs future climate scenario, respectively; b LS and f LS are the local scaled baseline and future values, respectively, and b LS is the design rainfall event by the storm intensity formula revised in 2015 based on the Chicago approach. The change factor was 1.16 in 2030.

Simulation-Optimization Framework
NSGA-II was used in this study because it is a stable, fast and accurate MOO algorithm [46], and it was widely used in LID practices layout optimization [29]. SWMM lacks in an optimization module, which makes it difficult to recognize the optimal solution from the many simulation results. In this simulation-optimization framework, we used Python to interfere the parameters and extract the simulation results from outputs, and the process was embedded in the workflow of NSGA-II. To make the solutions distributed uniformly and save computing time, the parameters of NSGA-II were set as follows: crossover probability = 0.95; mutation probability = 0.05; number of iterations = 300, and number of populations = 80. The simulation-optimization framework is shown in Figure 4.

Calibration and Validation
In this study, the calibration and validation for the model of hydrology module were implemented by comparing the simulated and observed runoff volume at the outfall during the rainfall event, while those of water quality module were implemented by comparing the simulated and observed TSS concentrations at the outfall during the rainfall events.
The rainfall event on 13 July 2019 was used to calibrate the model, and the rainfall event on 4 July 2019 was used to validate the results. Figure 5 shows the comparison of the simulated values and the observed values. The NSE values for calibration of the hydrology and water quality modules of the model were 0.89 and 0.82, respectively, while those of validation were 0.69 and 0.65, respectively. The NSE values of runoff volume for both calibration and validation are greater than those of TSS concentration, and the accuracy of the model for runoff volume simulation is greater than that for TSS concentration. In summary, the calibration and validation values for both hydrology module and water quality module are greater than 0.5, and it indicates that the accuracy of the model is acceptable, which can be applied to further following studies. The results of calibration for hydrology and water quality parameters of the model are shown in Table 3andTable 4.

Effect of VAB of Green Infrastructure on Optimizing
In this section, we mainly discussed the effects of different α coefficients on optimizing. The results are shown in Figures 6-8. The green and the gray in the figures (Figures 6, 7 and 8a) are optimizing green infrastructure and gray infrastructure separately, while the green-gray is synchronously optimizing both green and gray infrastructure, and the color scale at the right of the graph is the ratio of gray infrastructure cost to total cost. Each point in Pareto solution sets obtained by multi-objective optimization is the optimal solution, which has the best benefit at a certain cost, and the solutions in the sets also contain a lot of information, such as the area of the green infrastructure in each sub-catchment and the volume of the gray infrastructure, which is helpful for decision makers to choose the most suitable solution intuitively based on water quality objectives and cost constraint. solution, which has the best benefit at a certain cost, and the solutions in the sets also contain a lot of information, such as the area of the green infrastructure in each sub-catchment and the volume of the gray infrastructure, which is helpful for decision makers to choose the most suitable solution intuitively based on water quality objectives and cost constraint. The optimization results when α = 1 are shown in Figure 6. The cost-effectiveness curve of gray infrastructure is significantly different from green infrastructure (Figure 6a). The algorithm will not identify the additional benefits of green infrastructure when α = 1, and gray infrastructure is a more cost-effective stormwater management measure at this point. When synchronously optimizing both green and gray infrastructure, the algorithm would select the gray infrastructure almost exclusively in the range of the solution set interval. Therefore the gray infrastructure accounts for most of the cost weight of the green-gray infrastructure (Figure 6b). Due to the limited treatment capacity of gray infrastructure, its marginal benefit decreases gradually with the increasing cost, and the green infrastructure would be present slightly at TSS reduction rate above 60%. Although green infrastructure is a more promoted stormwater management measure, it is an inefficient measure for NPS control compared to gray infrastructure when its VAB were not considered. It means that gray infrastructure should be given priority as an NPS control measure when the investment is constrained.
When α = 3, the optimization results are shown in Figure 7. At this point, the green infrastructure shows a comparable performance to the gray infrastructure (Figure 7a). Since green and gray infrastructure treatment capacity is limited, their solution set curves tend to flatten gradually with the increasing cost, which is more pronounced for green infrastructure, and the finding is very similar to Leng et al. and Qiu et al. [25,31]. This means that the marginal return decreases gradually with increasing investment. The optimization results when α = 1 are shown in Figure 6. The cost-effectiveness curve of gray infrastructure is significantly different from green infrastructure (Figure 6a). The algorithm will not identify the additional benefits of green infrastructure when α = 1, and gray infrastructure is a more cost-effective stormwater management measure at this point. When synchronously optimizing both green and gray infrastructure, the algorithm would select the gray infrastructure almost exclusively in the range of the solution set interval. Therefore the gray infrastructure accounts for most of the cost weight of the green-gray infrastructure (Figure 6b). Due to the limited treatment capacity of gray infrastructure, its marginal benefit decreases gradually with the increasing cost, and the green infrastructure would be present slightly at TSS reduction rate above 60%. Although green infrastructure is a more promoted stormwater management measure, it is an inefficient measure for NPS control compared to gray infrastructure when its VAB were not considered. It means that gray infrastructure should be given priority as an NPS control measure when the investment is constrained.
When α = 3, the optimization results are shown in Figure 7. At this point, the green infrastructure shows a comparable performance to the gray infrastructure (Figure 7a).
Since green and gray infrastructure treatment capacity is limited, their solution set curves tend to flatten gradually with the increasing cost, which is more pronounced for green infrastructure, and the finding is very similar to Leng et al. and Qiu et al. [25,31]. This means that the marginal return decreases gradually with increasing investment. It can favor decision makers to choose a reasonable range of investment and thus to avoid inefficiencies. In optimizing process, the green and gray infrastructure showed a close probability of being selected by the algorithm. For example, the cost share of gray infrastructure is approximately 0.36 when the TSS reduction rate is 50% (Figure 7b). However, the cost share of gray infrastructure will gradually increase with the increasing cost. This may be due to the weaker performance attenuation effect of gray infrastructure caused by increasing cost compared to green infrastructure so that the gray infrastructure accounts for an increasing proportion when the cost increases. In addition, the synchronous optimizing green-gray infrastructure can slow down the tendency of the curve to flatten compared to the single infrastructure, which makes the solution set curve be more linear. In other words, if green and gray infrastructure have similar cost-effectiveness, green-gray infrastructure will show the co-benefit, delaying the tendency of their marginal benefit to decline. This result is similar to Zhou et al. [36]. However, it is worth noting that the co-benefit is not obvious when the cost is not too high. For example, the TSS reduction rate is 50.6%, 48.3% and 50.8% for green, gray and green-gray infrastructure when the cost is USD 3 million, respectively, while the TSS reduction rate is 56.7%, 55.5% and 58.6% for green, gray and green-gray infrastructure when the cost is USD3.5 million, respectively (Figure 7). The co-benefit increases the possibility that watersheds reach the maximum NPS control rate and will gradually appear with the increasing cost.
The optimization results when α = 5 are shown in Figure 8. At this point, there are enough VAB to be considered, and obviously, green infrastructure is a more cost-effective NPS control measure than gray infrastructure (Figure 8a). Gray infrastructure would not be selected by the algorithm for green-gray infrastructure when the TSS reduction rate is below 51% (Figure 8b). However, since green infrastructure has higher attenuation performance caused by increasing cost, consistent with the discussion above, a better solution at this point as the green infrastructure is that the gray infrastructure would still be considered when the TSS reduction rate is above 51%. In summary, when α = 5, decision makers could choose green infrastructure as a principal stormwater management measure. Only when a higher water quality objective needs to be reached should a supplement of a certain scale of gray infrastructure on the basis of constructed green infrastructure be considered. It can favor decision makers to choose a reasonable range of investment and thus to avoid inefficiencies. In optimizing process, the green and gray infrastructure showed a close probability of being selected by the algorithm. For example, the cost share of gray infrastructure is approximately 0.36 when the TSS reduction rate is 50% (Figure 7b). However, the cost share of gray infrastructure will gradually increase with the increasing cost. This may be due to the weaker performance attenuation effect of gray infrastructure caused by increasing cost compared to green infrastructure so that the gray infrastructure accounts for an increasing proportion when the cost increases. In addition, the synchronous optimizing green-gray infrastructure can slow down the tendency of the curve to flatten compared to the single infrastructure, which makes the solution set curve be more linear. In other words, if green and gray infrastructure have similar cost-effectiveness, green-gray infrastructure will show the co-benefit, delaying the tendency of their marginal benefit to decline. This result is similar to Zhou et al. [36]. However, it is worth noting that the co-benefit is not obvious when the cost is not too high. For example, the TSS reduction rate is 50.6%, 48.3% and 50.8% for green, gray and green-gray infrastructure when the cost is USD 3 million, respectively, while the TSS reduction rate is 56.7%, 55.5% and 58.6% for green, gray and green-gray infrastructure when the cost is USD3.5 million, respectively (Figure 7). The co-benefit increases the possibility that watersheds reach the maximum NPS control rate and will gradually appear with the increasing cost.
The optimization results when α = 5 are shown in Figure 8. At this point, there are enough VAB to be considered, and obviously, green infrastructure is a more cost-effective NPS control measure than gray infrastructure (Figure 8a). Gray infrastructure would not be selected by the algorithm for green-gray infrastructure when the TSS reduction rate is below 51% (Figure 8b). However, since green infrastructure has higher attenuation performance caused by increasing cost, consistent with the discussion above, a better solution at this point as the green infrastructure is that the gray infrastructure would still be considered when the TSS reduction rate is above 51%. In summary, when α = 5, decision makers could choose green infrastructure as a principal stormwater management measure. Only when a higher water quality objective needs to be reached should a supplement of a certain scale of gray infrastructure on the basis of constructed green infrastructure be considered.

Optimizing Green and Gray Infrastructure under Future Uncertainties
In this section, we chose α = 5 as the VAB of the green infrastructure for the following analysis of green and gray infrastructure under future uncertainties, for green infrastructure has a comparable cost-effectiveness to gray infrastructure at this additional benefits level as discussed in the last section.

Impact of Water Quality under Future Uncertainties
The TSS load discharged from the outfalls under different scenarios are shown in Table 5. Compared to C0U0, C30U0, C0U30 and C30U30, the TSS load increases by 2.71%, 5.26% and 9.92%, respectively. The change of impervious ratio and the increased pollutant load on the surface under urbanization scenario have a greater impact on TSS load than the higher wash-off effect caused by the increased rainfall intensity under climate change scenario. In addition, urbanization and climate change have a synergistic rather than a linear effect on TSS load, the finding is similar to previous studies [1]. This means that the change in land use and climate in the future poses a serious threat to the water ecosystem by producing more NPS.  Figure 9 shows the optimization results under the four scenarios. By observing the Pareto solution sets, we found that both green and gray infrastructure are more sensitive to climate change than urbanization (see Figure 9a  For example, at USD 4 million, the TSS reduction rates for green and gray infrastructure are approximately 59.6% and 59.8% under C0U0, respectively, but 56.7% and 53.8% under C30U0, respectively. This is due to the green and gray infrastructure which is mainly relied on to reduce the runoff volume and thus to control the TSS load. Under C30U0, the green and gray infrastructure were filled up earlier by the increased runoff volume generated from the higher intensity rainfall, and those have a greater impact on storage tanks.
According to the optimization results under C0U30, although there is more runoff generated from increasing ratio of impervious surface under C0U30, the concentration of TSS in polluted runoff will also be higher with increasing TSS load on the surface caused by the change of land use, so the green and gray infrastructure can remove more TSS load at a certain runoff volume reduction rate compared to those under C30U0. In addition, it is worth noting that the Pareto solution curve of green infrastructure under C0U30 is almost unchanged. Besides the reasons discussed above, another reason is the phenomenon that green infrastructure was placed targetedly on critical areas with high TSS load resulting from urbanization. In other words, in optimizing process, the algorithm reallocated cost and placed green infrastructure at the source, that is, new critical areas with high TSS load generated from the change of land use type. However, gray infrastructure cannot be placed on the new areas caused by urbanization with high TSS load. The finding elucidated that gray infrastructure is more sensitive under urbanization compared to green infrastructure. In practice, storage tanks are a huge project and difficult to continue expanded with the urbanization process, while LIDs are a small project placed at source, and decision makers can continually adjust their planning schemes to suit the process of urbanization. Figure 9c depicts the variation in ratio of gray infrastructure cost in green-gray synchronous optimization results. For example, when the cost is USD 3 million, the ratio of gray infrastructure cost is 0.39, 0.35, 0.33 and 0.27 under C0U0, C30U0, C0U30 and C30U30, respectively. Compared to C0U0, gray infrastructure is less selected in other future scenarios. When the cost is below USD 2.93 million, the curve C0U30 and the curve C0U0 almost overlap, and green infrastructure is more selected compared to other future scenarios, while the curve C0U30 cannot maintain the same increasing trend as the curve C0U0 with the increasing cost. This can be explained by the differences in the sensitivity of green and gray infrastructure under future uncertainties which can affect the synchronous optimization results of green and gray infrastructure. In other words, the same construction planning of green-gray infrastructure under C0U0 is not the optimal solution under other future scenarios. The finding indicates that to adapt the urbanization and climate change by adjusting the combination strategy of green and gray infrastructure can obtain a higher resilience in future uncertainties.
As a whole, the green-gray infrastructure shows comparable resilience in future uncertainties with the green infrastructure. For example, when the cost is USD 4 million, the TSS reduction rates for green and green-gray infrastructure are approximately 59.6% and 63.2% under C0U0, respectively, but 55.1% and 58.4% under C30U30. However, as discussed in the last section, the marginal benefit of green-gray infrastructure decreases more slowly with increasing cost than green infrastructure only. Therefore, green-gray infrastructure should be considered mainly if green infrastructure has a similar resilience to green-gray infrastructure under future uncertainties. Figure 10 shows the costs of maintaining the water quality levels to C0U0 under C30U30 (based on the optimization results). In the X-axis, −9.91% is the TSS load under C30U30 (increase of 9.91% compared to C0U0), and 0, 10, 30 and 55 are the same level as C0U0 and reductions by 10%, 30% and 55% compared to those of C0U0, respectively. The first values of cost are USD 0, meaning there are no control measures to be implemented under C30U30. The second values of cost mean that it will cost USD 0.635, USD 0.617 and USD 0.614 million to maintain the same water quality of C0U0 under C30U30 for gray, green and green-gray infrastructure, respectively, and at this level of water quality objective, the co-benefit of the green-gray infrastructure are less pronounced. The third values of cost are USD 1.211, USD 1.172 and USD 1.165 million for gray, green and green-gray infrastructure, respectively. The fourth values of cost are USD 2.482, USD 2.295 and USD 2.278 million for gray, green and green-gray infrastructure, respectively. The fifth values of cost are USD 4.730, USD 4.902 and USD 4.332 million for gray, green, and green-gray infrastructure, respectively, and at this level of water quality objective, gray infrastructure is a more cost-effective measure than green infrastructure. In addition, the co-benefit of green-gray infrastructure is relatively obvious. Figure 10 indicates that as the TSS load decreases, the cost increases markedly. However, the increase of cost with TSS load decreasing is various for green, gray, and green-gray infrastructure, and the cost-effectiveness of green-gray infrastructure is higher than that of green or gray infrastructure only. Therefore, when the low control objectives need to be reached, green infrastructure should be considered as the main measure, while gray infrastructure should be the main measure when the control objectives are relatively high.

Conclusions
In this study, a simulation-optimization framework was developed by integrated NSGA-II and SWMM to optimize green and gray infrastructure for NPS control under climate change and urbanization scenarios. Considering the VAB of green infrastructure, different values of α coefficients were applied to assess the impact on optimization results. The following conclusions were obtained:

1.
Gray infrastructure is a more cost-effective measure for NPS control, unless there are sufficient VAB of green infrastructure to be considered. When a high NPS control rate is reached, the marginal benefit of green infrastructure decreases faster than that of gray infrastructure, so when the water quality objective is relatively high, it is a less economical measure to adopt only green infrastructure for NPS control.

2.
When there are enough VAB (i.e., α = 3) for the green infrastructure, green-gray infrastructure will show co-benefits compared to green or gray infrastructure only, and the co-benefit will be more obvious with the increasing cost.

3.
Green and green-gray infrastructure have a greater resilience under future uncertainties, especially under urbanization scenario, while gray infrastructure is more sensitive under climate change scenario.
The results of this study can serve as a reference for the cost-effectiveness issues involved in watershed NPS control. However, the VAB of green infrastructure was not accurately determined in this study, and it is a relatively simplified approach to represent the VAB green infrastructure by α coefficients. In addition, the gray infrastructure optimization scheme is also a relative simplification, and only the volume of storage tanks is considered in optimizing. Changing the location of the storage tanks may have an impact on its performance. Therefore, in the further following exploration, a method to accurately evaluate the VAB of green infrastructure needs to be established, and both location and volume of storage tanks should be optimized.