Process-based modeling for ecosystem service provisioning: Non-linear responses to restoration efforts in a quarry lake under climate change

. In this study, we propose an ES modeling tool that integrates lake physics, ecology and service provisioning into a holistic modeling framework. We applied this model to a Dutch quarry lake, to evaluate how nine ESs respond to technological-based (phosphorus (P) reduction) and nature-based measures (wetland restoration). As climate change might be affecting the future effectiveness of restoration efforts, we also studied the climate change impacts on the outcome of restoration measures and provisioning of ESs, using climate scenarios for the Netherlands in 2050. Our results indicate that both phosphorus reduction and wetland restoration mitigated eutrophication symptoms, resulting in increased oxygen concentrations and water transparency, and decreased phytoplankton biomass. Delivery of most ESs was improved, including swimming, P retention, and macrophyte habitat, whereas the ES provisioning that required a more productive system was impaired (sport fishing and bird watching). However, our modeling results suggested hampered effectiveness of restoration measures upon exposure to future climate conditions, which may require intensification of restoration efforts in the future to meet restoration targets. Importantly, ESs provisioning showed non-linear responses to increasing intensity of restoration measures, indicating that effectiveness of restoration measures does not necessarily increase proportionally. In conclusion, the ecosystem service modeling framework proposed in this study, provides a holistic evaluation of lake restoration measures on ecosystem services provisioning, and can contribute to development of climate-robust management strategies.


Introduction
We have entered a human-dominated geological epoch, coined the Anthropocene (Lewis and Maslin, 2015), characterized by an increasing impact and over-utilization of ecosystems by humans.In the Anthropocene, the demand for almost all ecosystem services is on the rise (Millennium Ecosystem Assessment, 2005).Ecosystem services (ESs) are defined as direct and indirect contributions of ecosystems to human well-being (Carpenter et al., 2009).In recent years, this concept of ecosystem service has guided ecosystem management and restoration efforts, aiming to integrate social, economic and ecological perspectives (Kull et al., 2015;Seppelt et al., 2011;Valencia Torres et al., 2021).
Quantifying ecosystem services can be instrumental in recognizing the benefits humans receive from ecosystems, providing stronger arguments for ecological restoration (Grizzetti et al., 2019;Guerry et al., 2015).Conveying restoration impacts in terms of the loss or gain of ESs can facilitate effective communication of restoration outcomes to policy-makers and river basin authorities responsible for implementing restoration measures (Wortley et al., 2013).While modeling terrestrial ecosystem services often focuses on mapping ESs provisioning through spatial variations of catchment attributes (e.g., land use, topography, lithology) (Nelson and Daily, 2010), the dynamics of water quantity and quality necessitate a more explicit consideration in aquatic ecosystem service modeling (Grizzetti et al., 2016).
There is increasing evidence that freshwater ecosystem services provisioning is closely linked to the ecological quality (or ecological state) of different aquatic environments, including shallow lakes (Janssen et al., 2021), deep lakes (Seelen et al., 2021), rivers and coastal waters (Grizzetti et al., 2019).Based on data reported under the European Water Framework Directive, Grizzetti et al. (2019) demonstrated that higher provisioning of ESs is mostly correlated with better ecological states, particularly for regulating services (e.g., water purification, erosion retention, flood protection) and cultural services (e.g., recreation).However, current modeling tools for water-related services primarily focus on water quantity (Grizzetti et al., 2016), with limited integration of services closely related to water quality (Keeler et al., 2012).Water quality dynamics are mediated by complex interactions among a myriad of ecosystem processes, which are often oversimplified in large-scale modeling frameworks.For instance, one widely-used ecosystem service model, InVEST, simplifies by using nutrient loading as a proxy for determining the availability of lake-related ESs (Nelson et al., 2009;Polasky et al., 2011), assuming simple linear responses of ecosystems to nutrient loading.This approach contradicts the resistance theory of ecosystems (Gómez-Baggethun et al., 2011;Ibelings et al., 2007), which supports threshold-type ecosystem responses to pressures.Consequently, the assessment of management actions often relies on variables collected at the landscape scale (e.g., Burkhard et al., 2012;Hernández-Romero et al., 2022), which may be inaccurate due to the aforementioned nonlinear responses or ill-fitting when assessing the impacts of in-lake restoration measures (Lürling and Mucci, 2020).Keeler et al. (2012) proposed a conceptual framework linking ecological-related services with corresponding water quality variables based on a review of existing ES models, emphasizing the importance of this link in assessing management actions.Nevertheless, an integrated ES modeling framework closely tied to water quality dynamics is still lacking.
Successful lake management efforts should also account for the effects of climate change in addition to eutrophication control (Salk et al., 2022).Previous studies have suggested nonlinear ecosystem responses to climate change pressures (Burkett et al., 2005).Climate change can influence lake physics and ecology, e.g., by changing stratification patterns and mineralization rates.This in turn may lead to water quality deterioration (Moss et al., 2011) and impair the effectiveness of restoration measures (Cabrerizo et al., 2020;Zhan et al., 2021Zhan et al., , 2022)), further affecting aquatic ES provisioning.It remains an urgent question for water managers to determine how robust their current measures for restoring ES provisioning are in light of future climate conditions.Addressing this question requires a model capable of simulating climate change impacts on lake physics (e.g., thermal regimes) and ecological processes.
To address the challenges posed by various pressures on lakes, such as eutrophication and climate change, different restoration strategies have been developed and can be divided into two groups: technological solutions and nature-based solutions.Technological restoration strategies focus on employing engineered solutions to mitigate stressors affecting lakes.One prominent technological approach involves the controlled removal or inactivation of excess phosphorus from water bodies and sediments using techniques like geoengineering, or dredging (Zhan et al., 2022).Phosphorus reduction strategies are amenable to implementation in diverse lake types and meet the increasing demands for spaces of different types of land-uses, as technological-based measures often require relatively less space and act at shorter time scales (Kim et al., 2021).In contrast, nature-based solutions (NbS), exemplified by wetland restoration, emphasize the use of natural features and processes to restore lakes.Importantly, NbS need to provide benefits for human well-being and biodiversity.Wetlands are a hotspot for biodiversity and act as natural filters that absorb excess nutrients, sediment, and contaminants from water bodies (Sollie et al., 2008).There is also a growing demand for nature-based solutions, as they are suggested to be potential adaptation strategies to climate change (Seddon et al., 2020;van Leeuwen et al., 2021).Yet establishment of a dose-effect relationship is needed for the evidence-based implementation of these restoration measures (Seddon et al., 2020).The main objective of this paper is to assess how effective our current lake restoration measures are in restoring ecosystem services provisioning, for which we followed this work flow: 1).Develop a model framework linking lake physics, ecology and ecosystem services provisioning.2).Assess effectiveness of a nature-based restoration (wetland restoration) and a technological restoration measure (phosphorus reduction) in a eutrophic quarry lake.
3).Study the impacts of climate change on ES delivery, by running the model under two climate scenarios.

An integrative ecosystem service modeling framework
A freshwater one-dimensional lake physics model -FLake (Kirillin et al., 2011) was used to predict the vertical temperature profile in the water column, based on weather data or climate scenarios.FLake uses climate data (air temperature, wind speed, solar radiance, humidity and cloud cover) as input.The output of FLake including mixing and water temperatures is subsequently imposed into a lake ecosystem model (PCLake+) as boundary conditions reflecting climate-related impacts (Fig. 1).
PCLake is a process-based ecological model that was developed to simulate water quality and assess the trophic state of lakes based on ecological interactions (Janse, 2005).It models nutrient cycling including nitrogen and phosphorus and a simple food web consisting of three functional groups of phytoplankton (cyanobacteria, green algae and diatoms), zooplankton, and fish.Moreover, it involves a wetland module that exhibits a purification function of the lake water, through biogeochemical processes explained in Janse et al. (2001).PCLake+ is an expanded version of PCLake that allows for water column stratification to take place, to model deep lakes (Janssen et al., 2019).After that, an extension was done by Chang et al. (2020) to include more realistic cyanobacterial traits, which modelled surface cyanobacterial biomass accumulation in addition to the epilimnion and hypolimnion biomass.In this study, we further expanded PCLake+ with a module for ecosystem services provisioning.The most up-to-date PCLake+ version can be found at GitHub: https://github.com/pcmodel/PCModel.git.

Modeling of ecosystem services provisioning
We built an ecosystem service module into PCLake + to translate the water quality into a quantitative description of ecosystem services provisioning.We followed a framework for assessing ecosystem services proposed by Seelen et al. (2021), which links ecosystem state indicators with ecosystem service provisioning through a threshold approach.The threshold values reflect the values that certain water quality parameter in a lake must attain to support the provision of a given service.The threshold values were based on published peer-reviewed literature, a field campaign covering 51 quarry lakes in the south of the Netherlands, and expert judgment (Table 1; see Seelen et al., 2021 for the supporting materials).Per ecosystem service, different aspects of the water quality requirements of the service are considered.For instance, the service of swimming is only suitable when the lake has sufficient transparency, the Q.Zhan et al. cyanobacterial biomass is at a safe level, and there is adequate vegetation-free water column.In the ES module, the suitability of delivering each ES was expressed by an indicator function ranging between 0 and 1, with "1" representing a fully suitable provisioning, "0" representing an unsuitable provisioning, and values in between representing a moderate suitability.
In total, nine ESs are modelled with their water quality requirements summarized in Table 1.We followed the Common International Classification of Ecosystem Services (CICES; Haines-Young and Potschin, 2012), in which the ESs are divided into four different groups: provisioning (water, materials, energy and others), regulation and maintenance (remediation and regulation of the biophysical environment, flow regulation, regulation of the physic-chemical and biotic environment), cultural (physical or experiential use of ecosystems, intellectual representations of ecosystems), and abiotic (abiotic materials, energy, and space) services.We selected the ecosystem services that can be provided by quarry lakes following Seelen et al. (2022).Our final selection was constrained by the capacity of PCLake + to compute quantitative ecological state indicators.Some adjustments have been made to the ecosystem model that are described in detail in SI section 1.A brief summary of the main modifications: 1).Hypoxia inhibition effect on fish growth is included into PCLake+; 2).We ran our model using hourly meteorological input data to capture the diurnal variations in light, temperature, evaporation, and rainfall, as drivers of water quality dynamics.

Model application and validation
Our developed ES modeling framework was applied to a quarry lake located in the south of Netherlands -"Put aan de Omloop" (51 • 79′22.8″N, 4 • 95′15.2″E)."Put aan de Omloop" or Lake de Omloop is a typical quarry lake that was created as a result of sand mining  activities, showing up on topographic maps since 1969 (see https://www.topotijdreis.nl/).After creation, the lake has undergone a land-use induced eutrophication process, resulting in an increasing frequency of algae blooms.The lake is characterized by a relatively small surface area (59,370 m 2 ) with an average water depth of 7.7 m (see Figs. SI-1 for lake bathymetry).The lake sediment is identified as sand-type soil, covered with a fluffy layer of organic matter resulting from an accumulation of detritus originating from primary production.The water and nutrient budgets are summarized in Tables SI-1, with the detailed methodology on estimation of each source described in SI section 4. In short, this lake is isolated from any surface flows and fed by both groundwater and precipitation (respectively 47.8% and 52.2%), resulting in a long residence time of ca.1724 days.The background eutrophic level is high attributed to fertilizer application in the surrounding agriculture area, with especially high nitrogen loading (≈0.018 g N/m 2 /day), and phosphorus loading (≈0.45 mg P/m 2 /day).
Given the fact that this lake is mainly regulated by the ground water state which is relatively stable over time within a time window of a year, for simplification, we assumed static water inflows and external nutrient loading as the forcing data for PCLake+.Wind fetch was calculated based on the measured wind speed and direction using an approach introduced in Janssen et al. ( 2017).We reduced the wind speed by half to account for the windbreak effects by the forests surrounding the study site (Jeong and Lee, 2020).
The ecosystem model PCLake + used in this study was calibrated based on a generic lake dataset of primarily Dutch lakes (Janse et al., 2010).PCLake + has a large set of parameters (>250), making overfitting the model a risk when subjecting it to site-specific calibration when data is not abundantly present.Hence, we rely on the generic calibration for our study and only adjust boundary conditions of the lake (i.e., depth, hydraulic and nutrient loads, climate forcing, wind fetch, etc., see Tables SI-1).
We collected meteorological input data from the closest weather station of the Royal Netherlands Meteorological Institute (Herwijnen, https://www.knmidata.nl/)from 2011 as a reference year.FLake simulation suggests that the lake is a dimictic lake undergoing thermal summer stratification (i.e., mixing depth ∕ = lake depth), which is confirmed by the measured depth profiles of water quality variables in the summer of 2014 (see Figs. SI-2), and occasional winter inverse stratification (Fig. 2-a).Under a climate change scenario, the summer stratification was prolonged by more than two weeks, with an earlier onset in spring and a postponed termination in autumn.Winter stratification however becomes much rarer with warmer temperature.A detailed description of the current ecosystem states can be found in SI section 2.
We combined the field measurements from the period from 2003 to 2016 (see appendix B for dataset, combinedly shown in Fig. 2 c-d & Figs.SI-3), to valid the PCLake + performance.A suite of water quality data was measured intermittently with water samples taken in the epilimnion, at ~50 cm below the water surface.We ran PCLake + for a period of 30 years to reach equilibrium states which are no longer dependent on the initial states.The simulation results of the last year Q.Zhan et al. was used for validating if the model is able to capture the generic dynamics of the lake ecosystem.Without any calibration of the default model parameters, PCLake + showed an overall adequate performance in capturing the generic water quality dynamics in the study lake, especially regarding timing of the onset and offset of summer peaks (see Fig. 2, panels c-d).The statistical agreements between the simulation and observation data were evaluated on the coefficient of determination (R 2 ) and Root Mean Square Deviation (RMSE).To exemplify, O 2 : R 2 = 0.36, Root Mean RMSE = 2.04 mg/L; TN: R 2 = 0.40, RMSE = 0.71 mg/L; Tot-Chla: R 2 = 0.44, RMSE = 13.3 μg/L.

Scenario analysis
(1) Restoration scenarios: After validation of the model performances, we simulated different restoration scenarios that can potentially tackle the water quality deterioration in Lake De Omloop.Note that no real restoration has taken place in the lake to date.We evaluated two different types of restoration scenarios with respect to their impacts on the lake ecosystem states and subsequently ecosystem service provisioning.The restoration scenarios include a technology-based measure phosphorus reduction, and a nature-based measurewetland restoration.We undertook a bifurcation analysis of the validated model for the studied quarry lake, to study how the lake ecosystem responds to different levels of restoration intensities.For each level of restoration, the model output variables during summer period (day 150 to day 210) were averaged to represent the response of ecosystem states as well as ecosystem service provisioning.
Technology-based restoration scenarios exhibit variable effectiveness in P reduction (Zhan et al., 2022).We simulated the full spectrum of P reduction effectiveness ranging between 0% and 100% (with 7% increments), to investigate the impact of varying levels of effectiveness on ecological states and subsequently ES provisioning.A P reduction of 0% represents the initial condition without restoration treatment, while a P reduction of 100% represents that the lake is devoid of phosphorus loading.A complete removal of phosphorus loading is not common in practice but theoretically can be achieved by technological approaches, for instance, through using a filtration system with P-binder pumping water thoroughly and repeatedly.Note that we assume technology-based measures solely reduce P without affecting nitrogen (N).This assumption is based on the fact that P is often the primary target element for most chemical adsorbents, though some engineering measures, such as dredging, target both P and N (Zhan et al., 2022).
In addition, we investigate the effectiveness of nature-based restoration scenarios in the form of purifying wetland creation.This was carried out via expansion of the wetland fraction in the model, a parameter representing the size of the wetland area relative to that of the lake area ranging between 0 and 100%.Several ecological processes are present in the wetland module: transport and settling of suspended solids, denitrification, nutrient uptake by marsh vegetation (increasing nutrient retention), and improvement of habitat conditions for predatory fish.The substances and process descriptions (mineralization, settling, P adsorption, nitrification and denitrification) are analogous to those in the lake model, except that the water depth is much lower (default 0.5 m), settling velocities are higher due to the absence of wind action and resuspension is assumed to be zero.Phytoplankton is assumed not to grow in the shadow of the reed vegetation.Mixing between the water columns of the lake and the wetland is described by an exchange coefficient (representing both dispersive transport and transport due to water level changes) multiplied by the concentration difference.We explored the impacts of an increased coverage of wetland area ranging from 10% to 100% (relative to water surface area, with 5% increments).Theoretically, however, this value can go beyond 100%., representing a larger wetland area relative to that of the lake.
To assess the magnitude of restoration effects, we applied linear regression models to the scenario results.For ecosystem state indicators, we calculated a slope-to-intercept ratio (=slope/intercept, in %) as a standardized indicator metric for the magnitude of change in the ecosystem state variables in relation to their initial levels.For instance, a slope-to-intercept ratio of 100% or − 100% indicates that, with 100% increments in restoration intensity, the response variable is increased or decreased by 100%, respectively, compared to their initial states without restoration.As for ecosystem services provisioning that are scaled and ranging between 0 and 1, we convert the slope into percentage as an indicator of the magnitude of the response.We reported the adjusted R-square of the regression model that corrects for the sample size effect (Thompson, 2007).An adjusted R-square value of 1 represented a linear relationship, whereas a value smaller than 1 indicates a non-linear relationship, validated by a visual inspection.All statistical analyses were performed in R language (R Core Team, 2019).
(2) Climate change: To evaluate the impacts of climate change on ecosystem states and subsequent service delivery, we followed the predictions of future climates scenarios for the Netherlands by KNMI (Royal Netherlands Meteorological Institute Ministry of Infrastructure and the Environment; Attema et al., 2014).We implemented the climate scenario changes for the prediction of climate around 2050 under global warming and high changes in air circulation patterns, representing the "most extreme" scenario at that time.This scenario shows an increase in air temperature by 2.3 • C, an increase in precipitation by 5%, and an increase in solar radiation by 1.2%.These future climate conditions were first implemented in FLake for simulation of lake physics under climate change, from which the outputted mixing depth and water temperature vertical profile were used to force PCLake + for prediction of lake ecology.To assess the magnitude of difference in ecosystem responses between two climate scenarios, we calculated an effect size metric using Hedges's g s (Lakens, 2013).Effect size is regarded as a standardized metric which can be understood regardless of the scale of measured variables.For interpretation of effect sizes we followed adopted Hedges's g s thresholds: no evidence (|g s | < 0.2); weak (|g s | < 0.5); moderate (|g s | <0.8) and strong (|g s | ≥ 0.8) (Munthali et al., 2022;Thompson, 2007).We used R package effsize for the effect size calculations (Torchiano, 2016).For visualization of our results, we used a colour-blind-friendly colour palette following Wong (2011).
Most ESs responded non-linearly to the intensity of the restoration measures (see Fig. 5), with the effects on ES being overall more pronounced at intensified interventions.For macrophyte habitat, the restoration effects were marginal until P reduction reached 50% or wetland coverage was larger than 30%.For nitrogen retention, a response was only found after the wetland coverage reaches ca.30%.For phosphorus retention, P reduction showed increased effect with intensified interventions, whereas the wetland restoration effect declined with increased wetland coverage.Similar response patterns were found for swimming service under P reduction, where a positive response emerged only after P reduction of ca.50-70%.

Climate change impacts on ecosystem state indicators and ecosystem services provisioning
Under the 2050 Dutch climate change scenario (see dashed lines in Figs.3-4 for bifurcation analyses, and Fig. 6 for the statistics of Hedges's g s effect size test), eutrophication symptoms were reinforced with declined oxygen concentrations in both epilimnion (For P reduction: strong Hedges's g s ; For wetland restoration: strong Hedges's g s ) and hypolimnion (For P reduction: strong Hedges's g s ; For wetland restoration: strong Hedges's g s ), decreased Secchi disc depth (For P reduction: strong Hedges's g s ; For wetland restoration: strong Hedges's g s ).Though the cyanobacteria chlorophyll-a in the upper layer (surface + epilimnion) showed negligible increase under future climate conditions (For P reduction: no evidence Hedges's g s ; For wetland restoration: no evidence Hedges's g s ), the surface bloom showed a strong response (For P reduction: moderate Hedges's g s ; For wetland restoration: strong Hedges's g s ).Climate change showed a negative effect on macrophyte growth under P reduction scenarios (For P reduction: strong Hedges's g s ), whereas negligible effect was detected upon exposure to wetland restoration (no evidence Hedges's g s ).
The provisioning of almost all ESs was hampered under future climate conditions (see dashed lines in Fig. 5 for bifurcation analyses), with the effect sizes varying among ESs and upon restoration scenarios (Fig. 6).The service 'swimming' was vulnerable to climate change (For P reduction: strong Hedges's g s ; For wetland restoration: weak Hedges's g s ).In contrast, the service 'fishing' showed limited improvement in suitability under 2050 Dutch climate conditions (For P reduction: weak Hedges's g s ; For wetland restoration: weak Hedges's g s ).Relative to P reduction, wetland restoration was able to mitigate climate change impacts on macrophyte-related ESs, which includes macrophyte habitat (For P reduction: moderate Hedges's g s ; For wetland restoration: weak Hedges's g s ), the availability of macrophytes for thatching (For P reduction: strong Hedges's g s ; For wetland restoration: no evidence Hedges's g s ), and particle capture (For P reduction: strong Hedges's g s ; For wetland restoration: no evidence Hedges's g s ).

Discussion
Our modeling framework enabled us to evaluate how the provision of nine ecosystem services in a quarry lake was impacted by a technology-based (phosphorus reduction) and a nature-based (wetland restoration) restoration scenario, under current and future climate scenarios.To this end, we provided an ecosystem service modeling approach that comprehensively incorporates complex ecosystem processes, filling what had been, to the best of our knowledge, a previously existing gap.Our results indicated that both types of restoration scenarios could mitigate eutrophication symptoms.However, the effectiveness of both restoration measures did not linearly increase with the restoration intensity.The ESs that require good water quality were improved, including swimming, irrigation and macrophyte habitat, whereas services requiring more productive systems were hampered (sport fishing and bird watching).Overall, climate change showed negative impacts on the provisioning of ESs.

The effectiveness of restoration scenarios on ES provisioning
The technology-based solution simulated in this study focused solely on phosphorus reduction, as phosphorus is the commonly targeted element in geo-engineering approaches (e.g., through the application of lanthanum-modified bentonite).Solutions that also target N are available, such as dredging or some chemical amendments (e.g., a modified zeolite Z 2 G 1 ; Gibbs et al., 2011).In contrast, the nature-based solution, i. e., wetland restoration, in our simulation targeted both nitrogen and phosphorus simultaneously.Furthermore, it can also contribute to removal of organic materials, which not only consists of organic carbon, but also contains organic nutrients (Reinl et al., 2022).
The technology-based solution simulated was more effective at high degrees of intensity, suggesting a lagged system response to the reduction of phosphorus loading.Such delayed system response could be attributed to the relatively high background nutrient loading and eutrophic state of this lake, suggested by the high internal loading measured from the lake sediments (1.35 mg P/m 2 /day, see SI section 2).In lakes that are degraded for years, most of nutrients are locked in biological forms, which are less available for nutrient binding (Zhan et al., 2022).Moreover, the nutrients released from sediments will be brought into the upper water layer during wind mixing, which constantly charges the primary production (Schindler, 2006).
In contrast, wetland restoration showed an opposite relationship, i. e., decreased effectivity with larger wetland coverage.This is in contrast to previous studies on eutrophication control by wetland restoration suggesting that only a large percentage of wetland area led to a significant restoration effect (Janse et al., 2001;Sollie et al., 2008).However, these studies were carried out on shallow lakes.Deep lakes are inherently different due to their pronounced seasonal stratification and water column mixing (Wetzel, 2001).Deep lakes showed higher water column stability due to stratification (Crisman et al., 2005).A potential explanation of our contrasting results could be that deep systems tend to become less stable with increased wetland restoration, impairing the Fig. 6.Summary of effect size values using Hedges's g s .DO = dissolved oxygen, Cyano-Chla = cyanobacterial chlorophyll-a.Strength of evidence takes the form of no, small, medium and large, with positive Hedges's g s indicating an increase under future climate conditions and negative Hedges's g s indicating a decrease.Note that DO here stands for epilimnion DO concentration, and cyano-Chla stands for cyanobacterial chlorophyll-a in the upper layer (surface + epilimnion).Statistics for hypolimnion DO and cyanobacterial chlorophyll-a in the surface layer can be found in Tables SI-2.
wetland restoration effectiveness.
Overall, our results revealed non-linear impacts by restoration scenarios on the ecosystem state indicators as well as on the provisioning of ecosystem services.These results hint at the existence of a threshold relationship between restoration efforts and societal benefits in the form of ecosystem services (Iwasa et al., 2007).In other words, incremental increases in restoration efforts do not translate into proportionate enhancements in the water quality or the provisioning of ecosystem services.Our modeling outcomes underscore the presence of an optimum level of restoration efforts, that produces the desired outcomes.Hence, river basin authorities and policy makers should be aware of such non-linear responses in their ecosystems and their services, and make use of tools such as the one presented here to better understand the optimal effort needed to reach desired outcomes.

Climate change impacts on ES delivery
We used FLake to provide a more explicit description of mixing regime and vertical water temperature, which was then used to force PCLake + for predicting lake ecosystem dynamics under future climate conditions.The model performance was validated by comparing with observed water quality data under current climate condition.Note that we took a simplified approach to climate change in this study, focusing on direct impacts to the lake ecosystem itself.However, over the time span of restoration, compounding changes may also take place that are not considered in the model, from climatic impacts on the wider catchment, e.g., increased droughts, land use change, increasing runoff, to changing nutrient loads and load ratios.The predicted impacts of climate change on eutrophication symptoms, including higher surface cyanobacteria biomass and lower water transparency, align with findings from earlier studies on climate change effects on freshwater ecology (Jeppesen et al., 2020;Nielsen et al., 2014).FLake predicted a prolonged summer stratification in 2050 climate scenario compared to the current climate conditions, which is in line with previous modeling studies (Feldbauer et al., 2022).
Our results indicate that climate change, overall, may have negative impacts on eutrophication control and ES provisioning, with the effects being more pronounced on the direct ecosystem state indicators than the ES provisioning.This supports the conclusion from previous studies (Zhan et al., 2021(Zhan et al., , 2022) that an intensification of nutrient intervention measures could overcome the negative impacts by climate change.In other words, to acquire the same magnitude of eutrophication control or ES provisioning, higher intensities of restoration measures (both technical-based and nature-based solutions) will be required under future climate conditions.Moreover, the Hedges's g effect size tests on climate change impacts showed that wetland restoration had overall lower effect size values, which indicates that nature-based solutions as restoration approach may offer greater potential for climate-adaption and resilience of wetland ecosystems.

Strengths and future development possibilities of the ES modeling framework
By building an ecosystem service module into a modeling framework that incorporates complex physical and ecological processes, our approach was able to study the linkages between ecosystem states and ecosystem services in a more quantitative way, as a follow-up of the recommendations by previous studies (Janssen et al., 2021;Seelen et al., 2021).Using our approach, managers can make more quantitative estimations of ES provisioning under different restoration scenarios, their effectiveness, and how they might evolve in the future.
Our results suggest a conflict between good water quality (low primary production and biomass) and high fishery production, as claimed in previous studies (Matsuzaki et al., 2018;Seelen et al., 2021).By integrating multiple ESs into one framework, we were able to illustrate the contrasting requirements of ecosystem states by different services, thus studying possible trade-offs between ESs (Janssen et al., 2021).
Our ES modeling framework can pave the way for comprehensive cost-benefit assessment between restoration scenarios.The monetary cost of wetland restoration was demonstrated by field applications to be higher over an order of magnitude than in-lake measures (Huser et al., 2016).In addition, nature-based solutions such as wetland restoration have higher demands of space and time, which are usually limited in intensively populated urban areas (Cooke et al., 2018).Which measure to select requires a comprehensive evaluation of its effectiveness as well as economical cost.In this study, we studied the scenarios of two restoration measures individually in order to derive mechanistic understanding of their individual impacts.However, our modeling framework is capable of evaluating a mix of the two restoration measures to study their combined effects.
It is important to recognize that our modeling outcomes represent a knowledge-based estimation of how ecosystems might respond to the two restoration scenarios tested, built upon a set of assumptions such as an extended timeline for restoration to take effect.In practice, restoration effectiveness will depend on a combination of political, societal, economic and ecological factors.Our model offers managers and policy makers an evidence-based, first-order impact assessment from an ecological standpoint, highlighting potential benefits gained or lost in terms of ecosystem service provisioning due to an envisioned restoration scenario.Importantly, our model provides working hypotheses on restoration strategies that need to be validated by experimental approaches and/or observational data.
Our analysis is not exhaustive of all ecosystem services provided by freshwater lakes.It is limited by the ability to quantify ESs as well as the capabilities of PCLake+ of modeling the required ecosystem state indicators, being limited to water-based ecosystem services.In addition, the assessment of ES suitability often relies on people's perceptions and can vary across cultures (Pereira et al., 2020).For instance, people's perception of aquatic plants as nuisance is found to be dependent on their relation to the area, with visitors being less likely than residents to perceive macrophytes as a nuisance (Hussner et al., 2017).Seelen's et al. (2019) survey data shows Europeans greatly underestimate their personal water use.To achieve a wider range of ecosystem services in our modeling framework, we envision expansions and improvements to the current model, such as making greenhouse gas fluxes explicit (Santos et al., 2022) and utilizing our model across lake networks (Teurlincx et al., 2019;van Wijk et al., 2022).Our modeling framework is designed to facilitate the incorporation of new and adjusted ecosystem services and is easy to modify.Importantly, making it open source allows for the synthesis of multidisciplinary knowledge required for ecosystem service assessment.

Conclusion
Incorporating an ecosystem service module into a modeling framework that considers complex physical and ecological processes has allowed us to explore the linkages between ecosystem states and ecosystem services in a quantitative manner.Our study employed this framework to a quarry lake, and evaluated the impacts of two restoration scenarios on ecosystem state indicators and provisioning of ecosystem services (ESs) under current and 2050 climate scenarios, leading to several take-home messages: 1.Our scenario analyses revealed non-linear relationships between the level of restoration intensity and the resulting outcomes of ecosystem service provisioning.Phosphorus reduction scenarios demonstrated increasing effectiveness with higher intensities, while the intensive wetland construction tended to exhibit decreasing effectiveness in the studied deep lake.2. Both measures showed positive effects on ESs that requires good water quality, such as swimming, irrigation, and macrophyte habitat.However, they had negative effects on ESs that require more productive systems, such as fishing and bird watching.3. Climate change had adverse impacts on the effectiveness of restoration measures regarding ecosystem state indicators and ESs provisioning.To achieve the same level of eutrophication control and ES provisioning, greater intensities of restoration measures (both technology-based and nature-based solutions) will be necessary under future climate conditions.Notably, our results indicate that nature-based restoration may display greater resilience to climate change, as evidenced by their overall weaker climate change effects.
Our ecosystem service modeling framework equips managers with the instruments to quantitatively estimate ES provisioning under various restoration scenarios, assess their effectiveness, and anticipate how the restoration impacts may evolve in the future.In conclusion, this framework is a valuable resource for decision-makers seeking optimal restoration strategies while considering the challenges posed by climate change.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Fig. 1 .
Fig. 1.Model chain for ecosystem service modeling.Rectangles denote state variables, ovals denote models, hexagon denotes ecosystem service module, rounded rectangles denote input data, solid arrows denote model input or output, dashed arrows denote data input.(Flake in blue, PCLake+ in green, input in white, output in orange).

Fig. 2 .
Fig. 2. Model results validation.Panel a: predicted mixing depth dynamics by FLake under present climate and future climate.Panels (c-d): Comparison between observed and predicted water quality variables under current climate conditions.The predicted oxygen concentrations are plotted for both the epilimnion (_Epi), while a surface layer (_Surf) is plotted in the total nitrogen and Chlorophyll-a which models surface cyanobacterial accumulation.

Fig. 3 .
Fig. 3.The response of ecosystem state indicators to increasing intensity of a technical restoration (dark yellow) or a nature-based solution (green).The solid lines represent the current climate, while the dashed lines represent the future 2050 climate scenario.

Fig. 4 .
Fig. 4. The response of ecosystem state indicators to increasing intensity of a technical restoration (dark yellow) or a nature-based solution (green).The solid lines represent the current climate, while the dashed lines represent the future 2050 climate scenario.

Fig. 5 .
Fig. 5.The responses of ecosystem service provisioning to P reduction (left panels) or wetland restoration (right panels).The solid lines represent the current climate, while the dashed lines represent the future 2050 climate conditions.Icons are in the same colors as the corresponding services for illustrative purpose.

Table 1 List of Ecosystem Services, their corresponding ecosystem state indicators and threshold values being included in the modeling framework. Ecosystem
state indicators with a symbol (*) represent an adjustment from Seelen et al. (2022), with the decisions explained in SI section 2.Category Service CICES Code Ecosystem state indicators Threshold values Provisioning Professional fishing -fishponds 1.1.4.1 Steady state fish density (kg/ha) >100 (suitable), 10-100 (moderate), Q. Zhan et al.