The San Francisco Bay Estuary (SFBE) in California is the second largest estuary in the United States and comprises nine counties with a population of over 8 million people (United States Census Bureau, 2021). The estuary has a Mediterranean climate with mild rainy winters and cool foggy summers. The weather is largely influenced by the cold water of the California Current flowing toward the equator year-round and the Davidson Current bringing warmer water toward the pole during the winter months. The SFBE is made up of multiple regions with different environments – South San Francisco Bay, Central San Francisco Bay, San Pablo Bay, Suisun Bay, and the Delta (Figure 1). In this inverted estuary, the Sacramento and San Joaquin Rivers flow through the Delta into brackish Suisun Bay and through the narrow Carquinez Strait into the greater San Francisco Bay embayment (Cloern and Jassby, 2012). Tidal marine waters are exchanged through the Golden Gate to the Pacific Ocean. The tidal regime is mixed semidiurnal with an average diurnal range of 2.6 m in South San Francisco (SF) Bay, 1.8 m in San Pablo Bay, and 1.2 m in the Delta (NOAA tide stations 9414509, 9415252, 9415316; https://tidesandcurrents.noaa.gov). Despite its highly urbanized landscape (Nichols et al., 1986), the SFBE is home to an array of wildlife that use tidal wetlands, including state and federally listed species such as the Chinook salmon (Oncorhynchus tshawytscha), California Ridgway’s Rail (Rallus longirostris obsoletus), salt marsh harvest mouse (Reithrodontomys raviventris), Delta smelt (Hypomesus transpacificus), and longfin smelt (Spirinchus thaleichthys) (California Endangered Species Act 1970; Federal Endangered Species Act 1973; California Natural Diversity Database (CNDDB), 2023). It is important habitat for migratory and nesting birds, supporting more than 50% of the Pacific flyway waterfowl and millions of shorebirds (Takekawa et al., 2001). Tidal marshes in the more saline parts of the estuary generally lie above mean tide level (MTL) and are characterized by relatively flat elevation platforms occurring around mean high water (MHW) (Takekawa et al., 2013a). Extensive modifications of the estuary have resulted in 80% loss of historic tidal wetlands, but restoration efforts are ongoing (Marcus, 2000).

For this effort, fifty-one tidal marsh sites across five regions (sub-embayments) were surveyed in the Delta, Suisun Bay, San Pablo Bay, Central SF Bay, and South SF Bay (Figure 1). These sites varied in geographic size, sample size, geomorphic setting, and tide range (Table S1). Some sites were close in proximity, such as many of the South Bay sites. However, sites were separated by large channels or land features that delineated their physical boundaries.

Vegetation surveys spanned ten years, from July 2008 to October 2018 (Table S1). A total of 5,026 plots were surveyed. Most plots were positioned on transects along an elevation gradient, while non-transect plot placement was used at one site (Arrowhead Marsh) due to its small size and the desire to avoid endangered wildlife species. Most transects were evenly distributed across each site, where possible, to capture spatial variability along elevation and distance gradients. At each plot, percent cover of all plant species, bare ground, and litter was visually assessed within a 0.25 m² quadrat. Total plant cover in a plot could exceed 100 percent due to vegetation layering. Bare ground and litter cover was estimated as total area visible through the vegetation from above the plot. Vascular plant nomenclature followed Baldwin et al. (2012). Geographic position (in UTM) and elevation (North American Vertical Datum of 1988, NAVD88) were recorded at each plot. Sampling season was not expected to significantly influence species identification and cover estimates. The Mediterranean climate of the SFBE with short mild winters creates fairly unseasonal wetlands. Most of the wetland species are perennial and identifiable throughout the year.

Environmental gradients were selected based on previous studies that identified soil salinity and tidal inundation as the most influential drivers of tidal marsh vegetation structure (Engels and Jensen, 2009; Watson and Byrne, 2009; Janousek and Folger, 2014). Channel salinity was calculated at the site or multi-site scale (in PSU, Table S1) with data from locally deployed water sensors or from San Francisco Bay National Estuarine Research Reserve (NERR) sites (Takekawa et al., 2013a, Thorne et al., 2019; Table S2). Salinity dataloggers (Solinst and Odyssey) were programmed to record specific conductivity (in µS/cm) at 12-minute (Solinst) or 30-minute (Odyssey) intervals and were deployed for a minimum of one year. Annual mean salinity was calculated across the deployment period and converted to PSU. Salinity values were assigned to each marsh based on nearest available measurements, from either a NERR site or logger deployment. All salinity values represent marsh channel or creek conditions, and not soil porewater salinity.

A survey-grade Real Time Kinematic (RTK) global positioning systems (GPS) rover was used to measure location and elevation (± 1 cm horizontal, ± 2 cm vertical manufacturer-stated accuracy; Leica Geosystems Inc., Norcross, Georgia). Rover positions were received in real time from the Leica Smartnet system using a CDMA modem (http://www.leica-geosystems.us/en/index.htm). The WGS 84 ellipsoid model was used for horizontal positioning and NAVD88 for vertical positioning. Rover accuracy and precision were evaluated by measuring positions at local National Geodetic Survey benchmarks; all errors were within the stated rover error. Elevation data was converted to z* (z* = [NAVD88 − MTL]/[MHHW – MTL]), a unitless measure of elevation relative to the local tidal frame which accounts for variation in tidal range and allows for direct comparison across sites (Swanson et al., 2014).

Local tidal datums were calculated from multiple sources including NOAA tide stations (https://tidesandcurrents.noaa.gov/), deployed water level loggers, and NOAA’s VDATUM model (Parker et al., 2003; Table S2). Water loggers (Solinst) were deployed at select sites and programmed to record water level information at 6-minute intervals. Data were corrected for atmospheric pressure with barometric data from barometric dataloggers (Solinst; Takekawa et al., 2013a; Thorne et al., 2019). These water level values were used to calculate mean higher high water (MHHW) while mean tidal level (MTL) was estimated from NOAA’s VDATUM model (Parker et al., 2003). At some sites, we used tidal datums from nearby NOAA tide gages where local loggers were absent (Table S2).

Distance to channel was calculated as the distance (in meters) between each plot and the digitalized boundaries of all nearby channels and bays. Channels were digitized based on 2020 NAIP imagery supplemented with LiDAR when necessary. The centerline was digitized on small channels (~1 – 2.5 m wide), while both edges were digitized on large channels (> 2.5 m), bays, and rivers.

To understand the relationship between plant communities and environmental gradients at different scales, generalized linear latent variable models (GLLVMs; Skrondal and Rabe-Hesketh, 2004) were analyzed in R version 4.2.2 (R Core Team, 2022; package ‘GLLVM’; Niku et al., 2019; Niku et al., 2021) at the site and estuary levels. This type of model incorporates latent variables that quantify species response correlations and environment-trait interactions, known as fourth-corner terms (Skrondal and Rabe-Hesketh, 2004; Niku et al., 2021). At the site scale, GLLVMs were built with elevation (tidal position) and channel distance as environmental covariates. At the estuary scale, one model was built using all plot data. Site sampling intensity was accounted for by including row effects by site. Salinity was added as an environmental covariate to the estuary model as salinity data were collected at the site/multisite levels. All models were fit using negative binomial distributions. They identified species with strong associations, either negative or positive, between environmental covariates and species presence and abundance. Strong associations were defined by models with confidence intervals that did not overlap with zero. Associations were used to identify how elevation (tidal position), channel distance, and salinity influenced plant presence and abundance at the site and estuary scales. Species that changed associations between the site and estuary models were identified to examine scale-dependency.

Species associations from site and estuary abundance models were combined with percent cover values from field surveys to calculate site-specific cover of species with strong associations to each environmental variable. All species with strong associations (either negative or positive for each environmental covariate) were identified by site, and their mean percent cover values were added together to compute a total percent cover of strongly associated species. Sites with high total cover values indicate marsh communities most influenced by elevation, channel distance, or salinity.

To investigate the impacts of environment conditions on species interactions, species correlation matrices were compared between models with and without environmental predictors (elevation, channel distance, and salinity). Species positively correlated with each other could indicate facilitation while competition may be identified by negative correlations. Latent variables used in GLLVMs include correlation values across all model response variables used to create correlation matrices. Species relationship differences across environmental gradients were also identified by comparing correlation matrices from data subset along gradient thresholds.

To examine why species differ in their response to environmental drivers, species traits were incorporated into a “fourth corner” estuary model. Species traits included: California native/non-native status, wetland indicator status (upland – almost always occurs in upland habitat; facultative upland – usually occurs in upland habitat; facultative – equally likely in upland and wetland habitats; facultative wetland – usually occurs in wetlands; and obligate – almost always occurs in wetlands), and plant form (tree, shrub, forb/herb, and graminoid) (USDA, 2023; calflora.org; Baldwin et al., 2012). Strong environmental–trait interactions, or fourth corner terms, were identified from model coefficient plots.

Differences in relative species abundance were examined by creating rank abundance curves (RACs) for each site and region as well as environmental covariates. Species richness and species evenness were evaluated by comparing RACs within regions and across the estuary. RAC differences were compared by site and region using R (package ‘CODYN’; Hallett et al., 2020).

Eighty-one species were identified across all sites surveyed (Table S3). The number of species per site ranged from 1 to 31 with an average of 8.0 ± 6.5 species per site. Most plots were dominated by one or two plant species.

At the site scale, elevation/tidal position (z*) was strongly associated with plant cover for more than half of the species (58%; Figure 2A). Most of these species (70%) had a positive association with elevation – they were more abundant at higher elevations. Some species that were more abundant at lower elevations included Schoenoplectus americanus and Spartina foliosa. Distance to channel had a weaker relationship with species abundance with an average association of 52%. Most of these species (73%) had a negative association with distance to channel indicating they were more abundant closer to a tidal channel or bay edge, such as Spartina foliosa.

At the estuary scale, environmental covariates (elevation/tidal position, distance to channel, and salinity) explained about 49% of the total community variation (Figure 2A). After incorporating region as a variable, the covariates explained about 65% of the variation. Of the three environmental covariates we examined, elevation and distance to channel had the strongest association with species composition. Elevation was associated with species abundance for 64% of the species, distance to channel for 57%, and salinity for 51% (Figure 2A). Coefficient estimates varied between site and estuary scale models as seen in 9 abundant species of interest in the estuary (Figure 2B).

Assigning a single environmental association for each species across the estuary, as opposed to each site, resulted in associations that were either confirmed (no association to positive or negative), changed (negative to positive or vice versa), or lost (positive or negative to no association) (Figure 3). With an increase in sample size, the estuary model confirmed existing positive associations to elevation (tidal position) for species such as Atriplex prostrata and Frankenia salina (Figure 3A). At the estuary level, Typha species had a negative association to elevation (tidal position) but at the site level these species had no strong association or even a positive association to tidal position in the Delta (Figure 3A). Salicornia pacifica had a negative association and no association to tidal position at multiple marsh sites (7 and 25 respectively of 49 sites), yet the estuary model identified a strong positive association (Figure 3A).

Schoenoplectus americanus had a strong association to tidal position at most of its marsh sites but did not have a strong association at the estuary scale (Figure 3A).

Similar differences in species association to channel distance occurred between site and estuary models. The increased sample size confirmed negative associations for Grindelia stricta and Limonium californicum, but more species changed or lost associations between scales (Figure 3B). Schoenoplectus americanus had a strong positive association at 5 sites and no association at 10 sites (of 20 sites), yet a negative association at the estuary scale (Figure 3B). Jaumea carnosa had a positive association with channel distance at 7 of 41 sites, including all sites in the Delta and Suisun Bay where it was located, yet a negative association across the estuary. Atriplex prostata had a strong association to channel distance at 9 of 21 marsh sites but did not have a strong association at the estuary scale (Figure 3B).

The total percent cover of species with strong association to these variables differed greatly by site (Figures 4A, C, Table S4). Sites with the highest percent cover of species associated with elevation were in Suisun Bay and the Delta, while the sites with the lowest percent cover were in South SF Bay (Figure 4A). Total percent cover of species associated with channel distance varied across the regions with the highest cover in Suisun Bay (Figure 4C). Given that species associations varied considerably between site models, associations at the site scale did not always correspond to associations at the estuary scale. Spatially, there was a narrower range in percent cover values across sites based on the estuary-wide model (Figures 4B, D, E, Table S4). Most sites had moderate to high cover of species strongly associated with elevation and channel distance (Figures 4B, D). The greatest cover of species associated with salinity occurred at more inland sites associated with freshwater sources in Suisun Bay, the Delta, and the southeastern area of South SF Bay (Figure 4E). More species had strong associations with elevation and channel distance when modeled at the estuary scale compared to the site scale (Figure S5). Species associations were less frequent for channel salinity (Figure S5).

Without accounting for environmental conditions, multiple positive species correlations were identified (Figure 5A). When environmental variables were incorporated into the model, almost no species correlations were observed (Figure 5A). Species were most likely found together because they had similar habitat preferences. The environmental gradients explained nearly all the species interactions when modeled across the entire estuary. However, species co-occurrence patterns did vary along the salinity gradient, where both positive and negative correlations were more common adjacent to less saline channels (Figure 5B). Species co-occurrence patterns appeared unaffected by position within the tidal frame and proximity to channel edge (Figure 5B).

Multiple environment–trait interactions were observed between regions and along environmental gradients (Figure 6). As expected, trees were more abundant higher in the tidal frame and farther from channel edges whereas shrubs, forbs/herbs, and graminoids were found closer to channels. Only forbs/herbs showed interactions along the salinity gradient as they were less abundant near saline channels. Species wetland status had no impact on abundance across the elevation gradient. All species had varying degrees of negative association to channel distance regardless of wetland status. Only upland species showed an interaction with salinity as they were more abundant near saline channels.

The influence of all species traits (form, wetland status, and native status) on abundance varied between regions (Figure 6). Trees were more abundant in San Pablo Bay and less abundant in Suisun Bay whereas shrubs had the opposite relationship. Forbs/herbs were more abundant in San Pablo Bay and graminids more abundant in South SF Bay. Obligate wetland species were less abundant in Central and South SF Bay. Facultative species were more abundant solely in South SF Bay. Facultative upland species were less abundant at all regions except for the Delta. Native plant species were less abundant in Suisun Bay. Non-native species were most abundant in the Delta and South SF Bay and less abundant in San Pablo Bay.

RACs illustrated different levels of species dominance by site and region (Figure 7). Lower species evenness was observed at sites in San Pablo Bay, Central SF Bay, South SF Bay, and Delta. Most of these sites had similar species diversity, with Salicornia pacifica as the most abundant species (Figures 7B, D–F). The more inland sites at Suisun Bay exhibited higher evenness. Most sites in this fresher part of the estuary (Delta and Suisun Bay) did not have an individual species that accounted for over 50% of the proportional abundance (Figures 7A, B).

RACs by environmental gradients demonstrated a clear influence of flooding and salinity on species rank and abundance (Figure 8). Lower species evenness and similar dominant species (Salicornia pacifica and Schoenoplectus acutus) were measured at plots lower in the tidal frame (z* < 1.5; Figure 8A). Plots higher in the tidal frame (z* > 1.5) exhibited higher species evenness with variable dominant species (Festuca microstchys and Baccharis pilularis). There was no measurable difference in species rank or abundance by distance to channel (Figure 8B), Salicornia pacifica was always most abundant. The communities located near oligohaline channels with assumed lower soil salinity had more species evenness, while communities adjacent to more saline channels were dominated by Salicornia pacifica (Figure 8C).

Pairwise comparisons of site RACs revealed that sites differed most by species richness, calculated as the difference in richness between samples divided by the number of unique species in both samples (40%; Figure S6A). Species rank was the next most important measure of site difference (25%) followed by species number (22%) and evenness (15%; Figure S6A). Species number is calculated as the difference in species number between samples divided by the total number of species in both samples. While all sites were treated independently in the RAC comparison, sites within the same region would be expected to have similar RACs. The unequal number of sites by region could possibly affect which metric detected the most difference in RACs. However, when comparing differences by regional RACs, the order and magnitude of importance metrics were comparable (species richness – 49%, rank – 27%, species number – 26%, and evenness – 3%) (Figure S6B). Species richness and rank explained most of the variability in the plant communities across the estuary at the site and regional scale.

Within tidal marshes, plant functional traits have the potential to influence both marsh processes (e.g., accretion, channel formation) and species distributions (Schwarz et al., 2018). In our analysis, there was some variation between the environmental gradients and species trait relationships (Figure 6). Trait relationships varied between regions, possibly signaling additional environmental conditions that were not represented in the models or region-specific species responses. Tidal wetlands differ in channel structure, hydrology, elevation, and salinity which can result in variable relationships with common plant species (Watson and Byrne, 2009; Janousek and Folger, 2014). Understanding linkages between biogeomorphic features in estuaries and plant species composition is key to understanding vegetation community organization at small and large spatial scales (Van de Koppel et al., 2012). While all study sites were tidal with similar biogeomorphic features, the patterns of old dikes, berms, and levees throughout each site could have impacted plant communities at a much finer scale. Individually constructed earthen levees can impact ecosystem processes at coastal marshes long after abandonment (Hall et al., 2022). These “hidden levees” can cause a disconnect between channel distance, elevation, and inundation regime resulting in less flood tolerant vegetation and shallower organic rich soils landward of the abandoned levees (Hall et al., 2022). More research is needed on the impacts of past biogeomorphic features at these sites.

Our analysis showed that distance to channel and elevation relative to tidal flooding were the primary drivers for vegetation community structure and plant variability, but the magnitude of importance of these drivers varied by marsh site and scale. However, other factors related to flooding could be useful in determining stress gradient impacts, such as inundation frequency and duration. Future studies would benefit from more in-depth interaction measurements and species-specific responses to environmental stress gradients. Different halophyte plant species have varying tolerances to flooding; for example Spartina foliosa (California cordgrass) can tolerate extended flooding and Salicornia pacifica (pickleweed) is tolerant of higher salinities (Woo and Takekawa, 2012; Takekawa et al., 2013b; Gallego-Tévar et al., 2020). The spatial variability and distribution of plants by traits can be used to infer multiple habitat characteristics, but recent research suggests that there are limits to the extent to which traits can be related to ecosystem properties (Van der Plas et al., 2020).

This analysis identified elevation (tidal frame position), channel distance, and channel salinity as variables that strongly influenced abundance for more than half of the species at both the marsh site and estuary scales. Previous research found similarly strong responses of plant diversity to channel distance (Sanderson et al., 2001). The observed difference in magnitude of effect is supported by previous work which found plant diversity was dependent on scale (Jones et al., 2021; Korell et al., 2021), elevation (Janousek et al., 2019), and salinity (Bonin and Zedler, 2008). At the site scale, the highest total cover of species associated with elevation were in the Delta and Suisun Bay, possibly indicating more zonation in areas with lower tidal ranges and an increase in plant cover related to brackish and freshwater conditions; these marshes are dominated by grass and graminoid families such as Poacae (grasses), Cyperaceae (sedges), and Juncaceae (rushes). Sites with high cover associated with channel distance were spread across the estuary at both the marsh site and estuary scales. Many marshes have higher cover near channels due to the presence of Grindelia stricta (gumplant) and Salicornia pacifica (pickleweed). There are other additional factors with site-specific gradients, such as soil salinity, that may contribute to the importance of channel distance in species abundance. Sites with freshwater influence (the Delta, Suisun Bay, inland South SF Bay) had more species with abundance driven by channel salinity. Future studies should include soil salinity, rather than channel salinity, to capture site variability and improve site-level model coefficient strength.

As sea levels rise, marsh plant communities will shift from primarily high marsh species to low marsh species (Donnelly and Bertness, 2001; Wasson et al., 2013), significantly changing habitats. These types of plant species shifts have been projected under accelerating SLR scenarios, demonstrating drastic changes over the coming century with loss of relative elevation (Schile et al., 2014; Thorne et al., 2018; Buffington et al., 2021). At the site level, Miner Slough in the Delta and Coyote Creek in South SF Bay had the highest plant cover negatively associated with salinity. In other words, plant abundance was highest at low salinities throughout these sites, making them particularity vulnerable to saltwater intrusion from SLR. Miner Slough is a freshwater tidal marsh site located off the Sacramento River and the only freshwater Delta site in our study; therefore, a larger sample size of similar marshes is needed to determine if this pattern holds. Coyote Creek marsh is located on a large channel with freshwater input from an upland sewage treatment plant, and therefore, is a unique case study in SFBE.

Species were found to respond differently to environmental factors depending on the scale. For example, Schoenoplectus americanus, a bulrush native to CA, exhibited a strong positive association with elevation in the estuary scale model and in multiple, but not all, marsh site models. Strong negative associations to elevation were documented at the two Suisun Bay sites as well as half of the San Pablo Bay sites. Similar site-specific variations in plant vertical zonation have been previously observed along the Pacific coast (Janousek et al., 2019). As an obligate wetland species, S. americanus is expected to be more abundant closer to channels (as was seen in the estuary model), but some site models documented the opposite relationship (namely both Suisun Bay sites and a few sites in San Pablo Bay and South SF Bay). Species-environment associations that changed between site and estuary models indicate scale dependency for select species and environmental gradients (Figure 3). A combination of multiple environmental factors, rather than a single gradient, is likely responsible for marsh plant zonation (Eleuterius and Eleuterius, 1979; Silvestri et al., 2005). Other factors with the potential to influence plant diversity and composition at different scales include competition (Costa et al., 2003), microtopography (Courtwright and Findlay, 2011), soil salinity (He et al., 2011), herbivory (Rand, 2002; Elschot et al., 2017), geomorphic setting (Shipley, 2010), and local climate (De Leeuw et al., 1991). There could be fewer species capable of withstanding saline conditions, impacting diversity and rank abundance. Saltwater intrusion in some marshes has also been shown to strongly affect species composition rather than factors such as biomass production (Li and Pennings, 2019). More investigation is needed to understand ecosystem wide drivers of plant communities, which could be done with remote imagery (e.g., Li et al., 2005; Rosso et al., 2005; Adam et al., 2010) and long-term monitoring of abiotic drivers (Trowbridge et al., 2016).

The difference in plant associations from the site to estuary scales and overall trait interactions could be influenced by unequal samples within regions. South SF Bay had the most sites (35 sites; 1,567 plots) and San Pablo Bay had the most plots (11 sites; 2,864 plots) compared to Suisun Bay (2 sites; 176 plots), the Delta (2 sites; 88 plots), and Central SF Bay (1 site; 361 plots). The dataset also contained unequal covariate distributions with most plots existing close to channels, near MHHW (z* = 1), and adjacent to mesohaline channels (Figures S1–S4). Given the deficient quantity of freshwater tidal marsh sites in the estuary due to development, there was not enough power to confidently detect the loss of rare plants, something common to freshwater systems since these communities experience drastic changes in response to shifting environmental gradients (Neubauer, 2013; Li et al., 2022).

Species associations with environmental gradients were sometimes stronger (larger coefficients and greater species cover; Figures 2, 4) at the site scale compared to the estuary scale, but when averaged across all sites, the estuary scale model had more species associations (Figure S5). Even though species pools were fairly similar across all sites (10 species made up 90% of occurrences), the relative abundance was determined by prevailing environmental conditions specific to the site. Local changes in environmental drivers influence plant communities within a small regional inference space, as seen in the variable RACs (Figures 7, 8). Understanding how local conditions influence plant variability is crucial in land management and conservation. At landscape scales, changes in vegetation patterns have been documented by large disturbance events (e.g., fire; Taylor et al., 2021) and human impacts (e.g., human pressure; Malavasi et al., 2016). Ecological studies at small scales would require different interpretations than a large-scale study.

There has been a growing interest in how best to track the effects of our changing climate. Increased flooding due to SLR can lead to vegetation loss and conversion of marsh to unvegetated mud flat (Orson et al., 1985; Morris et al., 2002), as well as an increase in erosion, creation of new channels, and the expansion of existing channels (Moffett and Gorelick, 2016). Increasing inundation can affect salinity concentrations leading to changes in plant communities, biomass, and productivity (Engels and Jensen, 2009; Ryan and Boyer, 2012; Janousek and Mayo, 2013; Snedden et al., 2015; Janousek et al., 2016). Previous modeling in California marshes showed drastic changes in marsh vegetation under SLR scenarios (Thorne et al., 2018); however, these models are based on generalizations of plant tolerances to flooding and focus on main functional plant groups, therefore missing plant community or diversity concerns that may manifest prior to state change from vegetated to unvegetated systems.

SLR impacts to marsh plant diversity are dependent on species associations and position along environmental gradients. Species more abundant at higher marsh elevations will be negatively impacted by rising sea levels due to loss of habitat, vegetation vigor, and soil condition from increased flooding time (Reed and Cahoon, 1992). If a physical barrier prevents upland migration, the marsh will experience a phenomenon known as “coastal squeeze” (Pontee, 2013). Species with abundance influenced by channel proximity may be negatively impacted as SLR inundates species near water forcing migration from channel widening (Hartig et al., 2002), which narrows existing habitat for species far from water. SLR will also expose more marsh plants to saline water, with saltwater intrusion predicted to be more common in the future (Cloern et al., 2011).

Future SLR is anticipated to lead to flooding, erosion, and saltwater intrusion (Cloern et al., 2011; Thorne et al., 2018). Sweet et al. (2022) project an increase in sea level between 0.3 m and 2.5 m by 2100 for the SFBE, depending on realized levels of global greenhouse gas emissions. The results from this study can be used when considering an unpredictable future related to SLR. Here we identified species community and abundance relationships across the SFBE and investigated how they relate to environmental variables. With rising sea levels, plant communities will shift lower in the tidal frame and experience more frequent saltwater intrusion. This shift is likely to result in less species evenness as select species become dominant in these lower and saltier environments (Figure 8). This type of information can be used by ecologists and wetland managers to project how vegetation communities will change in response to SLR and develop monitoring strategies that are scale dependent. Slight changes in vegetation composition or associations could be early indicators of change within the abiotic conditions of the estuary (Kearney et al., 2011; Schepers et al., 2017; Kutcher et al., 2022). This approach has been successfully deployed to monitor marsh restoration and holds the same premises (Handa et al., 2002).

Large scale analyses for abiotic drivers may be biased by other marsh associations and mask individual drivers at the local scale. At the local scale, site-specific drivers may dominate the influence on changes for species that were associated with high elevation, low salinity, and interior locations (Figure 2). At the regional or estuary scale, species may re-order to look more like the saline San Francisco Bay, dominated by a few species, and less like Suisun Bay with high evenness and patchiness (Jones et al., 2021). This has strong functional implications for carbon and nitrogen cycling (Knox et al., 2018; Windham-Myers et al., 2018; Baldocchi, 2020; Russell et al., 2023), habitat provisioning for threatened and endangered species (smelt, rails), recreation, and aesthetics.