Multiple methods confirm wetland restoration improves ecosystem services

ABSTRACT Global wetland loss has reduced biodiversity and ecosystem services. These declines have inspired many landholders to restore wetlands, but the success of these efforts remains unclear, in part, because quantifying ecosystem services requires diverse methods. Here, we blend participatory mapping and surveys, field measurements and high-resolution models to track ecosystem services from restored wetlands on private land. We ask: 1) What ecosystem services do people perceive from restored wetlands? 2) What modelled/field measured ecosystem services were enhanced through restoration? and 3) How do field measured, modelled and perceived ecosystem services in restored wetlands interact? Participating landholders mapped their restoration project and shared their perceptions of ecosystem services. Next, we modelled ecosystem service changes using the Land Use Capability Indicator (LUCI) model and contrasted these to field measured ecosystem services for each wetland. Landholders perceived ~6.5 services from their restored wetlands. For modelled services, restoration significantly enhanced nitrogen and phosphorous retention. For field-measured services, restoration increased soil organic carbon by ~20%, soil permeability to water by ~27% and native plant species richness by ~15 species, while reducing plant-available phosphorous by ~23%. Correlating across methods revealed that reduced plant-available phosphorus and site age and size were associated with more perceived services, whereas an increase in plant species richness was not a good proxy for gains in measured, modelled or perceived services. Based on the diverse ecosystem services gained, demonstrated by multiple methods, we contend that private wetland restoration can be successful as well as leveraged to meet multiple management and policy objectives.


Introduction
Momentum for wetland restoration is growing, in part, because the consequences of wetland loss for ecosystem services are increasingly apparent. Since 1900, global estimates of wetland drainage range from 64% to 71%. This drainage has resulted in losses of many ecosystem services, the benefits people gain from ecosystems (Ausseil et al. 2008;Tomscha et al. 2019), usually in exchange for agricultural productivity. Relative to their area, wetlands are disproportionately important for ecosystem services (Barbier et al. 1997). Wetlands store carbon, trap sediment and retain excess nutrients (thus promoting better freshwater quality), and attenuate floods (Clarkson et al. 2013). To regain the diverse ecosystem services from wetlands, restoration and reconstruction will be key.
While efforts to restore wetlands are widespread (Peters et al. 2015), providing evidence of the ecosystem services regained from wetland restoration remains challenging. The restoration research community has long grappled with how to define and measure wetland restoration success (Ruiz-Jaen and Aide 2005). Approaches have included measuring species diversity, vegetation structure and ecological processes (Ruiz-Jaen and Aide 2005). Although frameworks are developing, no globally applicable methods exist to assess wetland ecosystem services via modelling and/or field measurement (Janse et al. 2019). Tracking restoration gains may encourage more landowners to participate in wetland restoration, underscored by the growing market of wetland mitigation banking (Palmer and Filoso 2009). Additionally, systematic assessment of progress could enable restoration policy to be target towards landscape contexts providing the most favourable outcomes. Finally, evidence of gains in ecosystem services may inspire more government, industry and philanthropic support for what might otherwise be perceived as a costly endeavour.
Ecosystem services are notoriously difficult to measure, because they represent the dynamic flow of benefits from ecosystems to people (Haines-Young and Potschin 2010). To fully represent this social-ecological system, both biophysical and social components of ecosystem services should be measured (Reyers et al. 2013). Furthermore, ecosystem services are often produced and received at different spatial and temporal scales (Holland et al. 2011;Bagstad et al. 2013), and they may be highly variable across space and time (Rau et al. 2018). These attributes mean ecosystem services are difficult to measure in the field. Often, direct measurements are not sensible or even possible. Nonetheless, field-measured indicators are key for estimates of service gains from restoration projects, and a range of biophysical indicators have been developed to facilitate field measurements (Sutherland et al. 2016b;Ziter and Turner 2018). For example, soil organic carbon is an indicator of soil carbon storage, while saturated hydraulic conductivity is an indicator of flood mitigation (Ziter and Turner 2018). Multiple methods, including social science, can complement field measurements to capture the range of ecosystem services gained through restoration. Interviews and surveys may be particularly useful for asking about temporally dynamic process in difficult to access sites. For example, asking people how often and where they hunt is more feasible than waiting for observations of hunters on private land.
Another challenge in studying restored wetlands is identifying their locations. For example, in New Zealand, no centralized database or maps of restored wetlands exist, partly because many wetlands are located on private land (GWRC 2003), where management decisions are considered confidential. Furthermore, wetlands are often small and therefore difficult to map (Gallant 2015). Consequently, costly high resolution, multispectral imagery is often required to adequately map wetlands. High-resolution data also facilitate the detection of soil moisture (Maxa and Bolstad 2009;Adam et al. 2010). Even when wetlands can be detected, distinguishing restored versus intact wetlands may not be feasible if their spectral signatures are indistinguishable. In addition, land-use history and some restoration treatments, such as predator control and fencing, are not visible in satellite imagery. One solution to the problem of locating restored wetlands is participatory mapping; where researchers ask local people to identify restored wetlands and describe associated management activities. This approach is easily paired with surveys and interviews to capture the intangible ecosystem services from wetland restoration.
To estimate ecosystem service delivery over large areas, modelling is necessary. However, heterogeneity in small-scale ecosystem processes is generally unresolved by models that run at the watershed or sub-watershed scale. Modelling changes in ecosystem services due to small restored wetlands is therefore challenging (Tomscha et al. 2019). The Land Use Capability Indicator (LUCI), an extension of Polyscape, is one of the few tools that models multiple ecosystem service indicators at high spatial resolutions (<5 m; Jackson et al. 2013). LUCI is designed to answer questions about land and water restoration, making it ideal for analysis of wetland ecosystem services. The framework can incorporate minute variations in topography and diverse approaches to land management which can affect patterns of nutrient flow at the farm scale (White et al. 2015). For example, LUCI was used to examine which farms would provide the greatest gains in ecosystem services if receiving agricultural subsidies in Wales (Sharps et al. 2017). Importantly, one challenge in robustly estimating ecosystem service gains from wetlands is that the landscape context of each wetland is an essential input. Nutrient exports from upstream areas may change the magnitude of nutrient retention provided by a restored wetland (Trepel and Palmeri 2002). Thus, a spatially explicit model that accounts for upstream contributing areas, such as LUCI, is critical for accurate estimates of ecosystem service delivery.
Here, we use a blend of in-person surveys, geospatial modelling, field measurements and participatory mapping to identify changes in ecosystem services through wetland restoration. Drawing on the strengths of each approach, we ask three main questions: 1) What ecosystem services do people perceive from restored wetlands? 2) What modelled/ measured ecosystem services were enhanced through restoration? 3) How do field measured, modelled and perceived ecosystem services in restored wetlands interact? Using multiple methods to estimate ecosystem service gains from restoration, we create a more comprehensive picture of the social-ecological system that responds to wetland restoration.

Study site
This study was conducted in the Wairarapa region, the southern part of the North Island, New Zealand. The Wairarapa is bounded by the Remutaka mountain range and the Pacific Coast ( Figure 1). The 5938 km 2 region is home to ~45,000 inhabitants, including Māori iwi (tribes) Ngāti Kahungunu ki Wairarapa and Rangitane o Wairarapa. Primary land uses in the region include dairy and beef cattle and sheep farming, viticulture, forestry and conservation. Wetland area in the Wairarapa has declined 98.7% from its precolonial extent, exceeding global loss averages of 71% and the New Zealand loss average of 90% (Ausseil et al. 2008). Of the existing wetlands, 75% are located on private property, and the majority do not exceed 3 ha (GWRC 2003). Influxes of phosphorus and nitrogen have led to eutrophication of waterbodies (Conley et al. 2009), including Lake Wairarapa, the largest lake in the region, which is now classified as supertrophic (GWRC 2012). Within the region, landholders, farmers, iwi and community groups are restoring wetlands on private and public property, often with the goal of improving water quality.

Measuring ecosystem services using multiple methods
We assessed ecosystem services using multiple approaches; participatory mapping and interviews, modelling, and field and lab-based measurements ( Figure 2). Each method addressed different ecosystem services and produced different data structures (e.g. binary, continuous). First, we conducted a participatory mapping and survey exercise with private landholders who had restored a wetland on their property, including questions about their use and values for their restored wetland. The resulting maps of restored wetlands were then used to calculate the upstream contributing area for each wetland. Next, we used LUCI to estimate changes in ecosystem services resulting from restoration. Finally, we conducted field surveys of plant diversity and characterised soil properties in restored wetlands, and  adjacent unrestored land, to estimate differences in ecosystem services achieved through restoration.

Participatory mapping, surveys and site selection
We recruited landholders in the Wairarapa with wetlands -defined as wet paddocks, wet fields and the terrestrial margins of lakes, streams, rivers or pondson their property. We sent flyers advertising our research to all registered farms in the Ruamahanga Basin and received eight direct responses. We also contacted local restoration groups, the regional council and others through our research networks. In total, we interviewed 28 landholders. Twenty-five of these landholders had restored one or more wetland on their property, while the remaining three were interested in doing so but had not yet applied restoration treatments.
Participatory mapping and surveys were completed using Maptionaire, an online map-based survey tool. All mapping and surveys were completed in person by S. Tomscha from June -October 2018. Participants elected to input data independently using a laptop or to be guided in data input by the researcher. The mapping captured the size and location of one restored wetland on each property. Where more than one wetland had been restored, participants selected either their largest or most important restoration project.
We then asked the landholders a series of structured and semi-structured questions related to their values and restoration practices. For example, we asked: Which of the following are uses of the wetland on your property? (see Appendix 1 for survey details). Landholders identified ecosystem services from a list of 22 services, resulting in binary data for 22 ecosystem services per site. The 22 ecosystem services available in the list were selected through an iterative process of input from practice interviews with restoration practitioners and Wairarapa residents.
Next, we selected locations for field sampling and modelling work from among the 28 wetlands. Each site met the following criteria: 1) the property must have a wetland that has received a restoration treatment, and 2) the wetland must exceed 20 × 20 m in size (to allow standard vegetation monitoring). Eight sites did not meet these criteria, and a further two were unable to provide site access. Thus, fieldwork modelling, and the evaluation of landholders' perceptions, was conducted on 18 sites.
For modelling and fieldwork, we applied a paired sampling design, whereby restored wetlands were contrasted to a locally similar, unrestored site (for fieldwork) or to a baseline agricultural condition (for modelling). For fieldwork, locally similar unrestored sites were often located immediately adjacent to the restored wetlands, separated only by a fence (16/18 farms). In all cases, farmers identified a location as similar as possible to their restored wetland site prior to the application of restoration treatments (i.e. similar elevation, slope, aspect, soils, vegetation, etc.).

LUCI modelling
LUCI is a spatial modelling tool that can be used to create maps and measures of ecosystem service provision (Sharps et al. 2017;Trodahl et al. 2017). We used LUCI to quantify ecosystem service changes from wetland restoration, producing a 5 m modelled output. We determined the upstream contributing watershed area for each wetland site, including the wetland polygon itself. Within this watershed, we generated quantitative estimates for three soilhydrological processes: sediment export, nitrogen exported and phosphorus exported; first under restored conditions and secondly under 'no restoration' (i.e. pasture conditions). Calculating the difference between these two scenarios allowed us to estimate the amount of sediment, N and P retained per annum, due to the restoration of the wetland.

Input data
Seven map-layers were obtained from publicly available databases and used as inputs to LUCI, with minor alterations (Table 1). To generate multiple sets of land cover data, we created two edited land cover files by overlaying restoration site maps with the LCBD (New Zealand Land Cover Database) data. In the first, restored sites were classed as either 'Indigenous Forest' (swamp) or 'Herbaceous Wetland Vegetation' (fen) depending on their terrain. In the second, we assumed a 'norestoration' scenario and manually reclassified all sites as 'High Producing Exotic Grassland', the LCBD category which best corresponded to pasture.

Nutrient exports (N and P)
Nutrient losses were calculated using a modified 'export coefficient' approach, where nutrient load inputs are calculated based on soil, climate, topography and land management (Trodahl et al. 2017). Here, total accumulated nutrient exports were calculated by including contributions from runoff at nearby higher elevation sites and accounting for reductions in load where nutrients are routed through mitigating features such as forests and wetlands (Trodahl et al. 2017;Taylor et al. 2018). Here, the model outputs for each watershed were focused on one grid cell of the Digital Elevation Model (DEM); the water-flow exit point of the wetland polygon. The change in nutrient exports at the watershed outflow point, which was measured as the difference between restored wetland and pastureland cover conditions, provided a measure of the nutrient mass retained annually as a result of the restoration projects.

Erosion/sediment retention
A quantitative measure of how sediment retention responded to wetland restoration was obtained using LUCI's implementation of the RUSLE (Revised Universal Soil Loss Equation) algorithm. The algorithm calculates the flow and volume of sediment export in tons km −2 yr −1 according to soil, land cover, topographical and management characteristics (Benavidez et al. 2018). To convert sediment loss to soil retention, loss values in tons yr −1 for pastureland cover were subtracted from those calculated for restored wetlands.

Field-measurements of ecosystem services
Restored wetlands varied in size (0.4 to 33.7 ha with a median of 2.5 ha). Time since restoration treatments were first applied ranged from 6 months and 42 years prior to sampling (median = 9). All restored wetlands were fenced and planted with native or exotic species. However, the combination of additional restoration treatments applied to each wetland varied. These treatments included invasive plant species removal (herbicide or physical removal), pest removal (trapping, baiting, hunting), earthworks, pond creation and/or drain blockage or installation.
Of the 'unrestored' wetland sites (i.e. paddocks) used for comparison, seventeen were used for grazing livestock, while one lay fallow.
Fieldwork was conducted in summer between December 2018 and early February 2019. Within each restored and unrestored wetland pair, we established a 400 m 2 (20 m × 20 m) plot for the sampling of plants and soils ( Figure 3; 36 plots total). Plots were oriented so that one edge ran parallel to the waterbody, while the perpendicular edge followed the increasing elevational gradient of the wetland. If no surface water was visible, plots were oriented parallel to the lowest elevation locally. Within each of the 18 restored and 18 unrestored paired plots, we recorded vascular plants to species. When complete identification was not possible, vascular plants were recorded either to genus or family, noting their status as native or non-native. We also recorded the number of bryophytes species present, but these were not identified to species. The total number of native vascular species sampled per plot was used as a measure of native species richness.
Soils were sampled by coring with a 5 cm PVC pipe to a depth of 10 cm along three parallel transects at, 0 m, 10 m and 20 m within the plot. Nine cores were taken to measure soil organic carbon (SOC) for each plot, totalling  324 across all wetlands (Figure 3). Two further soil samples were collected at the highest and lowest elevations of each plot and used for determination of soil Olsen P (P) and saturated hydraulic conductivity (K s ). Cores were placed on ice before refrigeration at 4°C. Laboratory analysis took place between March and August 2019.

Soil organic carbon (SOC)
We measured soil organic matter (SOM) content using the loss on ignition method (Dean 1974). Samples were dried at 105°C in an oven overnight and sieved to <2 mm. The volume of all material over 2 mm was recorded and subtracted from the total volume of the core. The <2 mm fraction was then weighed to allow for calculations of soil bulk density. Subsamples (5 g) were taken from the <2 mm fraction and placed in a muffle furnace for 4 hours at 550°C, cooled and reweighed (Wright et al. 2008). The mass loss of the soil sample was taken as a measure of SOM. We estimated SOC by dividing SOM by 1.98 when SOC was between 1% and 10%, and by 1.86 when SOC was and >10% (Pribyl 2010).

Phosphorous
Plant-available phosphorous, measured here as Olsen P, is an indicator of phosphorous likely to be lost in surface runoff (McDowell et al. 2001). For Olsen P analysis, soil samples were lyophilised for 48 hours, mechanically crushed and sieved to 2 mm. Coarse roots and other large organic material were discarded, and the volume of lithogenic material over 2 mm was measured. Olsen P was determined on 5 g subsamples by bicarbonate extraction (pH 8.5) and measured with a Molybdenum Blue colorimeter by Hills Laboratory Ltd. (Hamilton, New Zealand).

Saturated hydraulic conductivity
Measuring saturated hydraulic conductivity (K s ) over large areas can be expensive and impracticable. Estimations of K s from empirical relationships (such as pedotransfer functions) regressed against easily obtainable lab measurements, however, offer an acceptable level of insight. To estimate soil K s we used bulk density, particle size measurements and soil organic carbon (Equation 1). The high organic matter content in these soils, required use of a pedotransfer function developed for organic topsoils (Wösten et al. 1999 To find bulk density and soil organic matter measurements, the results from the loss on ignition method were averaged for each lower and upper transects of each plot, from the cores used for SOC measurements. For particle size analysis, 3 g subsamples were taken from the two cores at high and low elevations in the plot. These subsamples were lyophilised for 48 hours and then sieved to 2 mm. First, we removed organic material by chemical digestion with H 2 O 2 . Samples taken from transects containing <10% organic material were digested in 27% H 2 O 2 until the reaction ceased. Samples taken from transects containing >10% organic material were digested in 52% H 2 O 2 until the reaction ceased. The samples were then neutralised and lyophilised. Particles were disaggregated in 0.5% calgon solution in an ultrasonicator for 20 minutes. We determined the size of the sand, silt and clay fractions of subsamples using a laser particle sizer (Beckman Coulter LS13320), following USDA/FAO classifications (FAO 2006).

1) What ecosystem services do people perceive from restored wetlands?
To characterize the co-occurrence of perceived ecosystem services from restored wetlands, we calculated a Gower's distance matrix between sites based on the number of perceived ecosystem services they shared divided by the total number of services assessed (here 22) (Gower 1971, McCune et al. 2002. The cooccurrence patterns were then visualised in two dimensions using non-metric multidimensional scaling (NMS) ordination based on Gower's distance. NMS was run using PC-ORD version 7.08 (McCune and Mefford 2018) in the 'auto-pilot' mode, which performed 500 iterations (250 runs each of real and randomized data) with random starting configurations and assessed dimensionality by minimizing stress.
We then applied Ward's method of hierarchical cluster analysis of the distance matrix to examine how perceived ecosystem services co-occurred on sites. This identified four distinct clusters of 'perceived ecosystem services', which we coded as a categorical factor in subsequent analyses. The statistical significance of the perceived ecosystem services clusters was assessed by calculating chance-corrected withingroup agreement (A) using non-parametric Multiple Response Permutation Procedure (MRPP; Zimmerman and Goetz 1985) with Euclidean distance in PC-ORD. The chance-corrected within-group agreement (A) is a measure of within-group homogeneity compared with that expected by chance, where A = 1 corresponds to identical members within each given group (maximum effect of categorical factor), and where A ≤ 0 corresponds to within-group heterogeneity equal to or larger than that expected by chance (no effect of categorical factor (McCune et al. 2002)).

2) What modelled/measured ecosystem services were enhanced through restoration?
We standardised modelled estimates of phosphorus, nitrogen and sediment retention by the area of each restored wetland, and significant differences between restored and unrestored wetlands were determine via paired t-test. Data were log transformed so that they met the assumption of normality.
We tested field-measured differences in restored versus unrestored wetlands using either paired t-tests (for native plant species richness) or linear mixed-effects models (for all soil variables), because the latter had within-plot replication of samples. In the SOC, Olsen P and Ks models, distance from water's edge and restoration state (restored vs. unrestored) were additive fixed effects. Site identity and plot nested within site were included as random variables affecting the intercept only. Plots of residuals vs. fitted values were inspected visually for normality and homoscedasticity, and a nonconstant error variance test from the car package was used to check homoscedasticity among levels of categorical variables (Fox and Weisberg 2019). We applied ANOVA with a type 3 Wald chi-square test to partition the variance associated with fixed factors. To test if restoration had a significant effect on native species richness, a paired t-test was used. All response variables were non-normal and were therefore log-transformed to conform to the assumptions of the linear mixed-effects models and t-test. All analyses were conducted in R studio (R Core Team 2019), using the package lme4 for linear mixed-effects models (Bates et al. 2015).

3) How do field measured, modelled and perceived ecosystem services in restored wetlands interact?
To understand how measured, modelled and perceived ecosystem services of restored wetlands interact, first, we visualised how sites clustered according to the perceived ecosystem services that landholders held for their restored wetlands. To this end, we transformed the initial matrix so that it contained 18 sites x 22 perceived values. We then calculated a Gower's distance between sites based on the number of perceived ecosystem services shared divided by the total number of perceived ecosystem services assessed (here 22) (Gower 1971, McCune et al. 2002. Next, we used NMS ordination to visualize sites in 'perceived ecosystem service space'. This step resulted in a high proportion of variance explained by the first two ordination axes. Then, we calculated Pearson correlations among each of the individual modelled and measured ecosystem services and the significant ordination axes of perceived ecosystem services. We also calculated Pearson correlations among the significant ordination axes of perceived ecosystem services and the number of years since restoration began, restoration site size, upstream contributing areas, the ratio of restored area to upstream contributing area and the total count of perceived ecosystem services at each site. We visualized the significant relationships among measured or modelled ecosystem services and the ordination axes of perceived ecosystem services as joint-plots, which were overlaid on the NMS in PC-ORD. The angle and length of the joint-plots lines describe the direction and strength of these correlations.

1) What ecosystem services do people perceive from restored wetlands?
Landholders used their restored wetlands in different ways and had a variety of perspectives on the ecosystem services they provided. Of the 22 possible ecosystem services of their restored wetlands (Figure 4), landholders perceived a median of 6.5 (IQR = 5.25-9.75). The maximum number of ecosystem services perceived was 21, while the minimum number was two ( Figure 5). The most commonly cited ecosystem service of wetland restoration was wildlife habitat or biodiversity, and the second most commonly identified ecosystem services was aesthetics or beauty. No landholders reported gathering fibres as an ecosystem service.
The NMS ordination revealed that 88% of the variation in landholder-perceived ecosystem services of restored wetlands could be described by two NMS axes ( Figure 5). Moreover, we found statistical support by MRPP for four significant clusters of ecosystem services (A = 0.222, p ≪ 0.001) ( Figure 5). The majority (77%) of variation in landholder-perceived ecosystem services for wetlands encompassed the transition from very tangible, mostly provisioning services which clustered high on NMS axis 1 to abstract/intangible ecosystem services for aesthetics, relaxation, history and habitat for wildlife ranked low on axis 1. While axis 2 accounted for only 11% of the variation in landownerperceived ecosystem services it helped to differentiate two additional ecosystem service clusters. The third cluster, which ranked low on axis 2 includes ecosystem services that are potentially more individually enjoyed (e.g. family, photography and space for bee foraging) versus the final cluster of ecosystem services which provide collective services, potentially for those residing off the landowner's property (e.g. sediment and nutrient retention, education and hiking).

Modelled gains in ecosystem services
The difference between the amount of nitrogen, phosphorous or sediment exported in restored versus unrestored wetlands can be interpreted as the amount retained per annum through restoration. Upstream contributing areas of the restored wetlands varied considerably (4.2 ha to 2263.3 ha) (Table 2); therefore, we standardized our ecosystem service data by restored area to increase comparability across sites. Consistent gains in modelled ecosystem service were found with all sites benefiting from restoration for all three indicators.

Nitrogen retention
Biophysical models implemented in LUCI indicated lower total nitrogen exports in restored wetlands ( Figure 6). Restored wetlands exported a mean of 148 ± 43 kg N ha −1 yr −1 , while unrestored wetlands exported a mean of 166 ± 48 kg N ha −1 yr −1 , and this

Phosphorus
Olsen P, an indicator of downstream water quality, varied greatly by site. Soils in restored wetlands contained between 6 and 62.5 µg P cm-3 dry soil, with an average of 25 µg P cm-3 dry soil. Unrestored wetlands contained between 11 and 51.5 µg P cm-3 dry soil of Olsen P, with an average of 28 mg/L (Figure 7). Overall, the linear mixed effect model showed that restoration decreased Olsen P, by 22.95% ± 13% (χ 2 (1) = 4.36, p = 0.036). Distance from water did not affect Olsen P quantities (χ 2 (1) = 0.030, p = 0.862). Olsen P was higher in restored wetlands than in unrestored wetlands for 10 of 18 sites and had no effect at one site (Figure 8).

Saturated hydraulic conductivity
Saturated hydraulic conductivity, an indicator of flood mitigation, ranged from 0.41 to 2.36 mm hr −1 in restored and unrestored wetlands with lower  values having a lower capacity to attenuate floods. Pedofunction-derived saturated hydraulic conductivity is reliable to ~1.02-1.67 mm hr −1 (RMSE) (Tóth et al. 2015). As such, the relative uncertainty of each Ks value is as much as 100%; however, inclusion of organic matter tends to decrease RMSE. This means that between site comparisons in Ks are tentative; however, we have high confidence in the restoration comparison because of the similar soil characteristics between restored and unrestored wetlands within each site. The average K s of restored wetlands was 1.239 mm hr −1 , whilst unrestored wetlands average K s was 0.97 mm hr −1 (Figure 7). The linear mixed model showed that restoration increased a wetland soil's saturated hydrologic conductivity by 27.3% ± 11% (χ 2 (1) = 5.58, p = 0.018). Distance from water did not significantly change soil Ks (χ 2 (1) = 2.87, p = 0.090). Restoration improved K s relative to paired unrestored wetlands for 14 of 18 wetlands (Figure 8).

Native plant species richness
Restoration had a significant effect on species richness (t (17) = −8.72, p < 0.0001), with restored plots having a mean richness of 17.9 species, and  unrestored plots having a mean richness of 2.6 species (Figure 7). In one case, the restored wetland plot contained 37 more native plant species than the unrestored plot. Restoration increased the native plant richness in 17 of 18 wetlands (Figure 8).

Total gains in measured ecosystem service indicators
Of the four ecosystem service indicators measured in the field (SOC, K s , Olsen P and native plant species richness), wetland restoration increased at least two ecosystem services indicators for all sites, and restoration never resulted in a net loss in ecosystem services. For four of 18 sites, all four ecosystem services increased. A further nine sites showed a net gain in three of the four measured ecosystem services. Five sites were net neutral, showing gains in two ecosystem services and losses in two others (Figure 8).

3) How do field measured, modelled and perceived ecosystem services in restored wetlands interact?
The NMS ordination of the perceived ecosystem services by site showed that 79.2% of the variation in the data could be captured in two axes (Figure 9). Axis 1, which explained more than 60% of variability in perceived ecosystem services was strongly driven by the total number of perceived ecosystem services at each site, as revealed by the very strong Pearson correlation between axis 1 and the count of perceived ecosystem services (r = 0.972, p ≪ 0.001). We also found a weaker correlation between the number of years that a site had been restored and axis 1, although this was only statistically significant at P < 0.1 (r = 0.403, p = 0.097). Correlating the measured and modelled ecosystem services with Axis 1 of landownerperceived services revealed that it was significantly correlated with the change in plant-available P due to restoration (Table 3, Figure 9). We also found Axis 1 was significantly correlated with upstream contributing area and restored area (Table 3). Axis 2 generally corresponds to the split between cluster 2 and cluster 3 services; with landowners of sites p and n reporting all three of the cluster 2 ecosystem services, while sites i and k reported none of the three ecosystem services in cluster 2, but two or more of the ecosystem services in cluster 3. In fact, field measurements of SOC were significantly correlated to axis 2, along which sites rich in collective ecosystem services (e.g. sediment and nutrient retention and education and hiking) separate from sites rich in individually experienced benefits. The association between SOC and cluster 3 values is germane as soil C storage is a benefit that is collectively enjoyed, even if it may be not always be perceived by landowners engaged in wetland restoration.

Discussion
Wetland restoration enhanced ecosystem services on private land as demonstrated by field measurements, modelling and landowners' perceptions. Our approach quantified both the scope and magnitude of ecosystem service gains, using multiple methods to gain a holistic perspective of restoration outcomes. Although the median size of the restored wetlands was only 2.5 ha, all restored wetlands increased overall biodiversity or ecosystem services, which shows promise for landholders interested in small restoration projects. Nonetheless, as restored wetland size and their upstream areas increased, landholders identified a higher number of perceived services. These associations suggest that landholders who perceive more wetland services are inspired to restore larger wetlands. The strong positive correlation between perceived ecosystem services and upstream contributing areas suggests landholders may be aware of how landscape position affects the ecosystem services provided by their wetlands. From a biophysical perspective, these perceptions are supported by the significantly negative correlation between plantavailable P and wetland size, which suggests that larger restored wetlands may be more important for improving water quality. Green tiles indicate this ecosystem service increased following restoration at this site, while light purple tiles indicate that this ecosystem service decreased. White tiles mean this ecosystem service remained neutral. Soil organic carbon (SOC), saturated hydraulic conductivity (Ks) and Olsen P were variable by site, where as plant species richness increased at all but one site.
Landholders differed primarily by the total number of ecosystem services they perceived, rather than by their perception of individual or clusters of ecosystem services. These shared perceptions suggest broad commonalities among landholders engaged in wetland restoration. The most commonly perceived ecosystem service was biodiversity (17/18 landholders) followed by aesthetics (15/18 landholders). The perception of wetlands as beautiful signals landholder opinions of wetlands have shifted, as historically, many farmers perceived wetlands as unattractive wastelands (Nassauer 2004). In contrast, no landholders perceived gathering fibres as a service from their wetlands, suggesting privately restored wetlands are not providing this service. In New Zealand, harakeke (Phormium tenax) is a wetland plant traditionally used by Māori for weaving. People who would gather these fibres may not have access to these privately held wetlands, which underscores that public restoration projects are still key for maintaining this important cultural resource.
Ecosystem services often co-occur in clusters, sometimes called bundles. Typically, these clusters are defined across diverse land uses and stakeholders (e.g. Raudsepp-Hearne et al. 2010;Queiroz et al. 2015). Our approach is unique in that we clustered ecosystem services from one ecosystem (restored wetlands), among a narrow group of stakeholders (landholders with restored wetlands). Within this narrow group of landholders, perceived services separated into four clusters: tangible, intangible, collective and individual ecosystem services. While there are justifiably numerous frameworks to categorize ecosystem services (e.g. MA 2005; Costanza 2008), our results give insight into how current classifications align with clustered perceptions of ecosystem services. For example, our clusters showed recreational services were more tightly tied with provisioning services, rather than to broader cultural ecosystem services, despite them being generally categorized as cultural (MA 2005). In contrast, collective vs. individualistic services reflected the degree to which these ecosystem services are used exclusively by landholders, aligning well with ecosystem service excludability categorizations (Costanza 2008). Discerning clustered patterns among socio-cultural preferences for ecosystem services can facilitate better management for their holistic provision (Martín-López et al. 2012;Mouchet et al. 2014).
Ecosystem service interactions are often measured using correlations (e.g. Raudsepp-Hearne et al. 2010). However, assessments of interactions across measurement methods are rare (Mouchet et al. 2014). By assessing these relationships, we found interactions among perceived and measured services, as well as among perceived services and site attributes. Site attributes, such as the upstream contributing area, the area of wetland restored and restoration age are relatively easy to measure and all correlated with the major axis of variation in perceived services (axis 1). Likewise, the Figure 9. NMS ordination of sites in 'perceived ecosystem service' space. Axis 1 represents 60.5% of variability in perceived ecosystem services, while axis 2 shows 18.7% of this variation. Symbol size indicates the total count of perceived values per site. Significant correlations among axis one and two and measured and modelled ecosystem services indictors overlaid onto the ordination as vectors, dashed lines indicate P < 0.1, solid lines indicate P < 0.05. change in plant-available P with restoration was correlated to axis 1. These associations, together with the shared perceptions of landholders engaged in wetland restoration, suggest that people may prioritise the restoration of wetlands with large contributing areas. Over time, as trees replace pasture plant species, these restored wetlands support greater plant biomass, which enables greater uptake of phosphorus (Chapin et al. 2011). However, because our analysis of ecosystem service interactions is based on correlation, we cannot determine causality. In contrast to the suggestion that large wetlands and/or large wetland contributing area sizes are explicitly prioritised, it is possible that as landholders become increasingly aware of the ecosystem services their wetlands provide over time they become motivated to restore greater wetland areas. Most probably, some combination is true. We contend that perceived ecosystem services from wetlands, which are shaped by landscape context, evolve as the cultural context and biophysical ecosystem services of restored wetlands develop over time.

Relationships among biodiversity and ecosystem services in a restored wetland context
Native plant species richness increased in all but one of the restored wetlands. However, richness was not significantly correlated with any other measured or modelled service, nor was it correlated with perceived services. A lack of correlation among biodiversity and ecosystem services has been found in a variety of other systems, and relationships among biodiversity and ecosystem services remain complex and unclear (Ricketts et al. 2016). Furthermore, understanding these relationships is challenging, in part, because diverse methods are used to measure both biodiversity and ecosystem services (Harrison et al. 2014). Relationships among biodiversity vary widely depending on the biodiversity metrics and ecosystem services considered. Using spatial correlations, many studies have associated modelled ecosystem service to maps of species counts (Watson et al. 2020), but equally, biodiversity proxies, such as functional diversity and/or stand age, may also be correlated with ecosystem service delivery (Harrison et al. 2014). Our study represents a highly localized example of data collection for ecosystem services and biodiversity, e.g. measurements were taken in the same location, and thus have the benefit of not being subject to spurious correlations which may result when data are pooled across large spatial areas. While native plant species richness was uncorrelated to other service gains accrued through restoration in our study, its nearly universal increase suggests that wetland restoration on private land will, at a minimum, be effective in promoting biodiversity. All but one landholder perceived this to be true of their restored wetland. Table 3. Pearson correlation among the first two ordination axes of perceived values from the NDS shown in Figure 9 and measured and modelled ecosystem services in this study, significant correlations at P < 0.1 are underlined, significant correlations at P < 0.05 are in bold type. Asterisks indicate significance in models as follows: *p < 0.05, ** p < 0.01, *** p < 0.001, **** p <.0001.

Site attributes
Axis

Site differences may obscure interactions among modelled, measured and perceived ecosystem services
Apart from a strong positive correlation among modelled P and N retention, modelled and measured services were not significantly correlated. This lack of correlations suggests that none of these measurements are good proxies for one another, nor can we expect their changes to occur in tandem. Furthermore, despite the associations among perceived services with measured services and restoration age, we did not find an association between restoration age and measured or modelled services, contrasting results from other studies which highlight stand age as important for ecosystem services (Sutherland et al. 2016a). This lack of correlation may be because restoration age is not equivalent to stand age or successional stage in our study. The initial conditions of our wetlands varied among sites, as did the range of management techniques applied to achieve restoration. Some sites were initially degraded wetlands, while others were pastures. Furthermore, some restoration techniques are more disruptive, such as earthworks, while others are less intrusive, such as fencing. Exploring these patterns across wetlands with more similar initial conditions and more standardised array of management regimes may elucidate some of these contrasting patterns.

Leveraging multi-method approaches for future benefit
Multi-method approaches can be leveraged to create better ecosystem service models. Ecosystem service models are faster and cheaper than field measurements. Similarly, rapid assessments of ecosystem services are popular for extrapolating site-level gains in ecosystem services from restoration, including economic gains (Peh et al. 2014;McInnes and Everard 2017). However, rapid assessments and models rely on quantitative ecosystem service estimates that are already published or qualitative assessments from stakeholders or experts. Nonetheless, models are parameterized based on field measurements, and their accuracy depends on the availability and quality of these field measurement data. By collecting primary data through field measurements, our work facilitates the adoption of better models. For example, accurate flood mitigation modelling is hampered by a lack of soil saturated hydraulic conductivity measurements (Ebel and Martin 2017), and our saturated hydraulic conductivity data, which includes high-resolution measurement of its variability across small wetlands may be used to supplement flood models. Similarly, spatially explicit carbon models are often based on national or global scale soil tables, but our work now provides better estimates of mean SOC, as well as the spatial variability of these quantities. Furthermore, Olsen P is used as an input to many spatially explicit phosphorus models, but measurements of Olsen P in wetland-specific sites are rare. Future work could leverage our measured data to improve ecosystem service models.

Meeting needs for the greater good on private land
Acknowledging multiple ecosystem service gains from wetland restoration can illuminate interlinkages that bolster the role of wetlands in currently disparate policy objectives. For example, soils from restored wetlands may be leveraged to mitigate greenhouse gasses. In New Zealand, only afforestation has been recognised under the emissions trading scheme. Consequently, landholders are not credited for the carbon gains in restored wetland soils. Concurrently, we demonstrated that by restoring wetlands, N and P are less likely to leach into nearby lakes and streams, resulting in improved water quality downstream. Finally, we add to evidence that wetland restoration increased soil permeability, and thus, the capacity to abate floods (Ameli and Creed 2019), indicating that wetland restoration may provide a promising alternative to engineered flood protection. Taken together, wetland restoration may reduce landatmosphere carbon emissions, improve freshwater quality and provide a promising and cost-effective alternative to engineered solutions for flooding mitigation. Despite the solutions wetlands can provide, current policies overlook wetland restoration as a key tool.
Elevating the role of wetlands in policy is especially urgent as the current scale of wetland restoration does not match the magnitude of their continued loss in New Zealand (Robertson et al. 2019).

Conclusion
Wetland restoration on private land provided significant gains in SOC, native plant species richness and saturated hydraulic conductivity, and also achieved desirable declines in plant-available P. Restoration also led to significant gains in N and P retention and minor gains in sediment retention. Additionally, landowners perceived their restored wetlands to provide a median of six ecosystem services, which in some cases, were correlated with measured gains in services. Using multiple methods allowed us to capture a wide range of ecosystem services and how they interacted, providing comprehensive and quantitative estimates of the gains received when wetlands are restored on private lands. The overwhelming evidence of significant gains in multiple ecosystem services from wetland restoration upholds a need for continued restoration efforts.
Harnessing the potential of wetland restoration may be key for meeting a range of policy objectives in New Zealand and globally.