Variability in marsh migration potential determined by topographic rather than anthropogenic constraints in the Chesapeake Bay region

Sea level rise (SLR) and saltwater intrusion are driving inland shifts in coastal ecosystems. Here, we make high‐resolution (1 m) predictions of land conversion under future SLR scenarios in 81 watersheds surrounding Chesapeake Bay, United States, a hotspot for accelerated SLR and saltwater intrusion. We find that 1050–3748 km2 of marsh could be created by 2100, largely at the expense of forested wetlands. Predicted marsh migration exceeds total current tidal marsh area and is ~ 4× greater than historical observations. Anthropogenic land use in marsh migration areas is concentrated within a few watersheds and minimally impacts calculated metrics of marsh resilience. Despite regional marsh area maintenance, local ecosystem service replacement within vulnerable watersheds remains uncertain. However, our work suggests that topography rather than land use drives spatial variability in wetland vulnerability regionally, and that rural land conversion is needed to compensate for extensive areal losses on heavily developed coasts globally.


Study area
Chesapeake Bay is microtidal, with the highest mean tidal range of 0.9 m at the mouth, and the lowest mean tidal range of 0.3 m at Annapolis, Maryland (Xiong and Berger 2010).Like much of the Mid-Atlantic, Chesapeake Bay is a global hotspot of relative sea level rise due to regional subsidence coupled with weakening of the Gulf Stream (Engelhart et al. 2009;Sallenger et al. 2012).Relative rates of sea level rise in the 1930s were 1-3 mm/yr1 .In 2011, SLR rates were 4-11 mm/yr and accelerating (Ezer and Corlett 2012;Ezer and Atkinson 2015).Low marsh vegetation typically includes Spartina alterniflora, and high marsh vegetation includes Spartina patens, Distichlis spicata, and Juncus romerianus (Perry et al. 2001).Coastal terrestrial forests which exist upland to adjacent marshes within this low elevation coastal area are typically comprised of Pinus taeda and Juniperus virginiana (Perry et al. 2001).Freshwater forested wetlands are primarily dominated by Nyssa biflora and Acer rubrum (Noe et al. 2021).Bands of Phragmites australis often occur at the marsh-forest boundary and are a sign of ecosystem change as the invasive reed colonizes areas of disturbance before marsh vegetation can establish (Smith 2013).
Part of the Mid-Atlantic, the Chesapeake Bay region is dominated by forest, cultivated/pastureland, and developed land uses (Epanchin-Niell et al. 2017).These land uses are distributed heterogeneously in the low-lying areas immediately surrounding Chesapeake Bay.Rural land uses, such as forests and palustrine forested wetlands, are the dominate land use classes.Southeastern Virginia, the Hampton Roads region, holds the majority of developed land use near the coast and has the highest concentration of people at risk from inundation with sea level rise in the region (Hauer et al. 2016).The Eastern Shore of Maryland and Virginia contains widespread agricultural land uses, including soybean fields and poultry production (USDA 2019).
Interestingly, the Chesapeake Bay region well represents the United States coastline which is largely rural, with pockets of major development.Along the North Atlantic Coastal Plain, 96% of land use within 1-m of sea level rise is rural (Epanchin-Niell et al. 2017), which includes our forest and forested wetland categories.Under higher sea level scenarios, 74% of U.S. watersheds vulnerable to sea level rise consist of inundated land that is predominantly rural, non-tidal wetlands, while only 16% of watersheds have developed land uses predominantly at risk for inundation (Holmquist et al. 2021).There are high concentrations of the U.S. population in major urban centers along the coast, typically not associated with the Chesapeake Bay.However, the Chesapeake Bay region (Virginia, Maryland, Washington D.C.) has the 5 th largest population at risk from 0.9 m of sea level rise, behind only Florida, Louisiana, California, and New Jersey (Hauer et al. 2016); this is based on 2010 Census Block Groups and is projected to increase with population growth.The largely rural nature of the Chesapeake Bay with pockets of concentrated development makes this region generally representative of the U.S. coast as a whole, and well positioned for a discussion of marsh migration potential beyond the Mid-Atlantic.The United States has a much lower coastal population compared to other countries located in regions of the world where marshes are prevalent (CIESIN -Columbia University 2017).In northwestern Europe and Asia, extensive levee systems protect these coastal populations, and limit the ability of marshes to migrate laterally (Kabat et al. 2005;Ma et al. 2014;Schuerch et al. 2018).This makes understanding where marshes have the potential to migrate landward in the United States increasingly important for estimates of future global marsh extent.

Data sources
The marsh-forest boundary points were obtained from a previously delineated marshforest boundary dataset comprised of over 200,000 points at a 30-m resolution along marshforest boundaries throughout Chesapeake Bay (Molino et al. 2020(Molino et al. , 2021)).This dataset used the National Wetlands Inventory (NWI) estuarine intertidal emergent wetlands category for salt marsh extent and forest data from the Chesapeake Conservancy to determine the boundary using geographic information system (GIS) analyses (Chesapeake Conservancy 2018a; b; U.S. Fish and Wildlife Service 2018).To determine the elevation along this boundary, we used the U.S. Geological Survey (USGS) Coastal National Elevation Database (CoNED) Topobathymetric Digital Elevation Model, which has a 1-m resolution (Danielson and Tyler 2016).Current salt marsh extent within Chesapeake Bay was determined from the NWI shapefiles (U.S. Fish and Wildlife Service 2018).
We used two different SLR projections for this study.The first were general SLR estimates of 0.5, 1.0, 1.5, 2.0 and 2.5 m.The second sea level rise projections were based on the Sweet et al., 2017 projections provided at each National Oceanographic and Atmospheric Administration (NOAA) tide gauge (Sweet et al. 2017).We used the projections for tide gauges which fell within our study region and calculated a Delauney Triangulation between the points.We then took a spatially weighted average to determine a single SLR value for the entire Chesapeake Bay.This method was repeated at every 10-year time step for each NOAA sea level rise scenario we considered: global Low (0.3 m), Intermediate (1.0 m), and High (2.0 m) scenario estimates for 2100 (corresponding to 0.45, 1.22, and 2.66 m of SLR, respectively) (Sweet et al. 2017).
Land use of area converted under the general SLR scenarios was assessed using the Chesapeake Conservancy High Resolution (1 m) Land Use and Land Cover datasets (Chesapeake Conservancy 2018a; b).We merged the land use types into five categories: forest (value = 8), turf grass (values = 9, 11, 12, 13, 15), agriculture (values = 16, 17), impervious (values = 1, 2, 3), and other (values = 5, 6, 7, 10, 14) (Supplemental Table 2).A sixth category was also included, forested wetlands, which was calculated by the overlap in land use wetland classes (values = 5, 6, 7) and the land cover tree canopy class (value = 3).This category was included separately from the wetland and forest categories because the physiographic position of forested wetlands results in a higher sensitivity of this ecosystem to flooding and salinity stresses (Doyle et al. 2007).The "other" category includes mixed open and mixed impervious as well as marsh that is located at elevations above the median threshold value for each watershed.For a given threshold elevation, we assume all land below this value is already marsh, and all land above this value is not marsh.However, because the threshold is based on the median elevation of all the marsh-forest boundary points, there is land below the threshold that is not marsh (i.e.forests and forested wetlands), and there is land above the median that is already marsh.We quantified the forest and forested wetland areas below threshold elevations for each watershed and found that it was greater than the area of marsh currently above the threshold elevation.Therefore, we include marsh above the threshold elevation in our estimates of potential upland area converted as a conservative estimate for the forest and forested wetland area below the threshold elevation which we assume will convert to marsh under each future SLR scenario.

Analytical approach
We quantified the area between the current marsh-forest boundary and multiple increments of sea level rise in ArcMap 10.7 to obtain an estimate of potential marsh migration area for the entire Chesapeake Bay coastal plain (Enwright et al. 2016).We took the points from the previously delineated marsh-forest boundary dataset for the Chesapeake Bay (Molino et al. 2020(Molino et al. , 2021) ) and extracted the elevation of each point from the USGS CoNED topobathy.We then eliminated points which were located along boundaries with >4% slope to minimize error associated with offset in lateral position.Additionally, points with an elevation value of less than 0.1 m NAVD88 (North American Vertical Datum of 1988) were removed as these were often located in marsh channels or other bodies of water, and points with an elevation greater than 2.0 m NAVD88 were removed as these were mostly located in residential areas adjacent to marshes.Together, these steps helped minimize any potential error associated with gross misclassification of the marsh-forest boundary.The inaccuracy of the point locations was due to misrepresentation of the marsh-forest boundary delineated prior to this study.
We then calculated the median elevation of the marsh-forest boundary, hereafter referred to as the threshold elevation, for USGS HUC10 (Hydrologic Unit) watersheds (USGS 2020).This approach allows us to account for spatial variability in the processes that likely influence the threshold elevation (e.g.tidal range, salinity) and quantify marsh migration at the watershed scale.The points along the marsh-forest boundary were contained within 86 HUC10 watersheds along the Chesapeake Bay (Supporting Information Fig. S1) (Molino et al. 2022).Four watersheds were removed because they contained areas of abrupt elevation shifts (>1 m change) corresponding to where individual elevation datasets had been combined to create the CoNED topobathy.For 11 watersheds that contained less than 100 points each, the threshold elevation was determined by additionally including all points from a neighboring watershed with similar topobathymetric features.One watershed that contained less than ten points was removed as it had no neighboring watersheds from which to combine points.Thus, our analysis includes 81 of 86 HUC10 watersheds surrounding Chesapeake Bay where elevation of the marsh-forest boundary can be most precisely defined (Figure 1A).
To determine areas of future marsh migration, we added increments of sea level rise to the threshold elevation of each watershed and quantified the area between the original marshforest boundary and new marsh-upland boundary.To do this, we split the CoNED topobathy raster into 81 watersheds and set the unique threshold elevation for each watershed as the baseline.Below this threshold, all land was assumed to be either water or marsh.We then reclassified the raster values by adding the increment of sea level rise we were interested in to the baseline and quantifying the area in between the threshold elevation and the threshold elevation plus the SLR scenario.For example, a watershed with a threshold elevation of 0.32 m assumes that the land which would convert to marsh with 0.5 m of SLR would be land from 0.32 m to 0.82 m (0.32m + 0.50m) in elevation.We used the reclassify tool in ArcMap which quantified the raster cells whose elevation fell between these elevations for each SLR scenario.This was completed for each watershed for the general SLR scenarios as well as the Chesapeake Bay specific predictions derived from the NOAA data.Marsh migration areas for each watershed were combined for a bay-wide projection of potential marsh migration area (Molino et al. 2022).Like previous work, we only consider land available for marsh migration, and do not account for accretion or edge erosion (Enwright et al. 2016;Borchert et al. 2018;Holmquist et al. 2021).However, we are also not accounting for processes such as peat collapse which could accelerate marsh migration (Miller et al. 2021).Following Gesch (2012), we calculated uncertainty envelopes for area converted under the localized NOAA predictions by adding and subtracting the root mean square error (RMSE) of the CoNED dataset (20 cm) from the new marsh-upland elevation within each watershed (Danielson and Tyler, 2016).
Land use inundated under each SLR scenario was quantified by individual watershed.To do this, we extracted the land use within the area that was projected to convert to marsh under the different sea level rise scenarios determined during the reclassification process described above.We combined these numbers for a bay-wide assessment of land use within the potential marsh migration area.This allowed us to examine how local variations in land use influenced marsh migration area regionally.

Limitations
We assume that the threshold elevations determined from the marsh-forest boundary are representative of marsh-upland boundaries in general given that forested uplands make up more than half of the upland land use and that other upland land use types, such as agricultural land, are often separated from wetlands by a thin band of forest.Preliminary results from analyses of marsh-agriculture, forested wetlands, impervious, and turf boundaries suggest threshold elevations for other land uses are similar to that of the marsh-forest boundary (Supplemental Table 1).The marsh-agricultural boundary was significantly higher in elevation from the marshforest boundary, which we hypothesize is for two reasons.First, some agricultural areas are located on high elevation river bluffs.The very steep, but narrow change in elevation means that the boundary between these two land uses is at a very high elevation.The second and most common reason is that agricultural areas often have small earthen levees (~1 m in height) to protect them from tidal inundation.Without these anthropogenic features, the natural elevation between marsh-agriculture would likely be similar to the other boundary elevations.However, more work is needed to refine this method for other land uses.
It should be noted that we do not directly assess the effect of connectivity on our calculated marsh migration area.When using high resolution digital elevation models for inundation estimates, it is possible that some isolated low-lying areas will be counted in the area of potential marsh migration even though they are not hydrologically connected to other marsh or water.However, this issue is minimized because our approach to calculating potential marsh migration includes only land that is above a watershed's threshold elevation.For example, there is a quarry in a northern watershed of Chesapeake Bay that is excluded from the analysis because it has elevations below the watershed's threshold elevation.As described earlier, our approach explicitly calculates marsh migration by counting the number of raster cells with elevations between the threshold elevation and the threshold elevation plus an increment of sea level rise (i.e., migration area includes locations with elevations between 0.61 and 1.61 m for a watershed with a threshold elevation of 0.61 m and a 1 m increment of sea level rise).In that example, elevations between 0.61-1.61m NAVD would be counted as potential migration area, regardless of connectivity.However, previous inundation mapping with high resolution (6 m) DEMs in coastal North Carolina suggests that connectivity effects are minimal for SLR scenarios greater than about 0.4 m (Poulter and Halpin 2008).For a 1 m SLR scenario, there was less than 3% difference in cumulative area inundated between approaches that did not consider connectivity and those which considered 4-sided or 8-sided connectivity (Poulter and Halpin 2008).Six of our eight sea level rise scenarios use projections greater than 1 m which suggests that the overestimation of area from low-lying area upslope is minimal.Nevertheless, future steps quantifying this potential area would improve local estimates of marsh migration on a watershed scale.
Additionally, our quantification of land use is limited by the accuracy of the land use and land cover datasets.In particular, forest and forested wetlands can be confused by remote sensing techniques (McCarthy et al. 2018).As a result, we combine all forested areas when delineating the marsh-forest boundary as those measurements rely on precise delineation between ecosystem types.Future work could benefit from improvements in identifying between forest and forested wetlands.Supp.Table S4: Area of uplands converted by 2100 under NOAA Low, Intermediate, and High SLR predictions and five generalized sea level rise (SLR) scenarios.Uncertainty envelopes were calculated for the NOAA scenarios based on the CoNED topobathy root mean square error (20 cm).Blank cells indicate SLR scenarios for which uncertainty estimates were not conducted.
Supp.Table S5: Area of land use types (km 2 ) predicted to convert to marsh under each sea level rise (SLR) scenario.
Supp.Table S6: Percentage of agricultural, impervious, and developed (combined agriculture and impervious) land use within predicted marsh migration area for 1.0 m of sea level rise (SLR) for each HUC watershed.Adjusted 1:1 has the developed land use removed from the marsh migration area when calculating the resilience index.Please note that watersheds 19 and 23 have the same name -this is not a mistake, those are the names from the original source data.Instead we differentiate them by their unique HUC ID numbers.For the location of each watershed, see Supp.

Figure
S1.Supp.TableS7: Percentage of forest and forested wetland land use within predicted marsh migration area for 1.0 m of sea level rise (SLR) for each HUC watershed.Adjusted 1:1 has the forest and forested wetland land use removed from the marsh migration area when calculating the resilience index.Please note that watersheds 19 and 23 have the same name -this is not a mistake, those are the names from the original source data.Instead we differentiate them by their unique HUC ID numbers.For the location of each watershed, see Supp.FigureS1.