Shinnecock Bay is located on the south shore of Long Island, NY ( Figure 1A ). It is part of a series of interconnected barrier-island lagoons that make up the South Shore Estuary Reserve, including, from west to east, Hempstead Bay, South Oyster Bay, Great South Bay, Moriches Bay, and Shinnecock Bay, which is the easternmost bay. Shinnecock Bay is microtidal (tidal range = 1 m), covers a ~40 km² area and has a constriction in the middle that separates the bay into distinct eastern and western portions ( Figure 1B ). An inlet to the Atlantic Ocean in the eastern half of the bay provides strong tidal flushing to this portion of the bay that attenuates to the west.

A series of criteria were set for the establishment of adult hard clam spawner sanctuaries in Shinnecock Bay. Regions that had previously been studied and shown to support robust clam growth, survival, and/or condition were clear priority locations. Regions that experienced extended periods of high temperatures (>25°C) and extreme levels of the harmful brown tide alga, Aureococcus anophagefferens (> 10⁶ cells mL^-1) were avoided due to adverse effects of these conditions on M. mercenaria (Ansell, 1968; Hamwii and Haskin, 1969; Greenfield et al., 2005; Gobler and Sunda, 2012). Sediments composed of extremely coarse sediments or anoxic, sulfidic muds, which can be harmful to bivalves (Wang and Chapman, 1999), was avoided. An exceedingly thick covering of macroalgae, which can be found in some regions of the Bay (Young and Gobler, 2016), was avoided assuring clams would not be smothered in seaweeds that might promote anoxia (Valiela et al., 1992). Regions with high levels of channel whelk (Busycotypus canalicultatus), a well-known predator on the south shore of Long Island (Polyakov et al., 2007), were also avoided. Given that hard clams have a two-week larval stage (Carriker, 2001), sites with extended residence times, located away from the Shinnecock Inlet were selected to ensure the retention of larvae within the estuary. Finally, sites selected based on these scientific criteria required regulatory and public approval.

Within the spawner sanctuaries established in Weesuck Creek and Tiana Bay, adult clams were transplanted at high densities (mean ± SD = 27.1 ± 4.7 m^-2) onto numerous smaller plots (~1,850 m²) to increase fertilization efficiencies and the production of larvae, instead of being spread out at low densities across the entire area. In total, 30 GPS-delineated plots were established in Weesuck Creek, and 34 plots were established in Tiana Bay, spaced in a near equidistant manner from one another with the borders of each sanctuary area defined and GPS-delineated. Between 25,000 to 83,000 clams were planted in each plot yielding densities of 12.3 - 44.6 clams m^-2, with a mean planting amount and density of 50,300 clams and 27 clams m^-2, respectively. In total, 3,219,250 clams were planted from 2012 – 2019.

Sanctuary plots were planted with adult, market-sized (i.e., ≥ 25.4 mm shell width) clams purchased from seafood distributors. The clams were harvested from New York and Connecticut waters, primarily in the Long Island Sound where most of the current hard clam fishery exists for these states ( Figure 1A ). The clams were harvested by hand with clam rakes, sold to local distributors on the day of harvest, and stored under refrigeration by the distributors for up to three days prior to transport to vessels (R/V Peconic or R/V Shinnecock) for planting. To plant clams, the research vessel slowly (~ 1-2 knots) moved across the entire latitudinal and longitudinal confines of the sanctuary plot boundaries defined by the four corner GPS coordinates while clams were gently, evenly, and continually dispersed into the water from the deck of the pontoon vessels. Any visibly dead or broken individuals were removed to minimize the attraction of predators (Polyakov et al., 2007). Clams were planted during spring and fall when water temperatures were above 10°C to ensure clams would be active enough to dig into sediments. Planting ceased during summer months to avoid stressful temperatures (>25°C; Ansell, 1968; Hamwii and Haskin, 1969; Greenfield et al., 2005).

Three plots in each of the two sanctuaries were surveyed annually during summer to quantify hard clam abundance and size. At each plot, divers hand excavated 20 randomly placed quadrats (0.25 m²) to a depth of approximately 15 cm. The shell height, length, and width of recovered clams were measured with calipers to the nearest 1 mm.

As a proxy for hard clam gametogenesis, clam gonad condition was monitored twice per year, near the beginning (June) and end (September) of the hard clam spawning season from 2013 - 2019. In 2020, gonad condition was monitored more frequently during summer (~weekly in June and July) to generate a refined understanding of the timing of spawning. On each sample date, ten clams were sampled from each of three plots in each sanctuary that were also monitored for survivorship (n=30 clams per sanctuary). The clams were either hand excavated by SCUBA divers (2012-2019) or collected using a clam rake from a boat (2020). Gonad ripeness (GR) was monitored qualitatively using a visual ranking system that ranged from 0 (no gonad) to 4 (most ripe), with declines in gonad ripeness providing an indication of spawning (Doall et al., 2008). Rankings were based on the relative amount of gonad tissue observed within the dissected visceral body, as well as the relative number and appearance of gametes observed in a tissue smear through a compound microscope as described by Doall et al., 2008.

Benthic surveys of ~60 randomly generated sites outside of the sanctuary plots were conducted by SCUBA annually during late summer from 2012-2021 to track changes in the abundance and size distribution of clam populations throughout eastern and western Shinnecock Bay. All sites were at least 200 m apart. At each site, divers excavated two quadrats (1 m²) to a depth of 15cm using a suction dredge to vacuum all quadrat contents into a 3-mm mesh collection bag. The contents of collection bags were sieved through a 3-mm screen on the boat, and all hard clams and crabs (mud crabs, Dyspanopeus sayi, blue crabs, Callinectes sapidus, lady crabs, Ovalipes ocellatus, spider crabs, Libinia emarginata, green crabs, Carcinus maenas, broad-clawed hermit crabs, Pagurus pollicaris, and long-clawed hermit crabs, Pagurus longicarpus) greater than this size were counted and measured for shell length, height, and width with calipers. Pre-market size clams (<25.4 mm shell width) were binned into one of three size categories (shell width) based on size-age relationships established for clams in the hydrologically-connected lagoon, Great South Bay (Kraeuter et al., 2005): 1) < 9 mm (0-year-olds); 9-17 mm (1-year-olds); 18-25 mm (2-year-olds). Market size clams were binned into four size categories corresponding to the commercial size classes of littleneck (26-35 mm), topneck (36-38 mm), cherrystone (39-41 mm), and chowder (>41 mm). The size ranges for each commercial size category were adapted from Kraeuter et al. (2005), with the addition of the ‘topneck’ size category added between littlenecks and cherrystones to correspond with NYSDEC commercial landings data.

Digital PCR was used to detect and quantify the concentration of hard clam larvae in water samples collected via twice weekly cruises across Shinnecock Bay during June and July 2020, when hard clams are known to spawn on the south shore of Long Island (Doall et al., 2008). A total of 14 cruises were performed, two per week, from 6/10/2020 – 7/27/2020, sampling from seven stations spread across eastern and western Shinnecock Bay, including: Far East (40.86386, -72.44838), East Central (40.8659, -72.4718), Bridge East (40.85854, -72.49019), Bridge West (40.84375, -72.51406), Tiana (40.85523, -72.5381), Weesuck (40.83698, -72.5678), and Far West (40.81905, -72.58297) ( Figure 1 ). At each site, two replicate 20L water samples were obtained from ~0.5m below the surface and were concentrated onto a 20µm sieve, the contents of which were backwashed into two, 50mL centrifuge tubes vials that were immediately placed on ice and stored at -20°C within two hours of collection.

To prepare samples for digital PCR analysis, total cellular genomic DNA was extracted using the Qiagen Blood and Tissue Kit (Qiagen, Germantown, MD, USA). Samples were filtered onto 0.2 µm polycarbonate filters that were transferred to a 5 mL bead tube and treated with a lysis buffer and proteinase K to chemically lyse all cells and remove non-DNA organic contents. The lysate was transferred to a new tube and 50% v/v solution of ethanol was added to precipitate total genomic DNA. The mixture was subsequently passed through a mini-spin column via centrifugation (13000 g) to capture the DNA. The columns were further washed with a buffer containing ethanol, followed by elution of high-quality DNA with 100 µl nuclease-free water. The eluted samples were analyzed on a Nanodrop Spectrophotometer (Thermo Scientific, Waltham, MA, USA) to ensure nucleic acid recovery and quality. The purified DNA samples were stored at -80°C until digital polymerase chain reaction (dPCR) analysis.

Digital PCR analysis was conducted using the chip-based Applied Biosystems™ QuantStudio™ 3D Digital PCR System (Applied Biosystems, Foster City, CA, USA) to quantitatively identify the number of mitochondrial cytochrome c oxidase subunit I (COI) gene copies from M. mercenaria. To identify unique regions of the COI gene for this species, sequences for various bivalve species found across the US East Coast, including the closest relative Gemma gemma, were obtained from publicly available databases, NCBI GenBank and Barcode of Life Data Systems (BOLD; Ratnasingham and Hebert, 2007). The sequences were aligned, and unique regions were used to design forward (GGGCTCTTTCTCATTCTGGTAG), and reverse (AGCAACACCATAACACCCG) primers, and probe (TGGGTGGTGCATCTTCTATTTTGGCG). Each assay was validated and optimized using the dPCR system prior to sample analysis using synthetic double-stranded DNA fragments of the target genes as standards (gBlocks, Integrated DNA Technologies, Coralville, IA, USA). Lyophilized gBlocks were resuspended in 25 μL of IDTE buffer + 100 ng/μL polyA carrier (Roche Life Sciences, South San Francisco, CA, USA) used to increase the recovery of synthetic standards (Miyaoka et al., 2016), quantified using a Qubit 4 Fluorometer (Thermo-Fisher, Waltham, MA, USA), and serially diluted to prepare standards with final concentrations of 800 copies µL^-1. Optimization trials, testing gradients of annealing temperature, primer-probe concentrations and numbers of cycles were conducted to identify optimal thermocycling conditions for the assay.Digital PCR amplifications were performed in 14.5 µL reaction mixtures consisting of 7.25 µL of QuantaStudio 3D digital PCR Master mix v2 (2x stock solution), 0.725 µL Taq Man assay primer and probe mix (20x stock solution, see Table 5.1 for final concentrations), 1.525 µL nuclease free water and 5 µL sample DNA. All samples were originally analyzed using a maximum 5 µL of extracted DNA to achieve an on-chip concentration in the optimal range of 200-2000 c/µL. The dPCR reactions were loaded onto QuantStudio™ 3D Digital PCR Chip V2 chips containing 20,000 well partitioning with the QuantStudio™ 3D Digital PCR Chip loader (Applied Biosystems, Foster City, CA, USA), sealed with immersion fluid and the chip lid per the manufacturer’s instructions. All chip preparation was performed in less than one hour per manufacturer’s recommendations to minimize degradation. Loaded chips were amplified using a ProFlex™ 2x Flat PCR System thermocycler (Applied Biosystems, Foster City, CA, USA). Thermocycling conditions included initial denaturation at 95°C for 10mins, followed by 40 cycles of three steps: 95°C for 10mins, 95°C for 10mins and 95°C for 10mins. Amplified chips were brought to room temperature to prevent condensation before imaging on the QuantStudio™ 3D Digital PCR instrument (Applied Biosystems, Foster City, CA, USA). All samples were run as duplicates, along with a negative (nuclease free water) and positive (gBlock standards, 800 copies µL^-1 concentration) control.Imaging data derived from the QuantStudio™ 3D Digital PCR instrument were analyzed using the Applied Biosystems QuantStudio^® 3D AnalysisSuite™ cloud software. This software provided quality control steps on a per chip basis determining wells suitable for further analysis. In this study the default quality threshold of 0.5 was used for all chips. Chips were also manually inspected for equal distribution of positive wells across the chips and chip damage, such as large bubbles or evaporation, resulting in loss of readable wells in which chips were omitted and the sample was reanalyzed. Software derived fluorescence (call) thresholds delineating the unamplified wells (negative calls) and amplified wells (positive calls) were manually reviewed for each chip and adjusted to a common threshold per assay based on the ranges of the positive control and negative control clusters. Additionally, the spread of reads along the secondary assay (non-target dye) was manually reviewed in which wells identified as positive located largely outside the range of the positive control clusters on the secondary axis were identified as no amplification to reduce false positives. The negative and positive well counts were then converted to absolute quantification (copies µL^-1) by the software using Poisson statistics and corrected for dilution/concentration factors during sample collection (filtration), DNA extraction, and PCR reaction preparation. Sample concentrations are reported in M. mercenaria COI gene copies per liter.

Discrete water samples were collected weekly-to-monthly from May through October from selected locations in Shinnecock Bay from 2012-2021, with weekly sampling in summer (June – September) and monthly in spring and fall. At each station, temperature and salinity were measured at the surface and near the bottom using a YSI^® 556 sonde to confirm that the water column, which was typically 2 m, was well-mixed (Wilson et al., 1991). Water samples were collected in 20-L carboys at ~0.5 m at each station and returned to the laboratory for further analyses within 30 minutes. Triplicate chlorophyll a samples were collected on GF/F glass fiber filters (pore size = 0.7 µm), frozen, and analyzed by standard fluorometric methods (Parsons et al., 1984). Densities of the brown tide alga, Aureococcus anophagefferens, were quantified on a flow cytometer using a species-specific immuno-assay on water samples preserved in 1% glutaraldehyde (Stauffer et al., 2008).

A 20-year (2000-2020) dataset compiled by the Suffolk County (NY, USA) Department of Health Services (SCDHS) water quality monitoring program was analyzed to compare conditions before and after the creation of hard clam spawner sanctuaries. This program samples dozens of sites across the south shore of Long Island on a biweekly-to-monthly basis. Sampling occurred monthly and during the day without regard for tidal state. Sites are evenly distributed across Shinnecock Bay, typically separated by 1 – 3 km ( Figure S1 ). Chlorophyll a and nutrients were measured using standard oceanographic methods (Parsons et al., 1984), and Secchi disk depths were quantified. Water quality conditions were compared pre- and post-hard clam spawner sanctuaries and over time.

New York has lost 90% of its seagrass meadows since the 20^th century (NYSDEC, 2009). Here, classification of orthoimagery was performed to generate an estimate of eelgrass (Zostera marina) coverage change from 2013 to 2017. The imagery was obtained from the USDA’s National Agriculture Imagery Program (NAIP) which has collected high-resolution four-banded aerial images of New York state every other year since 2011. NAIP imagery has a 1-meter ground sample distance resolution, meaning a raster cell size of 1m by 1m. The NAIP imagery is collected to survey agricultural land during peak growing conditions (June – August) which fortuitously aligns with the peak eelgrass growing season on Long Island, NY. Since the intention of the imagery collection does not include surveying coastal waters, poor sea state conditions and sun angles can obscure the benthos. The 2011 and 2019 imagery were completely obscured and unusable, the 2021 imagery is not yet available, and therefore, the 2013 and 2017 orthoimagery were used for this analysis.

The western side of Shinnecock Bay, where spawner sanctuaries were established, was selected as the target region to assess eelgrass coverage change. Classification of eelgrass coverage change was partly constrained by the 2013 imagery as sun glare obscured the imagery on the western and easternmost edges of the bay, further reducing the viable area for analysis ( Figure S2 ). Depths beyond two meters were excluded as target regions as 99% of eelgrass in NYS is shallower than this depth (NYSDEC, 2009) and deeper regions become unusable for imaging of benthic structures.Natural color bands from the NAIP imagery (bands 1-3, red, green, blue) were used to classify bay bottom into four categories (shallow non-seagrass substrate, deep non-seagrass substrate, shallow seagrass, and deep seagrass) through a supervised object-based support vector machine via Esri’s ArcGIS Pro 2.9.1 Image Classification Wizard tool following specifications for object creation described by Davenport et al. (2017). Forty training samples were chosen for each class with a minimal pixel coverage of 20px by 20px. After model training and classification, the deep and shallow seagrass classes were combined, extracted, and converted to a polygon feature class to allow for geometry calculations. A union was performed on the 2013 and 2017 seagrass polygons. Areas of overlap in eelgrass coverage between the two years were assigned as no change for eelgrass coverage, areas with eelgrass only from the 2017 polygon were assigned as a gain in eelgrass coverage, and areas from solely the 2013 polygon were assigned as losses in eelgrass coverage. The total area was then calculated for each of the three assigned categories of eelgrass coverage. Ground-truthing was performed via comparison of model outcomes to diver surveys to affirm Z. marina coverage in the predicted locations.

A previously developed bio-optical model based on the absorption and scattering of light by specific components of coastal waters (Gallegos and Kenworthy, 1996; Gallegos and Biber, 2004 2008) was calibrated for Shinnecock Bay and the characteristics of the dominant seagrass in Shinnecock Bay, Zostera marina. The optical model is based on the properties of three optically active constituents: non-algal particulates (TSS), phytoplankton (chlorophyll a), and colored dissolved organic matter (CDOM). Chlorophyll a was quantified as described above, TSS was determined gravimetrically (Parsons, 2013) and CDOM was measured by filtering seawater through a 0.7µm glass fiber filter and measuring levels on a Trilogy Laboratory Fluorometer with a CDOM module with ultraviolet GUI selection. Thirteen sites across Shinnecock Bay ( Figure S3 ) were sampled during the day at peak water temperatures (August), a time that may be stressful for the seagrass, Z. marina (Hammer et al., 2018). Measurements made in August 2017 were compared to August 2021. Previous studies have demonstrated that light levels below a minimum physiological requirement (typically 15 – 40% of incident surface light) result in depth-limits to seagrass (Kenworthy and Fonseca, 1996; Steward et al., 2005). Therefore, the calibrated optical model was inverted to set threshold concentrations for water quality parameters to meet the water column light requirement for Z. marina at a given water depth (Gallegos and Moore, 2000).

The hydrodynamic model, Environmental Fluid Dynamics Code (EFDC; Hamrick, 1992) was constructed and optimized for Shinnecock Bay to estimate the movement of water and bivalve larvae in this system. The model solves the three-dimensional, time-dependent equations of motion for variable-density water. An equation of state, which determines water density as a function of temperature and salinity, is used to couple conservation of mass (salinity) and energy (temperature) equations to the hydrodynamic equations of motion. Vertical mixing due to turbulence is simulated using the Mellor-Yamada turbulence closure scheme (Mellor and Yamada, 1982) as modified by Galperin et al. (1988). Time-dependent boundary conditions in EFDC can be specified at inflow and outflow boundaries. A complete description of the EFDC model, including numerical algorithms is provided by Hamrick (1992) and Hamrick (1996).

For this study, the spatial extent of the model included Shinnecock Bay, the Atlantic Ocean, and Moriches Bay ( Figure S4 ). Inputs to the hydrodynamic model were obtained from NOAA tidal gauges in Moriches and Shinnecock Bays at their connections to the Atlantic Ocean and weather data from the Francis S. Grabreski Airport. Bathymetric data were obtained from the NOAA Coastal Relief Model (NOAA, 1999) and sounding surveys. The hydrodynamic model grid comprises 901 total cells. The aerial extent of the grid is about 40 km (25 miles) in the east-west direction, and 5 km (3 miles) in the north-south direction ( Figure S4 ). The average dimensions of the cells are 370m (east-west) and 200m (north-south). The bays are relatively shallow (1 – 4m), and, in the absence of any observed vertical stratification, were simulated as a single vertical layer. The hydrodynamic model was calibrated to water surface elevations and temperature measured by NOAA at the Ponquogue Point station in Shinnecock Bay and the Moriches CGS station in Moriches Bay. Temperature data were obtained from the Suffolk County monitoring stations in Shinnecock, Quantuck, and Moriches Bay (SCDHS, 1985-2021).In the first iteration of the model, a conservative tracer (i.e., dye) was released from each spawner sanctuary to visualize water dilution and movement on the time scale of the larval life cycle (two weeks; Carriker, 2001) using physical conditions (winds, tides, temperature, salinity) found for Shinnecock Bay during summer when M. mercenaria spawns. In a more refined iteration, the transport of neutrally buoyant particles that mimic clam larvae was calculated using velocity vectors based on a grid resolution of approximately 250 meters and a time interval of 5 minutes. Discrete particle locations were computed independently of the hydrodynamic model grid extents and used the velocity magnitude and direction from the nearest model cell. This type of analysis provides hypothetical particle tracks that may vary depending on grid size, computation time interval, release locations, release times, and complexity of circulation patterns. Particle movement was tracked using the tidal and wind conditions present during July 2020, when a high abundance of larvae was quantified in Shinnecock Bay. Particle movement was determined for both spawner sanctuaries and for eastern Shinnecock Bay, with one particle placed in each grid cell for each region (14 in Weesuck Creek, 17 in Tiana Bay, and 3 in eastern, central Shinnecock Bay).

The commercial landings of hard clams from Shinnecock Bay and near-by lagoonal estuaries were tracked using long term data sets publicly available from the New York State Department of Environmental Conservation (NYSDEC). In New York State, the legal harvest size limit for hard clams is one inch (25.4 mm) shell width, which represents three years of growth (Kraeuter et al., 2005). Landings are reported to the NYSDEC by shellfish distributors, including information on quantity, date of harvest, harvest area, and clam size. Clam size is classified into one of four commercial size categories, including, from smallest to largest, littleneck, topneck, cherrystone, and chowder, with cherrystone and chowder size clams often reported together as a “large mix”. Prior to 2012, the NYSDEC tracked only total hard clam landings of all size categories combined, with landings expressed in terms of bushels. From 2012-2020, more detailed landings data is available from the NYSDEC, including information on the number of individuals harvested in each size category. To convert between individuals and bushels, the NYSDEC uses the following conversion factors: one bushel equals 600 littlenecks, or 200 topnecks, or 150 cherrystones, or 100 chowders, or 120 large mix (cherrystone and chowder mixed). Harvest locations are reported according to specific areas delineated by the NYSDEC, and from 2008 to 2021 separate harvest areas were delineated for western (SS12) and eastern (SS13) Shinnecock Bay. Prior to 2008, landings from the entire Shinnecock Bay and the adjoining smaller lagoon Quantuck Bay were grouped into one harvest area.

To quantify the potential filtering effect of M. mercenaria populations in Shinnecock Bay, biological filtration rates of clams within eastern and western Shinnecock Bay were determined from 2014 - 2021. Filtration capacities were derived from allometric relationships established between M. mercenaria weight and clearance rates for larger clams (>16 mm) fed natural seston (Hibbert, 1977; Doering and Oviatt, 1986) and for smaller clams (<16 mm) fed natural seston (Grizzle et al., 2001). Application of these relationships to previously reported clearance rates of M. mercenaria feeding on natural seston from Shinnecock Bay (Wall et al., 2008) yielded rates consistent with those measured empirically in that study. Size-specific clearance rates from the allometric relationship were applied to the densities of each surveyed size class of clam surveyed in each side of the bay and the clams planted in sanctuaries each year. These basin-wide clearance rates were further applied to the total water volumes of eastern and western Shinnecock Bay as determined via GIS and the NOAA Coastal Relief Model (NOAA, 1999) to generate basin-wide turnover times.

Two areas within western Shinnecock Bay were chosen for the location of hard clam spawner sanctuaries: Tiana Bay and Weesuck Creek ( Figure 1B ). Prior studies documented high survivorship and condition index of clams planted in these regions (Padilla and Doall, 2012). Hydrodynamic models indicated that dye released from Weesuck Creek and Tiana Bay in western Shinnecock Bay dispersed to regions east and west of these sites over the two-week larval period of M. mercenaria ( Figure 2 ). While both sites had a large percentage of the dye retained in western Shinnecock Bay, dye from Tiana Bay mixed across nearly all of eastern Shinnecock Bay and dye released within Weesuck Creek mixed across a portion of eastern Shinnecock Bay ( Figure 2 ), suggesting larvae spawned in these sites would be spread across, and retained within, the Bay. SCUBA surveys found sanctuary regions to be comprised of appropriate sediment types (sandy mud), with moderate levels of macroalgae (<25% coverage), and an absence of channel whelk. Water quality data (SCDHS, 1985-2021) revealed these locations experienced more moderate temperatures (<25°C) and harmful brown tides (<10⁶ cells mL^-1) than western regions in this estuary. Following approval by the NYSDEC, Southampton Town Trustees and Southampton Bayman’s Association, precise areas within Tiana Bay and Weesuck Creek were designated as harvest-free sanctuaries in 2012 ( Figure 1B ).

Adult hard clam plantings began in 2012. By the spring of 2015, more than one million clams had been planted ( Figure S5A ). Two and three million clams were planted by 2016 and 2017, respectively, while a slower addition of tens to hundreds of thousands of clams were planted annually from 2018-2019 ( Figure S5A ). By the end of 2019, 1,540,200 clams had been planted into 30 plots in the Weesuck Creek spawner sanctuary, and 1,679,050 clams had been planted into 34 plots in the Tiana Bay sanctuary, at densities ranging from to 12.3 to 44.6 clams m^-2 (Figures S5B,C). Clams planted were littleneck size and larger (>25mm).

In Weesuck Creek and Tiana Bay, hard clam densities in spawner sanctuaries have remained consistent across years ( Figure 3A ). Almost all hard clams collected during surveys were large adults (>39 mm shell width), with low juvenile abundance indicating low recruitment within the sanctuary areas ( Figure S6 ). Therefore, the steady hard clam population size at the planting sites indicated high survivorship of the adults that were planted, rather than new recruitment replacing losses due to mortality. In addition, a decrease in mean gonadal ranks during the hard clam spawning season (summer) indicated annual spawning activity within spawner sanctuaries ( Figure 3B ). Combining years and sites, mean gonad ranks in late spring (2.85 ± 0.28; mean ± SD) were significantly higher than in late summer (2.24 ± 0.28) (t=4.93; p<0.001; paired two sample t-test for means; Figure 3B ), with summer declines in gonadal rank averaging 21%. To further resolve the timing of spawning during summer, weekly monitoring was performed in 2020 which revealed high GR from late May through early July, but a large decline in GR thereafter, from 3.9 in early July to 2.7 by mid-July in Weesuck Creek, and from 3.7 to 2.9 in Tiana Bay, changes indicative of spawning ( Figure 4A ).

Quantification of hard clam larvae across seven sites in Shinnecock Bay revealed low levels (< 1 x 10⁴ gene copies L^-1) through June and early July, with small peaks within the eastern bay on June 10 and at isolated locations in early July ( Figure 4B ). On July 17^th, the highest density of hard clam larvae to that date (5.8 x 10⁴ copies L^-1) was quantified within Tiana Bay with a smaller peak level found at the Far East station (2.3 x 10⁴ copies L^-1; Figure 4B ). The timing of this peak corresponded with the timing of the largest decrease in gonadal rank of clams in the spawner sanctuaries ( Figure 4A ). Three days later (July 20), densities similar to those found in Tiana Bay on July 17^th were quantified to the east at the sites east and west of the Ponquogue Bridge in the center of the bay (6.2 x 10⁴ and 5.0 x 10⁴ copies L^-1, respectively), with high levels remaining in Tiana Bay but also appearing in the center of Eastern Shinnecock Bay (East Central; ~3 x 10⁴ copies L^-1 for both locales; Figure 4B ). Visual examination of formalin-fixed plankton samples from net tows confirmed the presence of high concentrations of bivalve larvae at the Ponquogue Bridge on July 20. Hard clam larvae concentrations declined, but remained at moderate levels, through late July with 1.2 x 10⁴ copies L^-1 at the East Central site on July 24 and an increase in the Far Eastern site on the same date ( Figure 4B ).

To mimic the sudden increase and presumed release of larvae from the spawner sanctuaries on July 17^th, the hydrodynamic model was used to release neutrally buoyant particles based on tidal and wind conditions present on July 17, 2020. Hydrodynamic transport conditions near the mouth of spawner sanctuaries were relatively small in a north-south direction moving from the mouth out into Shinnecock Bay where velocities were relatively higher and predominantly in the east-west direction with the highest velocities occurring along the southern side of Shinnecock Bay. Over a two-week period, particles migrated east from Tiana Bay and Weesuck Creek, past the Ponquogue Bridge, and settled in eastern Shinnecock Bay, with the majority of Tiana Bay-released particles migrating to eastern Shinnecock Bay and Weesuck Creek-released particles displaying both eastern and western dispersal ( Figures 5A, B ). In contrast, particles released from eastern Shinnecock Bay were tidally transported into the Atlantic Ocean in four days ( Figure 5C ).

Benthic surveys of adult hard clams across the non-sanctuary zones of Shinnecock Bay revealed increases in densities from 2012 – 2021 with distinct patterns for the eastern and western Bay ( Figure 6A ). When benthic surveys began in 2012 in western Shinnecock Bay and 2015 in eastern Shinnecock Bay, clam densities were 1.2 and 0.4 clams m^-2, respectively, with populations being comprised primarily of large individuals > 36 mm ( Figure 6A ). In the subsequent years, both locations experienced steady increases in densities caused almost exclusively by smaller size classes of hard clams ( Figure 6A ). In the western half of the bay, clam density increased 50% between 2012 and 2021, with average densities of 1.8 and 1.9 clams m^-2 found in 2020 and 2021, respectively. In the eastern half of the bay, increases in clam densities were more dramatic, with increases occurring from 0.4 clams m^-2 in 2015 to 1.3 clams m^-2 in 2016, to 1.7 clams m^-2 in 2018, to 3.0 clams m^-2 in 2019, 3.5 clams m^-2 in 2020, and 7.3 clams m^-2 in 2021. In total, clam densities increased 18-fold from 2015 to 2021 in eastern Shinnecock Bay with a clustering of sites across the center of the eastern Shinnecock Bay displaying 6 – 36 clams m^-2 by 2021 ( Figure 6B ). The higher recruitment in eastern Shinnecock Bay has resulted in a dramatic shift in the population size structure ( Figure 6A ). For example, large clams (>36 mm; topnecks, cherrystones, and chowders) represented 60% of the population in the eastern Bay in 2015 but were <10% of the population by 2021 as small clams (<35 mm; littlenecks and juveniles) increased from 40% to 90% of the population ( Figure 6A ). The most dramatic change was among smallest size classes of clams (< 17mm) that were not detected in 2015 but comprised more than a one quarter of clams surveyed in 2020 and 80% in 2021 ( Figure 6A ). In western Shinnecock Bay, a smaller shift in population size structure was observed, with the dominance of small juvenile clams (<17 mm) increasing from 2% in 2012 to 25% in 2021.

Commercial landings provide another indicator of recent hard clam population changes in Shinnecock Bay. Annual hard clam landings from Shinnecock Bay have increased since the onset of the planting of adult clams in 2012 by more than sixteen-fold from 878 bushels in 2012 to 13,983 bushels in 2021 (p<0.0001; linear regression with time), with a modest rate of increase from 2012 – 2015 and a more rapid rise from 2016 – 2021 ( Figure 7A ). This increase has been stronger in the eastern half of Shinnecock Bay than the west ( Figure 7A ) and collectively the landings in 2021 were at a level (>14,000 bushels) not seen in Shinnecock Bay and surrounding estuaries since 1981 ( Figure 7B ). This sharp rise in landings has been regionally specific to Shinnecock Bay as landings in other, interconnected south shore Long Island lagoons have declined or had much smaller increases during the same time period ( Figure S7 ). In the neighboring Moriches Bay, for instance, landings declined by 42% from 2012 to 2021 (Figure S7). While landings in Great South Bay increased from 2015 – 2021, landings in Shinnecock Bay now exceed those in Great South Bay, which is approximately ten-times the aerial size of Shinnecock Bay ( Figure 1 ; Figure S7 ).

Among the increased landings within Shinnecock Bay, the increase has been disproportionately among smaller clams. From 2012 to 2020, annual landings of littleneck clams increased faster (26-fold increase) than annual landings of larger cherry and chowder clams (seven-fold increase), indicating that the total increase in landings was largely due to recent population increases rather than solely increased fishing pressure ( Figure 7C ). With the disproportionate increase in landings of littleneck clams, the percentage of littleneck clams in total landings increased from 45% in 2012 to 67% in 2020 ( Figures 7C, D ).

The combined clams planted, and clam densities measured across the Shinnecock Bay provides an estimate of the total number of hard clams in eastern and western Shinnecock Bay ( Figure 8A ). In western Shinnecock Bay, this number rose significantly with time (p<0.01) from 1.9 x 10⁷ clams in 2014 to 4.2 x 10⁷ clams in 2021 whereas in the eastern portion of the bay the population rose from 0.6 x 10⁷ clams in 2015 (first year of full survey) to 1.1 x 10⁸ clams in 2021 (p<0.01; Figure 8A ). Increased population sizes on each side of the bay translated into a more rapid ecosystem clearance rate and thus shorter biological turnover time of the bay over time when applying allometrically determined, size-specific clearance rates of hard clams to the in-situ population ( Figure 8B ). The western portion of Shinnecock Bay had a turnover time of three weeks in 2014 but this rate dropped to 11 days by 2017 and dropped to 10 days in 2020 ( Figure 8B ). Changes in the deeper, eastern portion of the bay were more dramatic with a biological turnover time of 104 days in 2015 that progressively dropped to 18 days by 2020 ( Figure 8B ). Applying similar calculations to the sub-bay regions hosting hard clam spawner sanctuaries indicates biological turnover rates in Weesuck Creek and Tiana Bay were initially two-to-three weeks and are now both ~ three days.

Since the commencement of hard clam restoration in Shinnecock Bay, there have been distinct changes in several water quality parameters. Brown tides caused by A. anophagefferens had been annual events on the south shore of Long Island including Shinnecock Bay from 1985-2016 ( Supplementary Table 1 ). From 2007 – 2016, brown tides with peak densities from 0.3 – 1.9 x 10⁶ cells mL^-1 had occurred annually in this system with particularly dense and sustained blooms detected in 2012 and 2016 (> 1.5 x 10⁶ cells mL^-1; more than six weeks; Figure 9 ). Since 2016, a brief, small brown tide bloom (~ 1 x 10⁵ cells mL^-1; two weeks; Figure 9 ) occurred in 2017, with no brown tides from 2018 – 2021 (<10⁴ cells mL^-1; Figure 9 ). Since 2012, there has been a significant decrease in the densities of A. anophagefferens and chlorophyll a concentrations measured in western Shinnecock Bay (stations 180 and 190; p=0.002 and p=0.0001; linear regressions with time) during the months of April through November, when clams are actively feeding (Figures S8A,B). There has also been a significant increase in the Secchi disc depth of western Shinnecock Bay (p=0.03; linear regression with time; Figure S8C ). In contrast to Shinnecock Bay, other interconnected south shore lagoons not experiencing restoration have experienced significant increases in chlorophyll a and harmful brown tides caused by A. anophagefferens, as well as a significant decrease in water clarity during the same time period (Gobler et al., 2019).

Digital aerial imagery of seagrass in the western portion of Shinnecock Bay revealed that, among the areas visible via aerial data, there was an expansion of aerial coverage from 2013 through 2017 of 5.75 x 10⁵ m², as well as loss of 1.01 x 10⁵ m², yielding a net gain of 4.74 x 10⁵ m² ( Figure 10A ). Gains in seagrass were primarily concentrated in the western extent of the surveyed region, but also occurred within and near the Tiana Bay spawner sanctuary ( Figure 10A ). Sampling of parameters used in the bio-optical model demonstrated that there has been a substantial increase in light availability in western Shinnecock Bay between 2017 and 2021 that has been driven by the reduction in chlorophyll a in this portion of the bay ( Figure 10B ; Table S2 ). Light levels within the far western region and specifically in the Weesuck Creek and Tiana Bay spawner sanctuaries have been estimated to increase by up to 10% while light levels in some regions to the east have decreased, potentially associated with the transport, deposition, and resuspension of sediments recently deployed to fortify the barrier island near the ocean inlet.

Multiple lines of evidence demonstrate that the significant increase in the hard clam population in Shinnecock Bay can be attributed to the establishment of hard clam spawner sanctuaries. The more than 3.2 million hard clams planted in western Shinnecock Bay displayed undetectable mortality and patterns of gonadal condition indicative of annual spawning during the past decade (Doall et al., 2008). Physical circulation models demonstrated that larvae in the western part of the bay would be transported to eastern Shinnecock Bay within the two-week duration of the larval stage of M. mercenaria (Castagna, 2001) and also indicated that larvae spawned in eastern Shinnecock Bay would be tidally advected out to the Atlantic Ocean in less than a week. The quantification of hard clam larvae in Shinnecock Bay demonstrated that peak larval densities were found within the Tiana Bay spawner sanctuary at the time hard clams were spawning, and then appeared eastward by the Ponquogue Bridge three days later, and then in eastern Shinnecock Bay days after that, a pattern fully consistent with that shown by hydrodynamic models.

The timing of increases in recruitment and landings following the establishment of spawner sanctuaries also suggests that the spawner sanctuaries contributed to these increases. The planting of adult clams began in 2012 and rapidly continued through 2017, totaling 2 million clams planted by 2016 and 3 million by 2017. Increases in clam densities, particularly in the eastern Bay, were largely due to increases in small clams (<25 mm; < 3 years) produced during the time frame in which the spawner sanctuaries were established. Increases in commercial landings began to accelerate in 2016, four years after the initial establishment of the spawner sanctuaries and was driven primarily by an increase in landings of small clams. As previously noted, it takes ~three years for hard clams to reach harvestable size in waters on the south shore of Long Island (Kraeuter et al., 2005). While increased landings can simply be a function of increased harvest pressure, the differential changes in the size of clams landed (i.e., 26-fold increase in landings of littleneck clams versus 7-fold increase in landings of cherrystone and chowder sized clams) was also supportive of the hypothesis that the increase in harvested clams since 2012 originated from population increases caused by higher recruitment. And while the increased landings could be related to broad-scale environmental conditions favoring the expansion of hard clam populations in the region, neighboring lagoonal estuaries across the south shore of Long Island displayed unchanged or decreased landings during the same period landings increased 16-fold in Shinnecock Bay.

Despite the location of hard clam spawner sanctuaries in the western half of Shinnecock Bay and despite the prediction that hard clam larvae from these locations would circulate to both eastern and western Shinnecock Bay, the trajectories of hard clam populations in these two regions were different. While both halves of the Bay saw increases in the densities of hard clams since restoration began, there was a 18-fold increase in the eastern bay but only a 50% increase in the western bay where a portion of larvae were predicted to have settled. Survival rates of M. mercenaria larvae are very low (<1%; Fegley, 2001) as this and other early life stages of clams are highly sensitive to environmental conditions including acidification (Talmage and Gobler, 2009; Talmage and Gobler, 2010; Talmage and Gobler, 2011) and harmful brown tide blooms (Greenfield and Lonsdale, 2002; Gobler and Sunda, 2012). From 2007-2017, maximal and average cell densities of the harmful brown tide alga (1.7 x 10⁶ and 6.97 x 10⁴ cells mL^-1, respectively) in western Shinnecock Bay were 11- and 24-fold times higher, respectively, than in eastern Shinnecock Bay (1.48x 10⁵ and 2.88 x 10³ cells mL^-1; SCDHS, 1985-2021) likely due to the strong tidal flushing with Atlantic Ocean water in the eastern Bay. The absolute densities of brown tide are impactful as early life stage clam feeding slows at 3.5 x 10⁴ cells mL^-1 and ceases at 2 x 10⁵ cells mL^-1 (Bricelj et al., 2001), levels found in the west during the summers of 2012 – 2016 but never seen in the east (SCDHS, 1985-2021). Moreover, juvenile clams experience near complete mortality when exposed to > 10⁵ cells mL^-1 of brown tide algae for several weeks (Greenfield and Lonsdale, 2002), a level common in the west during summers from 2007-2016, but never seen in the east. Similarly, cruise tracks measuring pH during the clam spawning period (June-July) have mapped normal pH levels in the eastern Bay (7.9 – 8.1; 400 – 500 µatm) but low pH (< 7.5) levels in the western Bay (Figure S9), that can significantly increase the mortality of larval hard clams (Talmage and Gobler, 2009; Talmage and Gobler, 2010; Talmage and Gobler, 2011). In addition, recent field studies surveying the abundance of larval bivalves across various lagoonal sites revealed that the combination of high temperature and low DO reduced larval abundance (Weinstock and Collin, 2021), and western Shinnecock Bay is warmer and has lower DO than the eastern bay during the larval spawning season (SCDHS, 1985-2021). Further, benthic surveys suggest that densities of crabs are higher in the western Bay compared to the east (Figure S10). Hence, higher predation rates on juvenile clams (Juanes, 1992) as well as environmentally driven lower survival rates likely have both contributed to a slower recovery of the hard clam population in the western Bay compared to eastern Shinnecock Bay. It is possible that as hard clam population densities slowly increase within the western Bay and turnover times continue to decline, water quality will continue to improve, making environmental conditions more favorable for the survival of early life stage clams. It is also possible that organic-rich, acidic muds that minimize survival of newly set hard clam are an obstacle for recovery (Green et al., 2004; Green et al., 2009) in which case buffering of sediments with crushed shell may be needed to promote restoration of clam populations in some western regions (Curtin et al., 2022).

For decades, it has been recognized that bivalves are ecosystem engineers that, when present at high densities, can shape the character of estuaries by filter feeding (Jackson et al., 2001; Lotze et al., 2006) by transferring organic matter from the water column into sediments, increasing denitrification rates (Wozny, 2021), light penetration, and sediment nutrient stores to the benefit of seagrasses like Z. marina (Carroll et al., 2008; Wall et al., 2008) that in turn serve as a critical habitat for shellfish and juvenile fish (Heck and Orth, 1980; McDevitt-Irwin et al., 2016). During the end of the twentieth century, Shinnecock Bay was emblematic of an estuary that suffered the consequences of the loss of hard clam populations as landings had declined by 99.4% from 1970 to 2011 (NYSDEC, 1970-2011), densities were < 1 clams m^-2 (Weiss et al., 2007), dense brown tides occurred annually, and seagrass populations had been eradicated from much of the western Bay (NYSDEC, 2009). The strategic placement of hard clam spawner sanctuaries in regions shown to facilitate robust gametogenesis (Doall et al., 2008; Padilla and Doall, 2012) with suitable water and sediment quality and physical circulation patterns conducive for larval retention assisted in repopulating Shinnecock Bay with hard clams as well as significantly increasing the influence of their filter feeding on the estuary. The shallow nature of Shinnecock Bay (mean depth ~2 m) helps maximize the influence of bivalve filtration on the water column compared with deeper, more voluminous estuaries (Gerritsen et al., 1994; Pomeroy et al., 2006). While biological filtration times were three weeks and three months in western and eastern Shinnecock Bay, respectively, prior to the commencement of this project, the stocking and regrowth of hard clams across the estuary reduced the biological turnover time of the bay to 10 and 18 days, respectively, while the turnover times within spawner sanctuaries were less than three days. Prior studies have suggested a biological turnover time of three days or less was ideal for establishing biological control of algal blooms (Cerrato et al., 2004; Wall et al., 2008). While this was not achieved outside of the spawner sanctuaries, we note that tidal flushing likely achieves this standard in the eastern bay and acts additively with clam filtration rates and other mortality processes (e.g. other bivalves, zooplankton, etc.) in the western bay to provide removal rates that have contributed toward significant reductions in the levels of phytoplankton biomass in western Shinnecock Bay (p<0.0001), significant increases in water clarity (p<0.001), and the eradication of harmful brown tides from 2017-2021. This is the first time there has been five consecutive years without brown tide in any estuary on the south shore of Suffolk County, NY, since these events began in 1985 (SCDHS, 1985-2021). During the same time period in which water quality improved in western Shinnecock Bay, regions of the south shore estuary lagoon system in NY, USA, not influenced by spawner sanctuaries in Shinnecock Bay nor ocean inlets (i.e., central Great South Bay) saw densities of brown tide, levels of chlorophyll a, or Secchi disc depths increase in part to regional increases in nitrogen loads (Gobler et al., 2019). During the past decade there have been intensifying watershed nitrogen loads to south shore estuaries (SCSWWP, 2020) accounting for the declining water quality in Great South Bay and leaving increases in clam population densities in Shinnecock Bay as the most parsimonious explanation for improved water quality there.

These improvements in water clarity in western Shinnecock Bay were coincident with both the restoration of hard clam populations across the estuary and the regrowth of more than 4.74 x 10⁵ m² of Z. marina meadows from 2013-2017. This time frame is coincident with the time period when the clam-induced turnover time of the western bay declined from 21 to 11 days and suggests the ecosystem-wide impacts of this restored, dense population of filter feeding hard clams contributed to this change. Bio-optical modeling also revealed distinct increases in light penetration in western Shinnecock Bay, specifically within the regions of hard clam spawner sanctuaries where biological turnover times decreased to three days. Given the bio-optical model also demonstrated the water column light requirement for Z. marina has been extended to 2.4 x 10⁵ m² of bottom in western Shinnecock from 2017 to 2021 (data not shown), it is anticipated that the availability of clear aerial images of western Shinnecock Bay from 2021 may reveal a further expansion of seagrass during the period of the bio-optical model, 2017-2021.

The restoration approach used for M. mercenaria in Shinnecock Bay was based on the premise that the population was recruitment limited and that improving fertilization success by densely co-planting millions of clams in individual sanctuaries with >25 individuals m^-2 would significantly improve reproductive output (Levitan, 1991; Mann and Evans, 2004). Similar approaches with other bivalves such as bay scallops (Argopecten irradians) have met with success in other regions of the US east coast (Peterson et al., 1996; Goldberg et al., 2000; Tettelbach et al., 2013; Tettelbach et al., 2015). The significant recovery of the hard clam population in Shinnecock Bay catalyzed the establishment of the Long Island Shellfish Restoration Program by the NYSDEC in 2018 using a highly similar approach to what has been reported here: creating small regions (~2,000m²) with high densities of hard clams (~25 m^-2) in no-harvest regions with water and sediment quality that is supportive of multiple life stages of hard clams and residence times that will allow larvae to settle in neighboring estuarine regions. Given the high mortality of juvenile hard clams (Juanes, 1992; Peterson et al., 1995; Gosselin and Qian, 1997; Munroe et al., 2015), it would seem that mass spawning events associated with spawner sanctuaries are required to overcome such mortality and is a preferrable restoration approach given the ability of each pair of adult hard clams to create millions of offspring annually (Kraeuter and Castagna, 2001). Moreover, the larger aerial coverage of larval settlement compared to that obtainable through manual planting of juvenile clams provides a greater likelihood of juvenile clams settling in regions hospitable for their survival. And while dense assemblages of bivalves in experimental and aquaculture settings that can lead to food limitation (Rheault and Rice, 1996; Wall et al., 2011), the planting densities used here were below such thresholds, a conclusion supported by the high survival and gametogenesis rates of clams within spawner sanctuaries.

Several attributes of the hard clam (M. mercenaria) made it an ideal choice for restoration of Shinnecock Bay. First, as a native species, there was a background, albeit diminished, population of individuals to build upon and there were not concerns associated with invasive species (e.g. zebra and quagga mussels in freshwater systems; Strayer, 2010). In addition, as a longer-lived bivalve (40+ yr), increases in M. mercenaria densities should be more stable than shorter-lived bivalves. For example, working with the shorter-lived bay scallop (Argopecten irradians; 2 yr), Tettelbach et al. (2013; Tettelbach et al., 2015) described a faster population recovery than was described here, but this shorter-lived population was less stable and has recently experienced sharp declines associated with climate change and other factors (Tomasetti, 2022).There are multiple attributes of shallow lagoons such as Shinnecock Bay that make them well-suited to the simultaneous restoration of bivalves and water quality. The shallow nature of this estuary (~ 2 m) maximizes benthic-pelagic coupling and, therefore, the ability of bivalves to exert maximal filtration pressure on this system as compared to deeper estuaries where such impact is minimized (Pomeroy et al., 2006). In addition, Shinnecock Bay is microtidal, with a typical tidal range of ~1 m, allowing bivalve filtration pressure to be maximal relative to tidal export of plankton (Banas et al., 2007) while also maximizing larval retention (Chen et al., 1997).Another key aspect of the repopulation of Shinnecock Bay with hard clams has been the use of spawner sanctuaries where clams cannot be harvested. No-take zones are a common and increasingly used restoration practice in marine ecosystems (Halpern, 2003; Lester and Halpern, 2008; Edgar et al., 2014; Gownaris et al., 2019; Grorud-Colvert et al., 2021). There have been multiple examples of hard clam populations being overfished in NY (Kraeuter et al., 2008) and elsewhere (Peterson, 2002). Given the extended time for clams to grow to harvestable size (at least three years; Kraeuter and Castagna, 2001), reductions in hard clam populations to recruitment limiting densities is a situation that may not be reversible without the introduction of new adult spawning stock (Peterson et al., 1996). The creation of spawner sanctuaries assures a continued supply of new recruits into ecosystems and, in some ways, mimics a situation present on Long Island in the middle of twentieth century when duck farms were prominent along south shore lagoonal estuaries (Great South Bay, Moriches Bay) and forced the closure of northern sections of bays to shellfishing due to fecal contamination (Ryther, 1954; Swanson et al., 2010). This allowed large adult clams to become established in closed regions that produced settlers that recruited and grew in the southern regions of these estuaries and contributed to the largest harvests of hard clams ever recorded in NY (McHugh, 1983). Through the latter half of the twentieth century, duck farms closed, fecal contamination was reduced, northern regions of estuaries were re-opened to shellfishing, and were subsequently overharvested (Swanson et al., 2010). With the loss of this spawning stock, the hard clam fisheries in these regions collapsed (Kraeuter et al., 2008). In some respects, the restoration approach used in Shinnecock Bay recreated the situation that supported the record setting harvests of the mid-twentieth centuries, with the difference being that today, spawner sanctuaries in Shinnecock Bay have been designated no-harvest zones for the purposes of restoration. This underscores the importance of political and stakeholder engagement, which is needed to establish and enforce such zones.

The longer life span of M. mercenaria (up to four decades; Kraeuter and Castagna, 2001) compared to other temperate bivalves (e.g., A. irradians, two years) makes it an ideal target for long-term restoration (Peterson, 2002). While this extended life span suggests M. mercenaria spawner sanctuaries could have a lasting impact on the Shinnecock Bay and similar coastal ecosystems, the rapid pace of climate change must be considered for restoration projects. The coastal waters of NY have been warming at two-to-three times global averages during summer since 1982 (Gobler et al., 2017; Griffith et al., 2019), and larval stages of M. mercenaria can experience enhanced mortality under elevated temperatures (28°C; Talmage and Gobler, 2011). This bivalve, however, has a geographic range that extends south into the Gulf of Mexico where individuals experience slowed growth associated with the warmest temperatures of the year (Kraeuter and Castagna, 2001). Hence, warming temperatures may alter the phenology of restored hard clam populations, with an increased prevalence of heat waves (IPCC, 2019) more likely to cause mortality in early life stage individuals (Talmage and Gobler, 2011; Stevens and Gobler, 2018) rather than simply progressive warming. Ocean acidification may be an even greater risk for restored hard clam populations as recent modeling estimated that the five-year risk of a 50% decline in M. mercenaria populations on the south shore of Long Island increases from 25% to 79% with an increase in pCO₂ levels from 400 μatm to 800 μatm (Grear et al., 2020). Revision of state and federally recommended water quality criteria to account for the multiple stressors associated with eutrophication and climate change may help improve the long-term viability of this fishery (Tomasetti and Gobler, 2020). Given that acidification in estuaries can be caused by eutrophication (Melzner et al., 2013; Wallace et al., 2014; Wallace and Gobler, 2021) regional efforts to combat excessive nutrient loading may also likely benefit hard clam spawner sanctuaries by lessening the intensity of harmful algal blooms (Heisler et al., 2008; Gobler et al., 2011) and the intensity of coastal acidification (Melzner et al., 2013; Wallace et al., 2014; Wallace and Gobler, 2021).

This study demonstrates the cascading, ecosystem-wide effects hard clam spawner sanctuaries can have on shallow lagoonal estuarine ecosystems. A decade of benthic and pelagic monitoring revealed that the stocking of no-harvest zones with dense assemblages of adult hard clams (27 m^-2) in areas with water and sediment quality supportive of gametogenesis and extended larval residence times caused a 18-fold increase in hard clam densities and a 16-fold increase in hard clam landings in Shinnecock Bay. Multiple lines of evidence linked the regrowth of clam populations to the establishment of spawner sanctuaries including the survival and spawning patterns of stocked clams, physical circulation patterns in the estuary, spatial-temporal patterns of spawning and clam larvae abundances, the timing of increases in recruitment and landings, the disproportionate increase in landings of smaller/younger clams, and differential landings within Shinnecock Bay compared to other regional estuaries. Higher clam densities caused biological filtration times of the bay to decrease from three weeks to 10 days in the western Bay contributing to concentrations of HABs and chlorophyll a significantly decreasing and water clarity and the extent of seagrass beds increasing. Given the trajectories of many bivalve populations across the globe (Jackson et al., 2001; Lotze et al., 2006; Beck et al., 2011; Ermgassen et al., 2012) and the recognition of the need for marine ecosystem restoration (Duarte et al., 2020), the approach taken here may be promising for recovering hard clam and other bivalve populations in estuaries around the globe. Finally, the use of protected sanctuaries in Shinnecock Bay and a long-lived bivalve species (M. mercenaria; 40+ yr) provides hope that effects of restoration presented here will persist over an extended period.