Modeling Wetland Functions: Is Space-to-Time Substitution of the Perimeter–Area Relationship Appropriate?

Abstract

Wetlands’ morphometric or shape properties, such as their area and perimeter, impact a multitude of ecosystem functions and services. However, current models used to quantify these functions often only use area as an independent variable, as the static area and perimeter of different wetlands have been found to be closely related. The study uses monthly inundation maps, derived from remote sensing data, to assess the temporal covariation of geographically isolated wetlands’ perimeter and surface area. The results show that using static representations of wetlands to evaluate temporal dynamic perimeter–area relationships can introduce significant discrepancies and that these discrepancies can be reduced if evaluations using static data are performed separately for each wetlandscape. This study concludes that models that use implicit area–perimeter relationships based on static wetland representations, as is usually the case, should be applied with caution. Additionally, it suggests that incorporating perimeter–area relationships from temporally dynamic data can improve estimates of wetland functions.

1. Introduction

Wetlands provide a wide range of invaluable ecosystem services, including contaminant buffering [1], flood protection [2,3], biodiversity support [4], and carbon sequestration [5]. Notably, several of these ecosystem services are strongly related to wetlands’ morphometric or shape attributes. For example, the removal rate of contaminants is related to wetlands’ inundation area and perimeter [6,7,8,9,10], the volume of methane emissions is affected by wetlands’ area [11,12,13,14], and the suitability of aquatic habitats for wetland bird communities is influenced by wetlands’ area and perimeter-to-area ratio [15,16,17,18]. Unsurprisingly, the effect of a wetland’s perimeter and area on its functions is a topic of active study. For example, several studies have focused on assessing habitat suitability [19,20,21,22,23,24,25,26,27,28,29,30,31,32] and erosion and water quality [33,34,35,36,37,38,39,40] vis à vis wetlands’ shape attributes.

The effect of a wetland’s perimeter and area on its function is of particular significance for geographically isolated wetlands (GIWs, hereafter). GIWs usually have a relatively small inundation area but higher perimeter length per unit area, characteristics that make them prominent control points for nutrient removal [41] and other biogeochemical processing at the landscape scale [42,43,44]. Even though these wetlands play an important role in modulating downstream water quality, they have been lost at higher rates [45], in part because they are often accorded limited protections [46], although variably so in different parts of the US [47].

Models quantifying wetland functions in relation to wetlands’ morphometric properties oftentimes only use area as an independent variable [7,44,48,49], in part because the static area and perimeter of different wetlands can be closely related [50]. However, it remains unclear if such a relation between static area and perimeter is valid temporally, as the area and perimeter of wetlands vary with time. Addressing this knowledge gap is essential, as a refined understanding of how wetland area and perimeter covary over time could significantly enhance the precision of models used to quantify ecosystem services. To this end, this study sheds light on the dynamic landscape patterns of GIWs. Specifically, we examine the intricate interplay between perimeter and inundation area, employing diverse shape metrics. The results of our investigation provide insights spanning wetlands of varying sizes, diverse wetness conditions, and an array of wetland landscapes. The goal is to unravel the nuanced relationships between wetland dimensions and their functions, paving the way for more accurate predictive models, and optimal designs of constructed and restored wetlands.

2. Materials and Methods

2.1. Metrics Used to Track the Perimeter–Area Covariation

2.1.1. Perimeter–Area (P:A) Ratio

The P:A ratio is a metric often used to capture a water body’s boundary configuration vis à vis its size. Increasing deviation from a perfectly circular shape results in a higher P:A ratio. In other words, an irregularly shaped wetland must be larger than a circular wetland to have the same P:A ratio value. This metric has been used for quantifying the suitability of water bodies to support biotic habitats [51,52] and to assess the potential of wetlands for sequestering contaminants in the boundary ecotone per unit wetland size [42]. A high P:A ratio is typically associated with wetlands that are hydrologically more dynamic and ecologically diverse, such as marshes and swamps. In contrast, a low P:A ratio is usually associated with wetlands that are more stable and less diverse, such as ponds and lakes. The P:A ratio is calculated as:

$$P:Aratio=\frac{P}{A}$$

(1)

where $P$ is the wetland perimeter in m and $A$ is the wetland area in m². The metric has the advantage of simple calculation, intuitive understanding, and representation. As the metric has a dimension of 1/L where L is a length dimension, its magnitude decreases monotonically with increasing area [42].

2.1.2. Shoreline Fractal Dimension

The shoreline fractal dimension (SFD) was originally developed to quantify the degree of irregularity in oceanic coastlines [53,54,55]. Recently, SFD has also been used to characterize wetlandscapes in US [56]. This metric is a measure of the complexity of the shoreline of a wetland, and is calculated using the following equation:

$$P=k{A}^{\frac{D}{2}}$$

(2)

where $k$ is a proportionality factor, and $D$ is the SFD. As the geometric shape changes from a circle to a line, $D$ changes from $D=1$ for a perfect circle with $P=k\sqrt{A}$ to $D=2$ for a line with $P=kA$. A higher SFD indicates a more complex wetland boundary with enhanced twists and indentations. Given that complex shorelines are known to provide more microhabitats, such as different water depths and water flow patterns, a higher SFD is usually considered to favor higher diversity of plant and animal species found in the wetland. Higher SFD for a wetland is also often interpreted to have larger reactive interface length per unit area, thus highlighting the potential for contaminant buffering. For the evaluation of D using static data of the area and perimeter, $D$ is calculated as the slope of the regression line between $\log P$ and $\log A$ of the wetlands under consideration. In order to evaluate the temporal dynamics of $D$, it is computed using the equation $D=2·\log P/\log A$, assuming $k$ as 1, for each time step. Notably, SFD, just like the P:A ratio, is also a scale-dependent metric.

2.1.3. Shoreline Irregularity

The shoreline irregularity (SI) is the ratio of the perimeter of wetland to the perimeter of a circle of the same area [57,58]. The SI is calculated using the following equation:

$$SI=\frac{P}{2\sqrt{\pi A}}$$

(3)

The smallest wetland area in remote sensing data corresponding to a wetland represented by a pixel has an SI of 1.13, while a perfect circle of the same area has an SI of 1. A crenulated wetland will have a larger SI. This metric has been used for evaluating human impacts on wetlands [45] or explaining wetlands’ functions within their landscapes [16,23,26]. Notably, unlike the two other metrics considered above, i.e., P:A ratio and SFD, SI is dimensionless in nature. However, it should be noted that a recent study highlighted that despite the dimensionless nature of this metric, data resolution may still alter its estimates [59].

2.2. Tracking the Dynamics of Perimeter-Area Covariation

To assess the dynamics of the shape attributes of GIWs, the three aforementioned metrics were evaluated for each time instant that the data existed. Next, variations in these metrics in relation to the wetness were assessed. Here, the wetness at a given time was quantified according to the fractional area of the wetland with respect to its maximal extent in the data period. The evaluation was performed separately for wetlands of different size ranges. For all wetlands with their time-averaged area falling within a size range, metrics were averaged for 20 wetness bins (i.e., 0–5%, 5–10%, etc.) to obtain the dynamical relation of shape metric with wetness.

2.3. Assessing the Representativeness of Dynamic Perimeter-Area Covariation Based on Static Wetland Data

To assess the applicability of perimeter–area relationships obtained from static wetland data for defining temporal dynamics of shape metrics, the static shape-based relationship between perimeter and area was first obtained. To this end, the area and perimeter for all wetlands at their respective maximum extent were derived. The reason for using the maximum extent was because it is consistent with the static relations derived using NWI or similar databases in previous studies. Notably, wetlands in NWI are generally developed based on vegetation, soil, and morphological descriptors that often respond to the maximum inundated area extent [56,60]. This is further confirmed by the 1:1 relation between NWI GIWs and the maximum area extent of remote sensing derived wetlands (Figure S1). Using the area and perimeter of the wetlands’ maximum extent based on remote sensing data, instead of evaluating them from static NWI data, also ensured that the methodology can be applied globally.

Once the area and perimeter of each wetland at its maximum extent were established, averages of the perimeter of all wetlands of identical sizes were obtained. The static area–perimeter relationship was then derived as the linear interpolation between each wetland pair ordered based on increasing area. To obtain the predicted perimeter, and consequently the three aforementioned shape metrics, at a given time, the area for that given time was substituted in the static relation. The differences between the predicted and observed perimeter or shape metric magnitude were quantified as percentage discrepancies. These differences were calculated as the ratio of abs(actual shape metric at any given time instant − estimated shape metric using static data relation)/(estimated shape metric using static data relation). Here, abs() indicates the absolute value. The evaluation was performed for each time step. Next, the discrepancies were time averaged.

2.4. Data

Spatio-temporal variations in wetland extent are usually mapped either by implementing distributed or spatially explicit hydrologic models [61,62,63,64] or post-processing of remote sensing data [65,66,67,68]. For mapping over larger areas (say, at continental scales) at high spatial-resolution, remote sensing data have proven to be particularly effective. For example, the Global Surface Water (GSW) [69] and Dynamic Surface Water Extent (DSWE) [70,71] maps have been used to study inundation dynamics at as fine as 30 m resolution [72,73,74,75] over large areas. Here, we used the monthly inundation maps of GIWs, as derived from remotely sensed data in [75], for ten wetlandscapes in the continental US (Figure 1). Using the GSW v1.0 inundation data from March 1985 to October 2015, the GIWs were derived within 1000 km² area rectangular regions within each wetlandscape. The selected regions were identified to ensure the presence of a sufficient number of wetlands in the study area and the frequent availability of remote sensing data. The method to derive the GIWs in the remotely sensed GSW data followed a two-step process. First, GIWs were identified in the National Wetlands Inventory (NWI) dataset following the rubric presented in [76]. This encompassed several steps including the exclusion of riverine, marine, and estuarine wetlands from ensuing analysis, as they are deemed to have surface connections with water bodies. From among the remaining palustrine and lacustrine wetlands within the study area, we identified the ones that have a higher likelihood of connections with other water bodies during wet periods. To accomplish this, proximity filters were used. Specifically, palustrine and lacustrine wetlands falling within buffer polygons with a radius of 10 m around rivers or lakes with size exceeding 80,000 m², or other water bodies such as reservoirs with area larger than 15,000 m², or water features such as bays/inlets, locks, levees, etc., were excluded. Water polygons from the National Hydrography Dataset (NHD) were used in the analysis. These steps yielded a map of GIWs in the NWI data. Notably, several of these GIWS did not include any wet pixels based on the GSW data. The reason for this discrepancy could be due to inherent uncertainties in these data, or mismatch between the periods when these data were generated. Next, GIWs in the GSW data were identified. To this end, the maximal extent of each wetland in the GSW data was mapped. Then, wetlands that included at least an NWI GIW and also did not touch other non-GIW water bodies in the NWI data were considered as GIWs. The ten wetlandscapes in which GIWs were derived included the Basin wetlands (BAW), Cypress domes of the Southeast (CYDS), Coastal plain wetlands (COPW), Delmarva bays (DELB), Nebraska Sandhills (NSAN), Playa lakes of Texas and High Plains (PLTHP), Prairie potholes of the Upper Midwest (PPUM), Vernal pools of California (VPOC), Vernal pools of Maine (VPOM), and the Pocosins (POC). This study used the derived monthly inundation maps of GIWs to quantify their dynamic shape attributes.

3. Results

3.1. Covariation of Perimeter and Area Based on Static Wetland Data

The covariation is evaluated using the best non-linear locally estimated scatterplot smoothing (LOESS) fit line of perimeter–area covariation metrics vis à vis wetland size (Figure 2). While performing this fit, all selected GIWs across the ten wetlandscapes are used. The results show a clear dependence of the P:A ratio, SFD, and SI on wetland size (Figure 2). For larger GIWs, the P:A ratio and SFD generally decrease. This is unsurprising, as both metrics are non-dimensionless, with a length dimension in their denominator, and are expected to decrease with area. The decreasing trend in P:A highlights that with increase in area, there is a proportional but sublinear increase in perimeter. Notably, the two metrics reduce at a faster rate when the wetlands are smaller in size. As SFD is obtained as the ratio of log-transformed perimeter and area, the decrease in it is much smaller than that in the P:A ratio. In contrast to the variations in P:A ratio and SFD, the only dimensionless metric here, i.e., the SI, increases with size. This indicates that larger GIWs have more complex edges or non-rounded shape compared to the smaller GIWs per unit area. It should be acknowledged that for the smallest wetlands with area that encompass only a few pixels (with pixel size = 30 $\times$ 30 m), SI could be uncharacteristically small, as the perimeter representation will likely be much smoother than reality. The aforementioned trends in the considered metrics’ in relation to wetland size are similar regardless of the wetlandscape under consideration (Figure S2), even though the wetlandscapes belonged to diverse hydroclimatic settings [75].

3.2. Dynamic Covariation of Perimeter and Area

Variation of the P vs. A relationship with time is tracked by evaluating the three metrics as a function of wetness. Covariation of inundation area and perimeter, as incorporated in the three metrics, shows that shape metrics vary significantly with wetness (Figure 3). In other words, as the inundation extent of wetlands or its wetness changes in time, so does the shape metrics. For example, for wetlands with a maximum size within the range 5000 m² to 20,000 m², the P:A ratio reduces from 0.155 to 0.05 (a 68.1% reduction), while SFD reduces from 1.447 to 1.379 (a 4.7% reduction). The expressed trends are valid for all considered wetland size ranges. Interestingly, SI is observed to show a non-monotonic variation, with SI first increasing with wetness and then decreasing subsequently. This trend is more apparent for bigger wetlands (size > 5000 m²). The non-monotonic relation indicates that for larger wetlands, both during extremely dry and wet periods, the perimeter is less crenulated. Figure S3 confirms that these variations of metrics with wetness are generally valid in most wetlandscapes.

3.3. Discrepancies between Observed Dynamic P vs. A Variation and One Estimated from Static Data

To assess the covariation of area and perimeter with time, as the wetness varies, similar to that obtained using the static data, the discrepancies between estimated shape metrics and the actual value at different times are evaluated (Figure 4). Discrepancies in SFD are found to be modest, ranging from 0 to 10%. In contrast, the discrepancies in the P:A ratio and SI estimates are found to be as large 80 and 120%, respectively. Notably, the discrepancies are usually higher for larger wetlands, irrespective of the metric.

We further investigate if static relations that are developed specifically for individual wetlandscapes can be used to reduce discrepancies in the prediction of temporally evolving shape metrics. To this end, averaged prediction discrepancies are re-evaluated for each wetlandscape separately by using corresponding relations derived based on static data in those landscapes. Among the 87 cases (10 wetlandscapes * 3 wetland size groups − 3, see Figure 5), 69 cases (79.31%) show that the errors in estimated shape metrics reduce when static relations are developed using local data of the relevant wetlandscape. The corresponding fractions for the P:A ratio, SFD, and SI are 82.76%, 79.31%, and 75.86% (24/29, 23/29, and 22/29), respectively. Smaller errors in estimated shape metrics when developed using local data can be partially attributed to the similarity of the morphological and bathymetric characteristics of wetlands within a wetlandscape, even when they vary in size. In fact, if all wetlands in a wetlandscape exhibit an identical shape, the relationship between perimeter and area across different wetlands, as determined from static data, will resemble that derived from dynamic data. The largely uniform influence of hydrometeorological factors on inundation dynamics across wetlands within a wetlandscape is expected to also contribute to reduced errors in the estimation of shape metrics when obtained using static relations within a wetlandscape. On average, the impacts of reductions are the highest for large wetlands (Table S1). Among the ten wetlandscapes, discrepancy reduction is the highest (lowest) in Pocosins (Basin wetlands), with an average reduction of 36.28%, 26.07%, and 27.3.7% (20.12%, −28.57%, and −15.97%) for the P:A ratio, SFD, and SI, respectively.

4. Discussion

This study quantifies three well-known shape metrics for GIWs, viz. P:A ratio, SFD, and SI, and assesses their variations with wetland size. The results show that shape metrics have strong dependence on wetland size and inundation status. While this is unsurprising for shape metrics that are non-dimensionless, but even the SI metric, which is dimensionless, exhibits the influence of inundation status, i.e., what fraction of the wetland is inundated with respect to its maximal area, showing an increasing followed by a subsequently decreasing trend with inundation fraction. The results indicate that larger GIWs have more complex edges or non-rounded shapes compared to the smaller GIWs per unit area. The information regarding the variation of these shape metrics with wetland size can be used to design or restore wetlands with the goal to optimize the ecosystem service yield.

Next, the study evaluated the extent to which the shape metrics vs. wetland size relationship derived using static data, as has been achieved in the past, is valid for temporal datasets of wetland inundation extent. Evaluations for three well-known shape metrics, viz. P:A ratio, SFD, and SI, indicated that using relations derived from static data to estimate shape metrics at varying inundation extents introduced large discrepancies (up to 80%, 10%, and 120% for P:A ratio, SFD, and SI, respectively). Notably, these discrepancies are small for smaller wetlands. Overall, the results highlight the need to develop and incorporate temporally dynamic shape metrics or perimeter–area relationships to improve wetland function modeling. This is further corroborated by a recent study [77] that highlighted that the consideration of transient wetland inundation dynamics can increase nitrogen retention estimates by up to 130%. Given that the prediction errors reduced markedly when static relations were derived for each individual wetlandscape, the results underscore the need to develop perimeter–area relationships that are wetlandscape-specific. It is important to note that the shape factors considered here do not capture all aspects of wetland functions. Different wetland functions, such as recharge, nutrient cycling, habitat suitability, or carbon sequestration, may be more or less sensitive to wetland area and perimeter dynamics, and this relationship may change over time depending on changes in meteorological conditions and wetland substrate properties. Hence, the extent to which the temporally dynamic perimeter–area relationship may improve estimates of wetland functions are expected to be dependent on the ecosystem function under consideration, the representativeness of the shape metrics used to capture these functions, and site-specific properties. Despite this, the findings of this study clearly emphasize the necessity of relinquishing the limitations of static data analyses and embracing a more dynamic approach to shape metric assessment in models used for quantification of wetland ecosystem services.