Estimating Combined Effects of Climate Change and Land Cover Change on Water Regulation Services of Urban Wetlands in Valdivia, Chile

Abstract

The relationship between cities and wetland cover varies across the globe, with some cities converting wetlands to low- and high-density urban cover and others preserving, conserving, or restoring wetlands, or constructing new ones. However, the scientific literature lacks studies relating changes in systemic flood risk in an urban stormwater management systems to changes in wetland cover. Furthermore, whether and how such relationships are affected by changing storm intensity under climate change is unknown. We present a case study on the effects of changes in urban wetland extent and storm intensity on flooding in an urban drainage system in Valdivia, Chile, under several co-produced future scenarios and historical trends of development. We used data derived from stakeholder workshops and historical landcover to determine four plausible scenarios of urban development, plus one business-as-usual scenario, in Valdivia through the year 2080. Additionally, we used historical precipitation data and downscaled climate data to estimate event rainfall from extreme storms in the year 2080. We found that system flood volume and time the system was flooded increased with increasing wetland loss and rainfall volume. Mean rate and hour of peak discharge were unaffected by wetland loss. Infiltration's relative role in reducing flooding diminished as wetland loss increased. Cities may still experience dangerous and/or unacceptable flooding even with extensive wetland coverage and will likely need to pair conservation with additional improvements in their stormwater management systems and contributing watersheds.

Introduction

Pluvial flooding is a major concern for residents of cities. Pluvial flooding is surface ponding or overland flow that occurs when rates of precipitation exceed the capacity of drainage systems and/or surfaces to remove it (Falconer et al., 2009). Pluvial floods can lead to loss of life, damage to property, and disruption of transportation networks (Chang et al., 2010; Douglas et al., 2010; Falconer et al., 2009; Yin et al., 2016). As a physical phenomenon, pluvial flooding results from interactions between rate of precipitation, urban stormwater management practices, and biophysical characteristics of the urban and peri-urban landscape (Westra et al., 2014). In many cities, one or all three of these interacting factors are changing in ways that may increase pluvial flood frequency, area, and damage. Even subdaily extreme rainfall has become more frequent and intense due to anthropogenic climate change (Westra et al., 2014; Wuebbles et al., 2014). Cities have historically prioritized mitigating the risks of fluvial and coastal flooding over pluvial flooding (Guerreiro et al., 2017). However, in recognition of pluvial flooding has in recent years garnered the attention of researchers and planners because understanding how to mitigate its causes and effects in urban areas is underdeveloped (Rosenzweig et al., 2018).

The conservation, restoration, and construction of wetlands have all been suggested as measures to mitigate the risk of various forms of flooding in many different ecosystem types. The ability of coastal wetlands to reduce coastal flooding has been explored in depth and in a diverse array of ecosystems (Arkema et al., 2013; Narayan et al., 2017; Nicholls et al., 1999; Rojas et al., 2019; Van Coppenolle & Temmerman, 2019, 2020). The effects of wetland presence on riverine flooding have received notable attention as well. Neri-Flores et al. (2019) modeled the capacity of wetland preservation to reduce riverine flooding caused by hurricane storm surges. Pomeroy et al. (2014) modeled how preserved inland wetlands can reduce riverine flooding driven by snowmelt. Yang et al. (2010) modeled how the restoration of wetlands in a Canadian prairie watershed can reduce peak river discharge and flooding. In a review of 28 modeling and empirical studies of the effects of wetlands on flow regimes in rivers, Kadykalo and Findlay (2016) found that wetlands generally reduced the frequency and magnitude of flooding, with one exception in a forest wetland system (Lundin, 1994). Historically, attributions of the positive water regulation services of wetlands have their bases in studies in non-urban riverine or coastal wetlands (Costanza et al., 1997; Millennium Ecosystem Assessment, 2005).

Only recently has research explored the abilities of inland urban wetlands to reduce urban pluvial flood risk, or how the incorporation of wetlands in an urban stormwater management system might alter the system's performance. The theory and practice of inland wetland restoration and construction in urban areas to reduce pluvial flood risk is relatively new in academia and among stormwater managers (Chan et al., 2018; Elmqvist et al., 2015), and modeling and empirical studies of the effects of wetland restoration and construction in urban areas are rare. Some cities have added inland wetlands to their portfolios of green stormwater infrastructure (GSI; otherwise referred to as a form of green infrastructure, urban ecological infrastructure (Childers et al., 2019), or, more broadly, nature-based solutions) or suggested that the construction, restoration, or incorporation of inland wetlands be included in sustainable urban drainage systems or low-impact development strategies to reduce pluvial flooding (Chan et al., 2018; Fletcher et al., 2015; Y. Li et al., 2020).

Wetlands may provide water-regulation services to cities through a variety of hydrologic processes. Depending on wetland morphology, wetland vegetation, environmental conditions, soil characteristics, water-table depth, and connectivity to drainage systems to which wetlands may be connected, wetlands may manage stormwater via some combination of impoundment (the temporary storage of water), infiltration (the removal of surface water via percolation into wetland soils), evapotranspiration (the removal of surface and soil water from the system via evaporation or plant-mediated transpiration), and conveyance (the movement of water through and out of the drainage system via passive flow; Bullock & Acreman, 2003). For many cities considering the use of wetland GSI, the key hydrologic functions of wetlands are those of detention and infiltration (Y. Li et al., 2020). Detention of stormwater in wetlands delays or reduces stormwater release to downstream waterways (Kadykalo & Findlay, 2016). Infiltration, facilitated by wetlands through their pervious soils, reduces the proportion of precipitation that converts to runoff (Fletcher et al., 2013). Widespread impervious cover in cities leads to high rates of conversion of precipitation to runoff, which in turn increases peak rates of discharge in drainage systems and can overwhelm the drainage system flood connected areas (Ogden et al., 2011).

Critically absent from the literature on the flood-mitigation services of wetlands are city-wide studies on how performance of the urban stormwater management system changes when urban wetlands are constructed, restored, or incorporated. Change in the value of water-regulation service of urban wetlands over is often estimated using simple land-use or land-cover change and look-up tables of water regulation service values according to regional wetland area (G. Li et al., 2022; Mukherjee et al., 2021). Such estimates assume water regulation services absent any details or consideration of the stormwater management system to which they are connected (C. Wang et al., 2018; Y. Wang et al., 2018; Zhang et al., 2020). But outside of an urban context, it is widely recognized that system-specific knowledge is necessary to accurately estimate effects of wetlands on the water regulation services that wetlands may provide (Acreman & Holden, 2013; Kadykalo & Findlay, 2016). Wetland dimensions, extent, antecedent storage conditions, rates of infiltration and evapotranspiration, and configuration within a stormwater management system are all likely to influence the performance of urban stormwater management systems.

While wetland GSI is often recommended to increase resilience against floods in cities under a changing climate (Stefanakis, 2019), its efficacy should not be taken for granted. Climate change will shift storm intensity and timing away from the conditions for which stormwater management systems, even those with wetland GSI, were designed, which are generally historical storms (ASCE/Environmental & Water Resources Institute, 2006). Sensitivity of the drainage system response to changes in precipitation intensity from climate change depends on, for example, the size of the contributing watershed and the size and configuration of wetland GSI within the system. Yet studies that espouse the benefits of wetland GSI for increasing resilience in the face of climate change rarely contextualize those benefits in terms of the scale of the flood risk that climate change poses, or examine how performance of systems with wetland GSI might also change with the climate.

In the present study, we modeled the coupled effects of inland wetland loss and impervious watershed expansion on stormwater management system performance under different scenarios of climate change. For the study system, Valdivia, Región los Ríos, Chile, we asked the following question: How does the loss of wetland GSI in an urban stormwater management system change the system's flood volume, peak discharge rate, and peak discharge timing? We hypothesized flood volume and rate of peak discharge would increase, and the hour of peak discharge would arrive earlier, with wetland loss. Additionally, we asked: how do the effects of wetland loss on flooding compare to the effects of changing rainfall during extreme storms? We hypothesized that there would be more systemic flooding, longer periods of flooding, and that peak discharge would be greater and arrive earlier due to increasing rainfall than by wetland loss. Finally, we asked: how much does infiltration contribute to flood reduction as wetland loss increases? We hypothesized that infiltration would contribute to lower flood volume and reduce flood duration under all extents of wetland loss.

Materials and Methods Study Site

Valdivia, Chile (area: 93.94 km²) is a city of approximately 166,000 people in the southern half of Chile, 850 km south of the capital Santiago, in the Región de los Ríos (Figure 1). Citizens and stormwater managers in Valdivia must contend with a high risk of pluvial flooding owing to high average annual precipitation, a long rainy season, the city's location 12 km inland from the Pacific Ocean, at the confluence of three rivers, and patterns of land development. Valdivia's ecosystem is classified as a temperate rainforest (Amigo & Ramirez, 1998; Hajek & Di Castri, 1975). Wetlands are a characteristic feature of Valdivia, covering 20.64 km² (22.7%) of the municipal area but are at risk from continued development.

Figure 1. Left: Location of study site, Valdivia, Chile (39.8336°S, 73.2154°W). Right: Valdivia's land cover based on spectral analysis of a 2010 orthophoto, and drainage system, as described in 2012 by the Chilean Ministry of Public Works.

Valdivia's average annual rainfall was approximately 1719.48 mm between 1990 and 2021 (Dirección General de Aeronáutica Civil, 2020), with pronounced droughts in the last decade. In 2015, rainfall in the region and snowpack in the Andés were low enough that the riverine potable water supply became too saline for treatment due to tidally forced saltwater intrusions from the nearby ocean, and the city was forced to pump groundwater for nearly all of its supply (Garcés-Vargas et al., 2020). In 2021, for the first time since the city began measuring precipitation at the nearby Pichoy Airport meteorological station in 1969, the city registered less than 1,000 mm of precipitation (Sepúlveda, 2021). This extreme departure from the prevailing rainfall patterns has added to concerns about sustainability and resilience in Valdivia under climate change (Garcés-Vargas et al., 2020; Sepúlveda, 2021).

Valdivia's stormwater management system is composed primarily of gray infrastructure components (e.g., pipes and canals), wetlands, and the rivers into which the system ultimately discharges. As of 2012, Valdivia's stormwater management system consists of roughly 245.7 km of drainage infrastructure, of which 41.2 km (16.8%) is wetland GSI. The origin of most of wetland cover in the city is a 1960 earthquake of magnitude 9.5, which caused up to 20 m of uplift in some areas and subsidence and rifting in others (Barrientos & Ward, 1990). Since the earthquake, the city has deliberately incorporated many of these wetlands into its stormwater management system (CMOP, 2012). In addition, the presence of wetlands in the city is owed in part to local conservation movements to maintain the cultural services of wetlands (Correa et al., 2018) and their function as habitat to charismatic species (e.g., Cygnus melancoryphus) tied to Valdivian identity (Silva et al., 2015).

General Approach

We used model estimates of future land cover change and estimates of future extreme rainfall as inputs to a 1-dimensional model of Valdivia's stormwater management system, and ultimately produced estimates of flood volume and flood location for a range of land cover and climate conditions (Figure 2). This process began by convening an in-person workshop in Valdivia, Chile to co-develop with practitioners the goals and objectives of four different scenarios of development for the city to achieve by the year 2080. We then combined historical data on land-cover change in Valdivia and scenario goals and objectives into rules governing land-cover change in the Dinamica EGO cellular automata-based model (Soares-Filho et al., 2002). The outputs of this model were five land-cover maps: one for each of the four scenarios developed in the workshop at the start of this process, along with an additional “business-as-usual” scenario estimating land-cover change in the absence of interventions to the status quo. We then used ArcGIS Pro (ESRI) to estimate changes in wetland volume and subcatchment area as a result of the changes in land cover areas in the five land cover maps. Separately, we used daily precipitation estimates from downscaled climate models to estimate rainfall of 100-year return period, 24-hr duration storms in the year 2080 under various climate conditions. Estimated changes in wetland volume and subcatchment area, as well as estimated changes to rainfall during extreme storms, were used to construct a 1-dimensional model of Valdivia's stormwater management system under various land-cover and climate configurations in the year 2080. This 1-dimensional model was then used to estimate flood characteristics that varied by land-cover and climate configurations.

Figure 2. Process diagram detailing convergent processes used in this study to produce estimates of flood characteristics under a range of land cover and climate conditions. The left branch represents work done to produce spatial estimates of land cover change by the year 2080 and to translate these changes to land cover to changes in the morphology of wetlands and watershed areas. The right branch represents work done to produce estimates of rainfall during 100-year interval, 24-hr duration storms.

Stormwater Management Model Characteristics and Calibration

In 2002, Chile's Ministry of Public Works (CMOP) commissioned the development a hydrologic model of the city's surface and stormwater management system flows using the Environmental Protection Agency's Stormwater Management Model (CMOP, 2012; EPA SWMM; Rossman, 2015). EPA SWMM is a 1-dimensional hydrologic model that converts rainfall to runoff for each subcatchment and routes this water through conduits and nodes. The model is commonly used to design and assess the performance of stormwater management systems in urban areas (Choo et al., 2021; Gülbaz et al., 2019; Iffland et al., 2021). Valdivia's stormwater management model (SWMM) was updated in 2012 to include system expansions and observational delineation of the city's urban subcatchments, among other updates and improvements (CMOP, 2012). The 2012 SWMM also included a tidal outfall curve to account for changing water levels in the rivers to which the stormwater management system interacts, peaking on hour two of simulation at 1.46 m above invert elevation, and on hour 14 lowering to 0.28 m above outfall invert elevations, and repeating every 12 hours until simulation completion. This curve was designed to represent an annual average difference in water levels at the outfalls under historical river and ocean-level conditions. This curve was conserved in our final models.

Valdivia's SWMM was calibrated using observed stormflow data from seven storms of different return periods, ranging from 0.67 to 24.52 years, in a sub-section of the larger SWMM (CMOP, 2012). The model was optimized to achieve similar rates of peak discharge and flood volume to those observed through manipulating parameters like Manning's roughness and rates of infiltration for pervious and impervious surfaces for the observed storms (Table 1). These calibrated values were conserved in our final models. The absolute differences between the simulated and observed flood volume and rate of peak discharge for each storm for the final values of these parameters range from 1% to 74% for flood volume and from 5% to 86% for peak discharge rates (Table 2; CMOP, 2012). Notably, the model was not calibrated using observed events with return periods greater than 24.5 years. Published reports of EPA SWMM models that estimated flooding for whole urban watersheds are uncommon, particularly those that estimate the effects of large magnitude storms (e.g., 10-year or greater) over long durations (e.g., 24-hr); however, for context, two studies examining the effects of storms of much lesser magnitude than we examined, but nonetheless in whole urban watersheds, reported relative errors between simulation and observation flood volumes between 5% and 20% (Wu et al., 2018) and between 1% and 100% (Barco et al., 2008) depending on the range of input storm magnitudes and the method of optimization.

Table 1 Calibrated Parameter Values Used in EPA SWMM Models for Valdivia (CMOP, 2012)

Parameter	Value
Manning's n, impervious (N-Imperv)	0.03
Manning's n, pervious (N-Perv)	0.09
Depression storage-impervious (Dstore-Imperv; mm)	1.25
Soil moisture retention, pervious (S-pervious; mm)	5
Percent with no depression storage (% Zero-Impervious; %)	80
Rate of infiltration, minimum	3
Rate of infiltration, maximum	4
Decay rate (seconds−1)	2
Drying time (day)	7
Evaporation (mm day−1)	2

Note. EPA SWMM parameter names and units in parentheses.

Table 2 Differences Between Simulated Model and Observational Flood Volume and Peak Discharge Rate for Storms of Different Return Periods (CMOP, 2012)

Storm return period (years)	Simulated flood volume (%)	Simulated peak rate of discharge (%)
0.67	+1	+20
0.88	−3	+5
0.94	−29	−36
1.78	+24	−15
2.40	+74	+86
6.56	−3	−5
24.5	+1	+20

Estimating Future Land Cover Scenarios and Wetland Dimensions

In May of 2017, the Urban Resilience to Extremes (UREx) Sustainability Research Network (SRN) hosted a workshop in Valdivia, Chile, to envision a series of long-term (2080) future scenarios and desirable future pathways of urban development. Participants in the workshop represented a diverse array of Valdivia's stakeholders, such as municipal and regional government employees, university professors, students, and members of community action groups. Participants collaborated to develop a suite of visions and strategies to undertake in order to achieve four unique, plausible scenarios for a future Valdivia: Inclusive City, Friendly City, Eco-Wetland City, and Resilient-to-Flood City. The scenario themes emerged from the concerns of the citizens of Valdivia and an analysis of Valdivia's governance documents as well as a publication from the Inter-American Development Bank (IDB, 2015). The visioning and scenario development process in the workshop followed methods described by Iwaniec et al. (2020).

The qualitative strategies of four scenarios—Inclusive, Flood Resilient, Friendly, and Eco-Wetland—developed in Valdivia's workshops were translated by the UREx SRN modeling team into quantitative spatial and temporal rules and introduced into cellular automata-based models of land-use/land cover (LULC). This phase represents an iterative process in which the modeling team gathered feedback from various stakeholders on the four co-produced scenarios, adjusted the quantitative rules based on that feedback, and released updated simulations. Paired with historical information on LULC changes (observed 1983 and 2010 LULC maps) in Valdivia, the cellular automata-based Dinamica Environment for Geoprocessing Objects GO model (Soares-Filho et al., 2001, 2002), hereafter Dinamica, generated predictions of LULC configuration in Valdivia in 2080 for each scenario, as well as for a “Business-as-usual” (BAU) scenario, which assumes LULC proceeded entirely according to historical patterns of development. Dinamica has been used to simulate LULC change in many studies (e.g., Gago-Silva et al., 2017; Kolb et al., 2013; Pathirana et al., 2014). Dinamica estimates LULC change quantity using a transition matrix obtained from the cross-tabulation of the observed LULC data. The transition matrix is then transformed into a Markovian Chain Probability Matrix, which computes the average percentage of each land class that changes to another class at each time-step (in our case, 1 year) which is the transition rate. Dinamica then spatially allocates the quantity of LULC change according to a transition rule with two components. The first component calculates transition probabilities of LULC-change global drivers (explanatory variables such as accessibility, elevation, and slope). The second component considers the influence of local neighbors on the transition of the LULC state of a cell. Dinamica adopts the Weights of Evidence method (Soares-Filho et al., 2002, 2004) to quantify the influence, or the weight for a set of explanatory variables, based on the occurrence of each LULC in specific ranges. Dinamica calculates the influence of local neighbors on each cell in the landscape using two complementary functions: Expander and Patcher, one to expand/contract previous LULC patches and one to generate new ones, as described in depth in Soares-Filho et al. (2002).

In all scenarios, wetland cover declined compared to the 2010 base map (Figure 3). However, co-developed scenarios showed lower wetland loss rates than the BAU scenario. Stakeholder proposals from the workshops played a significant role in determining the loss rate. For example, in the Inclusive scenario, for example, the proposal to “create a network of wetlands for connectivity within the city, and wetlands are protected and an important part of mitigating climate change impacts” by 2050 led us to introduce new wetland corridors and stop converting wetlands to other uses. Although some wetlands were lost (converted to other land uses-especially built-up) before 2050, the addition of new wetland corridors helped reduce overall loss over time. The Eco-wetland scenario did not include this specific role, resulting in a slightly higher wetland loss rate compared to the Inclusive scenario. Also in the Eco-wetland scenario, a proposal of declaring wetlands as protected zones and implementing a 100% prohibition of wetland filling by 2040 was essential for preserving more wetlands. However, some wetlands were still converted to other land uses before 2040 before the prohibition toggled on. Finally, many wetlands within the present-day and scenario land-cover maps are not included within the city's stormwater management model. As a result, the change in wetland area in the subset of wetlands in the SWMM differed from the change in wetland area for the whole city in the cellular automata-based models. In scenarios like Eco-wetland, where wetland cover overall was greater than in other scenarios like Friendly, much of its conserved or gained wetland cover was in the northwest and west where the SWMM model did not extend, while the wetland cover it lost was within the wetlands included in the SWMM model (Figure 3).

Figure 3. Land cover in the present day (2010) and under five scenarios of development by the year 2080. Wetland loss generally increases from left to right, and from top to bottom, compared to the present day. City-wide wetland loss for each scenario was: 9.72% in Inclusive, 13.3% in Resilient-to-flood, 18.3% in Friendly, 23.98% in Eco-wetland, and 37.3% in Business-as-usual compared to Present Day wetland coverage.

We determined change in wetland area by overlaying present-day land cover with scenario land cover in ArcGIS Pro (ESRI, 2019) and removed wetland area that existed in the present day that converted to low- or high-density urban land cover in the future scenarios. Conversion of wetland area to either form of urban land cover necessitates the in-filling and elevating of the former wetland's surface and reduces wetland storage capacity. In contrast, conversion of wetland area to either pasture/green or forest land cover types does not necessitate in-filling or affect storage capacity.

We then calculated wetland volume and change in wetland volume that resulted from change in wetland area. A 2019 contour map (1-m vertical resolution) of Valdivia was converted into a triangulated irregular network (TIN), which characterized the three-dimensional topography of the landscape. For each of the wetlands in the SWMM, for the present day and each scenario, wetland boundaries were used to generate pseudo-surfaces every 0.25 m from the base of each wetland to their lowest bank, and the volume of the TIN underneath the pseudo-surface was calculated using the Surface Volume tool in ArcGIS Pro.

The wetlands included in the SWMM were modeled as one of two elements in EPA SWMM: storage units or conduits. Wetlands with single inflows from subwatersheds, and wetlands that were spatially isolated from connecting wetlands, were generally modeled as single storage units with shape and volume determined by the previous step (Figure 4). Wetlands with multiple inflows, and that were only separated from other wetlands by short pipe segments under roadways, were generally modeled as a series of conduits linked by nodes (Figure 4). Modeling wetlands as storage units or as a series of conduits and nodes affects flow timing in the model, as a parcel of water moves in and out of a storage unit instantaneously but requires time to move through a conduit, but it is nonetheless accepted practice to model wetlands as storage units (Knighton et al., 2016).

Figure 4. Example wetlands in Valdivia illustrating differences in the construction of storage unit and conduit wetlands in EPA SWMM. The wetland on the left receives water from a single subcatchment and was modeled as a storage unit. The wetland on the right receives water from multiple subcatchments and was modeled as a series of conduits and nodes.

Owing to a high natural water table, proximity to three rivers, and high annual rainfall, the model developers assumed no infiltration in Valdivia's wetlands. While this may be an acceptable assumption during the rainy season (June–September) when the water table is particularly high, our own observations indicated substantial potential for infiltration in Valdivia's wetland soils during the summer months (December–February) when temperature and insolation are high and months may pass without rain. In Section 2.7, we attempted to account for this potential for infiltration in an experimental model subsection. Initial water levels in the calibrated model were set to zero, which may reflect summer conditions but not winter conditions. Data on groundwater inputs to wetlands were not available for this investigation, though our field data collected for as-yet unpublished research indicated primarily unidirectional flow from wetlands outward to the city's rivers. While wetlands in Valdivia are typically depressional they nonetheless are perched higher than river water levels, even at high tide.

Changes to wetland volume, as calculated in the previous step in ArcGIS Pro, were translated to the SWMM by conserving bank elevation, depth, and length, but, in the case that the wetland was modeled as a conduit its length, by altering cross-sectional width (referred to in EPA SWMM as station), such that the overall wetland volume was the same between ArcGIS and EPA SWMM. Subcatchment areas in the SWMM were increased by the amount of wetland area lost to low- and high-density urban land cover between the present day and the scenarios. In the case that a wetland was only connected to a single subcatchment, all lost wetland area was added to the subcatchment. In the case that multiple subcatchments were connected to a wetland, the subcatchments expanded according to the amount of nearby wetland lost. No other subcatchment properties, such as imperviousness or rates of infiltration, were changed, as it was assumed that new low- and high-density urban subcatchment area would be roughly the same as the present-day low- and high-density subcatchment area.

Downscaling Climate Models to Valdivia, Chile

We employed asynchronous regional regression models to downscale precipitation estimates from atmosphere-ocean general circulation models to Valdivia, Chile (Stoner et al., 2012). Input data were historical observational data on rainfall from the Pichoy Airport meteorological station, located roughly 32 km (22 miles) from Valdivia's centroid. This station has the most consistent and longest rainfall record of any station either within or around the city. These downscaled models produced estimates of daily precipitation for the years 1969–2080. Additional information on the downscaling methods can be found in Supporting Information S1.

Estimating Rainfall Volume of a 100-Year Return Period, 24-hr Storm

Estimated rainfall of historical and future 100-year return interval, 24-hr duration storms were derived from the generalized extreme value (GEV) distribution. The GEV distribution is commonly employed for modeling extremes in rainfall such as extreme events of various return periods (Bella et al., 2020; Reiss & Thomas, 2007). It is the combination of three extreme value distributions (Gumbel, Fréchet, and Weibull distributions), and can be represented by the following equation: [Image Omitted. See PDF]where, G(z) is the probability that the monthly precipitation will be greater than or equal to z mm, μ is the location parameter, σ is the scale parameter, and ξ is the shape parameter. The return period, T, of a rainfall amount greater than z mm can be represented as: [Image Omitted. See PDF]where P(exceedance) is the probability of event exceeding rainfall amount z. Thus, the return period can be related to the GEV via these intermediate equations: [Image Omitted. See PDF] [Image Omitted. See PDF] [Image Omitted. See PDF]

Using Equation 5, 1200 months for T, and the μ, σ, and ξ parameters from fitting the GEV model to the data, we estimated the precipitation expected to fall in a 100-year return period, 24-hr duration storm. Observational and modeled daily precipitation data were first used to determine monthly maximum rainfall. These monthly extremes were then grouped into periods of 30 years (e.g., 1986–2015) and fitted to a GEV distribution using the extRemes package in R (Figure 5; Gilleland & Katz, 2016; R Core Team, 2021).

Figure 5. Mean rainfall (solid lines) with 95% confidence intervals (transparent bands) of 100-year return interval, 24-hr duration storms over time for (top) the RCP 4.5 warming scenario and (bottom) the RCP 8.5 warming scenario. The thick black line and bands ending in 2021 represent event rainfall derived from observed data collected at the Pichoy Airport meteorological station. The thick dark purple line and band extending to 2080 represent the ensemble average event rainfall and 95% CI, respectively, for the generalized extreme value model estimates. Other colored lines and bands represent estimates for individual climate models. Gray dotted lines at years 2021 and 2080 mark the end of the observational and scenario modeling periods, respectively. The RMSE of the ensemble average for the mean event rainfall from 1998 to 2015 was 19.0 and 17.5 for the RCP 4.5 and RCP 8.5 climate scenarios, respectively, when compared to the observational period. The RMSE of ensemble average for the mean event rainfall from 1998 to 2021 were 16.9 and 16.0 for the RCP 4.5 and RCP 8.5 climate scenarios, respectively.

From the GEV distributions, we estimated the ensemble average rainfall amount expected during a 100-year return period, 24-hr duration storm (i.e., event rainfall), along with the 95% confidence interval, for the observational period, and the RCP 4.5 and RCP 8.5 climate scenarios (Table 3). The estimates for mean event rainfall in 2080 were 262.9 ± 1.3 mm under the RCP 4.5 climate scenario and 297.1 ± 1.3 mm under the RCP 8.5 climate scenario (Figure 5). These estimates represented increases of 15.9% and 31.0% in event rainfall by the year 2080 for the RCP 4.5 and RCP 8.5 climate scenarios, respectively, compared to the event rainfall for the year 2021 (Figure 5).

Table 3 Estimated Rainfall From 100-Year Return Period, 24-Hour Duration Storms (Event Rainfall)

	Event rainfall (mm)
	2021	2080, RCP 4.5	2080, RCP 8.5
Lower 95% CI	147.4	155.9	160.4
Mean	226.8	262.9	297.1
Upper 95% CI	306.2	369.9	433.9

Note. The estimates for 2021 were derived from the observational record, while the estimates for 2080 were derived from downscaled model estimates of daily precipitation for RCP 4.5 and RCP 8.5 scenarios of global warming.

From the GEV models, rainfall amount during extreme storms and extreme storm frequency were both estimated to increase by 2080 compared to 2021. From author experience monitoring wetland levels in 2019, when annual rainfall (1071.3 mm) was lower than the 30-year average (1715.2 mm), wetland stage was often at or near 0 m for much of the year. A drier climate and less frequent storms in 2080 would likely reduce surface and soil water stores in the region, leading to low wetlands stage, perhaps at or near 0 m except during the rainy season. Wetland stage in the SWMMs was then set at 0 m at the start of modeling to reflect drier surface storage conditions.

Modeled Hydrologic Processes Using EPA SWMM

EPA SWMM allows modelers to simulate infiltration in storage units by toggling an infiltration function, but this function is not available for conduits. In our SWMM for Valdivia, the majority of wetlands are conduits and would not allow for infiltration due to this model limitation. In our model simulating whole-city flood risk, we did not allow for infiltration in either storage units or conduits in order to consistently model the same hydrologic processes—storage and conveyance—in both model elements.

However, in a separate model we modified a subwatershed that contained one updrain and one downdrain wetland, where both were modeled as storage units, to allow for infiltration (Figure 6). Storm inputs for this subwatershed were the mean, lower CI, and upper CI event rainfall for the present day and 2080 from the RCP 8.5 climate model (Table 3). The soil moisture deficit in this wetland and park was set to 7.8%, which corresponds to the deficit in the driest soils we observed in a separate field study of wetland soil moisture content in Valdivia, and which could reflect the soil moisture deficit of Valdivia's wetlands in 2080. Rates of wetland loss in the two wetlands in this subsystem model were greater than the rate of wetland loss in the larger urban system (Figure 6).

Figure 6. The simple subwatershed used for exploring the contributions of infiltration on flood mitigation in the subsystem model. The subsystem consists of two storage unit wetlands and six storm drain nodes and it drains into a more complex downdrain watershed after the most southern node. Subsystem wetland loss for each scenario was: 30% in Friendly, 35% in Inclusive, 66% in Resilient-to-flood, 76% in Eco-wetland, and 97% in Business-as-usual compared to Present Day wetland coverage.

In the whole-city system and subsystem models, flood volume was calculated as the sum of flood volumes occurring in all nodes, including wetlands modeled as storage units and nodes within wetlands modeled as conduits, as tabulated in the EPA SWMM run reports. Our reported values for flooding then only consider drainage system flooding occurring outward from drainage system components and not the runoff produced by subcatchments.

Flood duration was calculated by summing the duration of flooding in each node. Mean rate of peak discharge (m³/s) and mean hour of peak discharge (from hour 0 (simulation start) to hour 48 (simulation end)) were calculated as the mean values of rate of peak discharge and hour of peak discharge, respectively, for all components in the system or subsystem.

Results

Flood volume increased linearly with wetland loss for all warming scenarios, thus more event rainfall was converted to flooding as wetland loss increased for any given rainfall amount (R² > 0.97 and p < 0.01 for all inputs; Figure 7a). Confidence intervals (95%) of mean values of flood volume did not overlap between warming scenarios; however, for rainfall representing the lower and upper 95% CI of each warming scenario, flood volume overlapped between emission scenarios. The mean slope of the trendlines of change in flooding per loss of wetland area was 11.0 × 10³/% of wetland loss. Flood duration also increased linearly with wetland loss for all warming scenarios (R² > 0.96; p < 0.01 for all inputs; Figure 7b). Mean hour of peak discharge increased linearly (R² = 0.77; p < 0.05) with wetland loss under mean storm rainfall in 2021 but not under other storm rainfall amounts (Figure 7c). Mean rate of peak discharge did not significantly change with wetland loss under any storm rainfall input (Figure 7d).

Figure 7. (a) Flood volume, (b) flood duration, (c) mean hour of peak discharge, and (d) mean rate of peak discharge versus wetland loss of whole-city drainage system. Points represent ratio values given mean values of rainfall as input for a given warming scenario. Ribbons indicate range of ratio values given lower and upper 95% CI values of rainfall as inputs for a given warming scenario. Presence of trendlines indicates significant linear relationship (all R2 > 0.97; p [less than] 0.01).

Flood volume increased linearly with rainfall in all scenarios (i.e., all amounts of wetland loss), indicating that more event rainfall was converted to flooding as rainfall amounts increased for any given amount of wetland loss (R² > 0.98 and p < 0.05 for all inputs; Figure 8). Confidence intervals (95%) overlapped between all scenarios, though maximum and minimum values for the confidence intervals increased with wetland loss from the present-day scenario to the business-as-usual scenario. The mean of the slopes of the trendlines was 27.0 × 10³ m³/mm, indicating the rate at which flood volume increased with each mm unit increase of rainfall (Figure 8).

Figure 8. Flood volume versus event rainfall for all land-cover scenarios (i.e., amounts of wetland loss). Input rainfall amounts are those included in Table 3. Ribbons indicate 95% confidence intervals of ratio values for given scenarios. Presence of trendlines indicates significant relationship (p [less than] 0.05). All R2 values were ≥0.98 and were significant at the p [less than] 0.01 level.

In the subwatershed, where flooding was compared between models that did and did not allow infiltration in the lone, updrain wetland of the subsystem, flood volume increased with wetland loss under all warming scenarios for models allowing infiltration (p < 0.05), similar to the models that did not allow infiltration (p < 0.05; Figure 9a). The effect of infiltration on flood volume diminished as wetland loss increased, as there was less area to allow for infiltration. In this subsystem, when infiltration was allowed and disallowed, flooding occurred at the updrain wetland storage unit as well as in downdrain nodes until nodes connected to the downdrain wetland (Figure 9a). The duration of flooding in the subsystem significantly increased with wetland loss when infiltration was allowed but not when infiltration was not allowed (p < 0.05; Figure 9b). The time the updrain wetland flooded also increased with increasing wetland loss (p < 0.05; Figure 9c). Under some rainfall conditions, particularly rainfall amounts representing the lower 95% CI of the different climate scenarios, the updrain wetland produced no flooding at all; but as wetland loss increased, the updrain wetland flooded for longer.

Figure 9. (a) Flood volume, (b) duration of subsystem flooding, and (c) duration of upper wetland flooding for different values of wetland loss in a sub-watershed that was modeled to allow infiltration. Wetland loss in this subwatershed occurred to a greater extent than it did in the whole system, but losses are still representative of the present day, four scenario visions, and the business-as-usual scenario. Presence of trendlines indicates significant linear relationship (p [less than] 0.05). R2 values for flood volume versus wetland loss were all ≥0.89 and p values were all [less than]0.01. R2 values for duration of subsystem flooding (when infiltration was allowed) versus wetland loss were all >0.89 and p values were all [less than]0.05. R2 values for duration of upper wetland flooding versus wetland loss were >0.84 and p values were [less than]0.05.

Discussion

Climate change is projected to alter rainfall intensity and timing across the globe (IPCC, 2021). In many cities, these changes will increase the load on stormwater systems— many of which are already underperforming—and render others inadequate. Cities are also facing increasing flood risk due to historical patterns of development and urban designs that proliferate impervious surfaces, which are known to exacerbate pluvial flood risk (Morita, 2014). The combined effects of these phenomena will be a worsening of existing pluvial flood zones or the creation of new ones. To evade such a future, stormwater managers must consider a range of future rainfall and land-cover scenarios and develop stormwater management systems that account for uncertainty.

From the results of this study, we identify arguments in favor of conservation, stewardship, and even perhaps construction, of wetlands to reduce future flooding in cities facing increasing storm intensity. We found that system flood volume and the duration of flooding increased with projected wetland loss under the present-day climate and in the two warming scenarios. These phenomena were present even for the upper confidence intervals of precipitation, which in the case of RCP 8.5 represented a 278% increase above the mean rainfall for a 100-year, 24-hr duration storm event in 2021. We did not find an exhaustion point for the mitigating effects of wetlands for any rainfall input. Further, models did not exhibit an exhaustion point even at the upper bound of wetland loss under the business-as-usual scenario, indicating that the presence of urban wetlands may still provide beneficial flood mitigation even when wetlands are scarce.

Urban wetlands may reduce flooding by removing water through infiltration, but the effect is limited by wetland extent in the system and real-world conditions not captured in our model. We showed in a subwatershed model that the contribution of infiltration toward the overall flood mitigation abilities of wetlands decreased as wetland loss increased, indicating that when wetlands are smaller the contribution of infiltration toward flood mitigation may play a more minor role in flood mitigation compared with other hydrologic processes, such as peak discharge rate reduction and surface storage (Fletcher et al., 2013; Kadykalo & Findlay, 2016). Additionally, because model infiltration proceeds in the EPA SWMM model until the run ends rather than up to soil saturation or soil layer depth, the reduction of flood volume attributed to infiltration in this study may be greater than what occurs in the real world. Further, if extreme storms occurred during wet-season conditions when soils are already waterlogged and the water table is high, infiltration may be negligible. Thus, depending on wetland extent and antecedent weather conditions, the emphasis on the infiltration function of inland urban wetlands in literature on wetland ecosystem services (Chan et al., 2018; Gülbaz & Kazezyilmaz-Alhan, 2014) may bias attention toward a relatively minor source of the flood mitigation services of wetlands.

We generally did not find significant relationships between wetland loss and rates or time of peak of discharge except under the present-day climate scenario. These findings differed with previous studies that found that wetlands reduced rate of peak discharge and delayed arrival of peak discharge in non-urban settings (Kadykalo & Findlay, 2016; Yang et al., 2010). As conduit wetlands in our model were lost, they effectively became more channelized as a result of how we altered channel width and discharge rates were increased, so we would expect the hour of peak discharge to arrive earlier. We posit that the lack of significant change to hour of peak discharge may be explained by other factors such as the location of wetlands with lost areas within the system and system complexity, which may reduce or nullify the effect of wetland loss by shifting flow metrics in other parts of the drainage system.

We then find qualified support for urban inland wetlands reducing flood risk in cities facing future increases in rainfall intensity. Inland urban wetlands in Valdivia reduced systemic flooding under present-day and future conditions where rainfall of 100-year storms was double that of the present day. However, the difference in the rainfall that was converted to flood water among wetland extents was much smaller than the difference in flooding between present-day and future extreme storms. Additionally, the stormwater management system performed worse overall as event rainfall increased, indicating that the flood-regulation services of wetlands are not static even when system dimensions are static. Under real-world conditions in the future in Valdivia, where estimates indicate lower annual rainfall volume, rivers levels may be lower, which may reduce tidal backflow to the SWMM below what we modeled. This may provide some relief to system-wide flooding, particularly in areas near river outfalls, but flooding in the parts of the system that are more inland would likely be close to our estimates.

Further, while wetland cover reduced system flood volume, wetlands in our models were sometimes themselves sources of flooding, depending on wetland extent and rainfall inputs, as we found in the subwatershed model. Although much of the literature presents wetlands as entirely beneficial water regulators (Costanza et al., 1997; C. Wang et al., 2018; Zhang et al., 2020), we argue that this flood generation represents a wetland disservice. We thus are in the minority of studies acknowledging the flood potential of wetlands in cities (Assefa et al., 2021; Lundin, 1994) or urban wetland disservices in general. Additionally, the likelihood of wetlands being sources of nuisance flooding and the amounts of nuisance flooding they may produce likely depend on wetland configuration within the system.

We suggest that cities can produce greater reductions in flooding through the strategic placement of wetlands in their network, diversion of more urban subwatersheds to wetlands, or perhaps by engineering wetlands (e.g., adding channels to increase flow during a storm and to lower water levels preceding a storm or adding gates to hold and release water as needed). Wetlands in Valdivia were mostly formed by an earthquake, and consequently they were placed and shaped in ways that may not be as effective as planners and engineers may be able to achieve. The effects of wetland loss and/or increasing storm intensity on flooding in other cities may be reduced beyond what we demonstrated in our system if our suggestions are followed. However, while these suggestions may augment the ability of wetlands to reduce flooding, our finding of rainfall's outsized effect on flooding nonetheless leads us to recommend urban wetlands as but one GSI tool in a diverse GSI portfolio to manage stormwater under a changing climate. Finally, urban wetlands, even in forms engineered to maximize water regulation, may provide a host of other ecosystem services not afforded by other forms of GSI (Davidson et al., 2019; Millennium Ecosystem Assessment, 2005; Wong et al., 2017; Xu et al., 2020), such as habitat for increased biodiversity, local climate regulation, aquifer recharge, aesthetics, educational opportunities, and spiritual connection (Correa et al., 2018; Millennium Ecosystem Assessment, 2005).

In addition to altering storm intensity, climate change may undermine or enhance the water regulation services of urban wetlands by altering antecedent conditions of the urban system. We may expect that an overall drier future for a city like Valdivia where annual precipitation is expected to decrease may, in turn, increase the available storage of wetlands soils and surfaces. Infiltration could remove more water from the system, and additional surface storage would likely reduce systemic flooding. Conversely, a wetter future for a city like Valdivia may reduce available wetland storage and thus flood mitigation. In wet seasons, wetlands may be identified as sites of nuisance flooding even if their effect throughout the whole year is to reduce net systemic flooding.

Our study is unique in estimating the changes in ecosystem services of wetlands under projections of urban development because we modeled the effects of land-cover change in scenarios that were developed by stakeholders in the city rather than simply by a team of scientists (Liu et al., 2021). In particular, we used the sustainable future scenarios framework (Iwaniec et al., 2020), which were based on positive visions of the future as an alternative to ones that simply occur according to historic precedent. Flood exposure in the scenarios was less than when development continued only according to historical trend, but further, the results of this work are intended to feed into a new round of stakeholder decision-making, where they can evaluate the risks and tradeoffs in the different scenarios and revise the scenario objectives. For example, based on our findings, Valdivian citizens may want to enforce wetland protections earlier in the Eco-wetland and Resilient-to-flood scenarios to ensure greater conservation of wetlands, or even work to expand wetland areas to counteract historical trends of wetland loss. Cities develop not simply by inertia but rather by human agency, and resilience research should work to account for and promote agency in the process.

Finally, we review the utility and limitations of computational modeling in assessing the effectiveness of urban GSI for reducing pluvial flood risk in a future impacted by climate change. In our work, we struggled with relating changes in wetland cover to changes in the 1-dimensional conduit conception of stormwater management systems, and we could not model with spatial accuracy the effects of rainfall on wetland storage units or conduits. We did not account for seasonal- or event-based rates of evapotranspiration, which may be another major source of removal of water from our urban system. Further, our model was limited in its original calibration by lack of consideration of storms with return periods greater than 24 years, which others have argued is a major source of inaccuracy predicting flood volumes and rates of peak discharge (Barco et al., 2008). Long-term data collection is necessary to increase the odds of collecting system response to storms of large magnitudes, and this is complicated by estimates of increasing magnitude of future storms. The city also lacked data on water table location and soil conditions to simulate region- and season-appropriate levels of infiltration. In our experience, cities rarely employ multi-dimensional models that incorporate both land cover and drainage systems or have such data on soil and evapotranspiration. Valdivia's unmodified model did not include infiltration, which, according to our subsystem model, may nullify the relationships between wetland loss and system flooding. To develop strategies to reduce the threat of flooding that future extreme storms evidently pose, such models and their supporting data, are critical to researchers and cities alike for building, maintaining, and improving resilient urban spaces in the face of climate change.

Conclusion

In this study, we simulated using EPA SWMM the effects of future climate change and land cover change on a city-wide urban stormwater management system containing a high proportion of wetland cover. We manipulated model rainfall to represent extreme storms (100-year return period/24-hr duration) in the present-day and by 2080 given the projections of downscaled regional climate models. We used scenario workshop products in cellular automata-based models to estimate future land cover in four intentional development scenarios and one business-as-usual scenario where landcover change proceeded as it had historically. Flood volumes and durations increased with wetland loss and with more extreme storms, though scenario workshop visions typically conserved more wetland cover than did the business-as-usual scenario, thereby proving useful in mitigating some flooding. We also found that the stormwater management system performed worse with more extreme storms under all land cover scenarios, highlighting the threat that climate change alone may pose for cities—even ones with extensive wetland cover. We demonstrated the diminishing contributions of infiltration of wetlands to flood volume and flood duration, indicating that mechanisms related to flow and storage in urban wetlands may contribute more than infiltration to water regulation perhaps unless urban wetland cover is expansive. Finally, we demonstrated a case where a wetland was a source of flooding under extreme storms and climate change. We found that while the conservation of urban wetland cover may reduce pluvial flood risk, we nonetheless caution that conserving or even expanding urban wetland cover is unlikely to resolve the threat posed by the effect of climate change on future extreme storms. More monitoring data are needed for urban stormwater management systems with wetland cover, particularly during extreme storms, and hydrologic pathways in wetlands need to be added to future models to improve understanding of contributions of urban wetlands to pluvial flood risk reduction.

Acknowledgments

This study was supported by the Instituto de Ecología & Biodiversidad at the Universidad Austral de Chile; Agencia Nacional de Investigación y Desarrollo de Chile (ANID/BASAL Grant FB210006); UREx SRN (National Science Foundation Grant 1444755); IRES: Interdisciplinary student research on urban resilience in Latin America (National Science Foundation Grant 1658731); NATURA (National Science Foundation Grant 1927468); and the Earth Sciences Division of the National Science Foundation (award no. 2204589).

Conflict of Interest

The authors declare no conflicts of interest relevant to this study.

Data Availability Statement

The data on which this article is based (rasters, shapefiles, EPA SWMM models, and .csv files) and R code used herein are available at Sauer et al. (2024).