Tidal Flat Depositional Response to Neotectonic Cyclic Uplift and Subsidence (1–2 m) as Superimposed on Latest-Holocene Net Sea Level Rise (1.0 m/ka) in a Large Shallow Mesotidal Wave-Dominated Estuary, Willapa Bay, Washington, USA

The late-Holocene record of tidal flat deposition in the large shallow Willapa Bay estuary (43 km in length), located in the Columbia River Littoral Cell (CRLC) system (160 km length), was investigated with new vibracores (n=30) and gouge cores (n=8), reaching 2–5 m depth subsurface. Reversing up-core trends of muddy sand to peaty mud deposits in marginal tidal flat settings demonstrate episodic submergence events resulting from cyclic tectonic uplift and subsidence (1–2 m) in the Cascadia subduction zone. These short-term reversals are superimposed on longer-term trends of overall sediment coarsening-up, which represent the transgression of higher-energy sandy tidal flats over pre-existing lower-energy tidal flat mud and peaty mud deposits in late-Holocene time. Fining-up trends associated with channel lateral migration and accretionary bank deposition occurred only infrequently in the broad intertidal flats of Willapa Bay. Vibracores and gouge cores were dated by 14C (n=16) and paleo-subsidence event contacts (n=17). Vibracore median probability 14C ages ranged from 0 to 6,992 yr BP and averaged 2,174 yr BP. Dated sample ages and corresponding depths of tidal flat deposits yield net sedimentation rates of 0.9–1.2 m ka-1, depending on the averaging methods used. Net sedimentation rates in the intertidal flat settings (~1.0 m ka-1) are comparable to the rate of net sea level rise (~1.0 m ka-1), as based on dated paleo-tidal marsh deposits in Willapa Bay. Reported modern inputs of river sand (total=1.77x104 m3 yr-1), from the three small rivers that flow into Willapa Bay, fall well short of the estimated increasing accommodation space (1.9x105 m3 yr-1) in the intertidal (MLLW-MHHW) setting (1.9x108 m2 surface area) during the last 3 ka, or 3.0 m of sea level rise. The under-supply of tributary sand permitted the influx of littoral sand (1.1x105 m3 yr-1) into Willapa Bay, as based on the net sedimentation rate (~1.0 m ka-1) and textural composition (average 60 % littoral sand) in analyzed core sections (n=179). The long-term littoral sand sink in Willapa Bay’s intertidal setting (55 % of total estuary area) is estimated to be about 5 % of the Columbia River supply of sand to the CRLC system, and about 30% relative to the littoral sand accumulated in barrier spits and beach plains during late-Holocene time. A 2.0 m rise in future sea level could yield a littoral sand sink of 2.2x108 m3 in the Willapa Bay intertidal setting, resulting in an equivalent shoreline retreat of 600 m along a 50 km distance of the barrier spit and beach plains that are located adjacent to the Willapa Bay tidal inlet. Willapa Bay serves as proxy for potential littoral sand sinks in other shallow mesotidal estuary-barrier-beach systems around the world following future global sea level rise.


Introduction
Modern tidal flat and channel deposits in Willapa Bay (Figure 1) have been the subjects of many studies due, in part, to their large size and accessibilities from bay shorelines and fishery ports (Hill & Chin, 1979;Clifton & Phillips, 1980;Wiberg et al., 2013;Nittrouer et al., 2013).The Willapa Bay study area also provides sequences of tidal flat and channel deposits in uplifted latest-Pleistocene terraces (Clifton, 1983;Gingras et al., 1999;Hubbard, et al., 2002).These late-Pleistocene analog deposits are widely exposed in retreating bay cliffs along the eastern shores of the large estuary.Comparisons between the modern and late-Pleistocene tidal deposits lead to some early conceptual models of deposit lateral trends and vertical sequences, as developed from migrations of tidal channels (Clifton & Phillips, 1980;Smith, 1988;Clifton et al., 1989;Clifton, 1994) and tidal inlets (Smith et al., 1999).However, no vibracoring or 14 C dating of the late-Holocene basin fill had been widely undertaken in the broad intertidal flats of Willapa Bay.As a result, it was not known whether the late-Holocene depositional sequences in Willapa Bay were indeed controlled by lateral channel migrations.In the mid-1980s Willapa Bay was investigated for evidence of tidal marsh subsidence records that could confirm great earthquakes in the central part of the Cascadia subduction zone (Atwater, 1987).As many as seven episodes of abrupt marsh submergence are recorded by their burials by bay mud sequences in Willapa Bay (Atwater, 1997), and in adjacent estuaries (Figure 1) of the central Cascadia margin (Atwater et al., 2003;Peterson & Cruikshank, 2014).The tidal marsh responses to the cyclic interseismic uplift and coseismic subsidence (1-2 m) with recurrence intervals of several hundred years are well documented in Willapa Bay during latest-Holocene time (3.4 to 0.3 ka) (Atwater & Hemphill-Haley, 1997;Atwater et al., 2003).Although exhumed wetlands had been observed in some eroding eastern tidal flat shorelines (Atwater, 1997) no focused studies of tidal flat response to the tectonically-forced sea level changes had been conducted in Willapa Bay.It was not known whether and/or how the intertidal flats might have responded to the episodic coseismic subsidence and interseismic uplift cycles.
Willapa Bay is located immediately north of the Columbia River mouth in the central part of the Columbia River Littoral Cell (CRLC) system (Figure 1).The CRLC has been intensively investigated for ocean shoreline response to changes in sand supply due to concerns about historic impacts from 1) jetty construction at the Columbia and Grays Harbor bay mouths, 2) impoundments of the Columbia River tributaries, and 3) offshore dredge disposal of Columbia River mouth sand deposits (Gelfenbaum & Kaminsky, 2010).Because of its relatively large size (~350 km 2 ) Willapa Bay has been thought to serve as an import littoral sand sink in the CRLC system.Lacking large river tributaries, the modern deposits of Willapa Bay are dominated by littoral sand (Andrews, 1965) supplied by the adjacent Columbia River (Ballard, 1964;Baker et al., 2010;Peterson et al., 2010a).However, long-term rates of littoral sand accumulations in Willapa Bay had not been established.The relative importance of Willapa Bay as a littoral sand sink in the CRLC system during latest-Holocene time was unknown.
In this reconnaissance study we analyzed sediment composition, grain-size distributions, and radiocarbon ages from 31 vibracores (2-5 m depth) in the late-Holocene tidal flats, accretionary banks and tidal inlet shoals in Willapa Bay.These analyses establish records of intertidal flat sediment accumulation under conditions of tectonic cyclic subsidence/uplift as superimposed on net sea level rise during late-Holocene time.Fining-up accretionary bank deposits were found adjacent to modern tidal channels, as expected.However, the marginal intertidal flats are dominated by 1) repeated fining-up (sand to mud) sequences from tectonic cyclic uplift/subsidence events and 2) overall coarsening-up sequences (peaty mud to sand) from shoreward transgression of sand flats over finer-grained estuarine deposits.Influxes of littoral sand into Willapa Bay are associated with both short-term coseismic submergence and long-term transgressive submergence of the broad tidal flats.These relations are of potential importance to evaluating littoral sand sinks in other large shallow mesotidal wave-dominated estuaries and lagoons around the world.Such inshore sand sinks could rapidly develop following predicted future sea level rise from ongoing global warming (DeConto & Pollard, 2016;Mengel et al., 2016).

Depositional Settings
Willapa Bay is a large wave-dominated estuary, which has been largely in-filled by sediment, leading to shallow depths of 55 percent intertidal surface area, between Mean Lower Low Water (MLLW=157.9km 2 ) and mean higher high water (MHHW=347.4km 2 ) (Andrews 1965).Unlike some incised valley estuaries (Dalrymple et al., 1992) there is no deep-water central basin in Willapa Bay (Figure 2).The broad tidal flats in the central reaches of Willapa Bay are dissected by two major tidal channels, including the Willapa River Channel and Nahcotta Channel, ranging from as much as ~30 m axial depth at the tidal inlet in the lower estuarine reaches to as little as ~5 m depth in fluvial-tidal floodplains in the upper estuarine reaches.The tides in Willapa Bay are mixed-semidiurnal, with a mean daily tidal range (MLLW-MHHW) of 3.1 m in the central bay reaches (Andrews 1965).Maximum predicted tidal ranges in the central reaches of Willapa Bay reach 4.0 m (Flater & Pentcheff, 2017).The three small rivers that enter Willapa Bay, including the North, Willapa and Naselle Rivers (Figure 2), yield relatively small annual suspended sediment discharges, totaling 131x10 3 mt yr -1 (Karlin, 1980).A large tidal prism in Willapa Bay results in a tidally-dominated hydrography for the shallow estuary.The estuary hydrographic ratio for Willapa Bay is H r =232 for the MHHW-MMLW tidal discharge relative to the mean normalized (six hour) fluvial discharge (Peterson et al., 1984).Willapa Bay is protected from high-energy ocean waves, H s1/3 = 7-13 m winter storm wave heights (Ruggiero et al., 1997;Tillotsen & Komar, 1997), by the Long Beach Peninsula, a progradational barrier spit (Smith et a., 1999;Meyers et al., 1996;Vanderburgh et al., 2010).
Gravely sand bottom deposits occur in the uppermost tidal reaches of the Willapa and Naselle Rivers (Figure 2), where river sand and mud likely contribute to the channel bank deposits.However, the tidal inlet shoals, subtidal channels, and broad intertidal flats in Willapa Bay are dominated by littoral sand (Andrews, 1965;Luepke & Clifton, 1983) as supplied by sand discharged to the littoral zone from the adjacent Columbia River (Figure 1) (Peterson et al., 2013;Peterson et al., 2016).Onshore winds blow across Willapa Bay from the southwest (winter) and northwest (summer), reaching sustained velocities of ≥8 m s -1 during winter storms (Byrnes & Li, 1998).
The Willapa Bay harbor entrance was not stabilized by harbor mouth jetties, though limited sand dredging was undertaken between the 1940s and 1960s in the tidal inlet and the Willapa River Channel (Hands et al., 2000).Large non-native oysters were first introduced to Willapa Bay from the United States East Coast (1896) and then from Japan (1928) (Clyde Sayce, Pers. Comm., 2003).For the purposes of this article the deepest burial depths of non-native (large) oysters provide a useful early-historic time marker (early 1900s) in intertidal deposits in the central reaches of Willapa Bay.Early historic navigation charts of Willapa Bay show a northward migration of the tidal inlet (~ 2.0 km) since 1887, but no significant lateral migrations (<0.5 km distance) of the Willapa or Nahcotta Channels in the central reaches of Willapa Bay since 1911 (USCGS, 1911).

Modern and Ancient Analogs for Lateral Trends and Vertical Sequences
Modern subtidal channels in the central and lower reaches of Willapa Bay are floored by sand (75-99 % by weight), but intertidal flats transition from sand (80-98% sand) in the lower reaches, to mud (70-95% mud) in the highest (most protected) fluvial tidal reaches (Andrews, 1965).Tidal marshes colonize mud flats at ≥0.5 m MTL in protected tidal creek settings of Willapa Bay (Peterson et al., 2000).Brackish-water tidal marshes 1±0.5 m MTL transition to fresh-water supratidal wetlands 2±0.5 m MTL in protected tidal creek valleys in Willapa Bay.Modern deposits and late-Pleistocene terrace deposits in Willapa Bay (Clifton & Phillips, 1980;Smith, 1988;Clifton et al., 1989) were used to develop vertical sequence models that could result from lateral tidal channel migration (Figure 3A).Such models were developed for different estuarine settings that vary in relative depositional energies from ocean waves, tidal currents, and local-wind waves.However, the overall pattern of channel migration and deposition is one of upwards grain-size fining from subtidal channel sand to intertidal sandy mud or peaty mud deposits.Under assumed conditions of shallow basin filling by tidal channel migrations in Willapa Bay (Clifton and Phillips, 1980), the late-Holocene tidal flat depositional records in the central reaches of Willapa Bay were initially expected by these authors to be dominated by amalgamated fining-up sequences.
Figure 3. Two conceptual models of tidal flat deposition in Willapa Bay Part A. Generalized depositional settings and vertical sequences associated with lateral tidal channel migration, and vertical accretion of sand, mud, and marsh flats in Willapa Bay, as summarized from Clifton & Phillips (1980), Smith (1988), Clifton et al. (1989), andPeterson et al. (2000).Part B. Generalized neotectonic cyclic sea-level curve (dotted saw-tooth line) superimposed on averaged relative sea level rise (dashed line) during latest-Holocene time (~3-0 ka).The averaged sea level curve is taken to be ~1.0 m below the dated tidal marsh deposits from vibracore site W13 (see Figures 6 and 10 below) and the nearby modern tidal marsh (BC in Figure 2).The dated coseismic subsidence events in Willapa Bay are taken from Atwater & Hemphill-Haley (1997).Generalized ranges of cyclic submergence/emergence events of relative sea level change (1.0±0.5 m and >1.5 m) are from Peterson et al. (2000) and Atwater et al. (2003).

Cyclic Coseismic Subsidence/Interseismic Uplift
Willapa Bay and the adjacent estuaries of the Columbia River and Grays Harbor (Figure 1) occur in a tectonically active subduction zone (Cascadia) that results in tectonic cycles of coseismic subsidence and interseismic uplift (Figure 3B).The coseismic subsidence events are associated with great earthquakes (Mw 8.5±0.5) and nearfield paleotsunami runups (10±5 m projected ocean shoreline runups) in the central Cascadia margin (Atwater, 1987;Darienzo & Peterson, 1995;Peterson et al, 2015).The coseismic subsidence events in Willapa Bay are recorded by sequences of wetland deposits that were abruptly buried by tidal flat deposits (Atwater, 1997).Similar wetland submergence/burial/emergence records are known from the historic coseismic subsidence (~ 2 m) that occurred in the tidal flats of Cook Inlet from the 1964 Gulf of Alaska rupture (Mw 9.2) (Ovenshine et al., 1976;Atwater et al., 2001).Macrofossil paleo-tidal indicators have been used to estimate the magnitude of past subsidence and uplift cycles in Willapa Bay (Peterson et al., 2000;Atwater et al., 2003).Estimated coseismic subsidence ranges are ~1.0 m for the smaller subsidence events and >1.5 m for the larger subsidence events in Willapa Bay, during latest Holocene time.Similar ranges have been estimated for the upper reaches of Grays Harbor and the middle reaches of the Columbia River estuary, as based on the macrofossil paleo-tidal indicators in those estuaries (Peterson and Cruikshank, 2014).The latest-Holocene coseismic subsidence events in Willapa Bay are dated to approximately 0.3, 1.1, 1.3, 1.7, 2.5 2.8, and 3.4 ka (Atwater & Hemphill-Haley, 1997;Atwater et al., 2003).The coastal subsidence/uplift cycles are largely elastic, so they did not affect long-term rates of relative sea level rise (~1.0 m ka -1 ) in Willapa Bay during late-Holocene time.It was not known whether or how the tidal flat deposits in Willapa Bay might have responded to such short-term changes in relative sea level.The short-term tectonic cycles of sea level change (1±0.5 and > 1.5 m) are at the scale of changing deposit facies types, as based on tidal elevation ranges, in the central reaches of Willapa Bay (Figure 3A).

Holocene Transgressive Bay Filling and Extension of the Barrier Spit
The pre-Holocene incised valley(s) of Willapa Bay are relatively shallow, reaching maximum depths of about -40 m elevation MTL (Figure 4).Two possible narrow-valleys extend seaward under the modern barrier spit and one wider-valley extends under the modern bay mouth, located seaward of the Willapa River (Wolf et al., 1998;Smith et al., 1999;Vanderburgh et al., 2010).By comparison, the pre-Holocene incised valleys of the nearby Grays Harbor and the Columbia River estuaries (Figure 1) reach depths of -70 m and -120 m elevation, respectively, at their modern bay mouth locations (Peterson & Phipps, 1992;Twichell & Cross, 2001, Baker et al., 2010;Twichell et al., 2010).Though the Holocene fill in Willapa Bay is relatively thin (< 40 m thickness), it does exceed the maximum axial depths of the modern tidal channels (10-20 m water depth) in the central reaches of Willapa Bay (Figure 2).Late-Holocene tidal channels located between the tidal inlet and Long Island should not have been restricted in lateral migration by deeper indurated Pleistocene deposits or bedrock.
Filling of the shallow incised valley(s) of Willapa Bay (Figure 4) was likely promoted, during middle Holocene time (9-5 ka), by littoral sand supplied though barrier shoal channels that existed under what is now the Long Beach barrier spit (Vanderburgh et al., 2010).Modern wave-swept barrier shoals at the north end of the spit (Figure 2) might serve as examples of a proto-barrier-shoal complex that protected the interior tidal flats of Willapa Bay during the early part of late-Holocene time (5-3 ka).The presence of beach retreat scarps along the bayside of the Long Beach Peninsula indicates that the marine transgressive ravinement surface extended to at least the eastern shoreline of the present barrier spit (Vanderburgh et al., 2010).Northward extension and consolidation of the progradational barrier spit throughout late-Holocene time is recorded by the preservation of abandoned foredune ridges and intervening coseismic beach retreat scarps (Meyers et al., 1996;Peterson et al., 2010b).The northward extension of the contiguous spit began at about 5.0 ka at profile PION and reached about 80 percent of its modern length by about 2.8 ka at profile PARK (Figure 4) (Peterson et al., 2010b).

Methods
In this reconnaissance study of Willapa Bay, continuous vibracores (7.5 cm diameter and 2-5 m length) were used to sample intertidal flats and channel accretionary banks in the central reaches of the estuary (Figure 5).The vibracores were measured for total penetration, then sectioned and drained before 1) core splitting, 2) scaled photography, and 3) subsampling for grain size distributions and radiocarbon samples.Vibracore site positions were established by GPS and site elevations were measured in the field relative to timed predicted tide levels (resolution 0.5 m MTL).Lower-intertidal site elevations (above MLLW) were later refined by Lidar surveys collected during a spring low tide in 2002 (NOAA, 2017), yielding site elevations estimated to be 0.1 m in resolution relative the NAVD88 datum.The NAVD88 datum is nearly equal to the local Mean Low Water (MLW) datum or about one meter below the local MTL datum.For the purposes of this study the Lidar elevations relative to the NAVD88 datum were adjusted by 1.0 m (subtraction) to approximate elevations relative to the local MTL datum.Vibracore sites were selected on the bases of 1) representative intertidal depositional settings, 2) relatively even spacing between sites, and 3) two denser sampling localities including sites KI11, W2, W3, and W4, as well as sites W13, W14, W15, and W16 to constrain local depositional variability in the marginal tidal flats (Table 1, Figure 6).A total of 30 vibracore sites were occupied including 18 intertidal flats, seven intertidal accretionary banks, two subtidal accretionary banks, two tidal inlet shoals, and one progradational tidal marsh.One dated vibracore (KI11) is included from a previous study of the Tokeland/Kindred Island spits (Figure 2) (Morton et al., 2000).Some of the vibracores lost sediment upon pull-out, leaving core top gaps.Vibracores that did not contain top gaps were analyzed for sediment length relative to measured total penetration to calculate sediment compaction (% shortening).Measured compaction ranged from 2 to 16 % and averaged 8 % for 18 measured vibracores.Core lengths are not expanded in this article, so basal radiocarbon age samples could be ~10 % deeper than reported.Figure 6.Map of core sites in Willapa Bay Vibracores (W) were collected from intertidal flats, channel accretionary banks, tidal inlet shoals, and one tidal marsh (Table 1).Vibracore position coordinates, surface elevations, and lengths (rounded to 0.1 m) are presented in Table 1.Gouge core paleotsunami sites (T) were collected from modern low tidal marshes located immediately adjacent to the modern tidal flats.A rooted stump was sampled for 14 C analysis from an exhumed remnant ghost forest (G1) in a modern tidal flat.5161980 426920 ITF -0.1 0 Notes: Vibracore sites (W), paleotsunami sites (T), and exhumed Ghost Forest (G), including position (UTM Sector 10 T) in meters (m) northing (UTM_N) and meters easting (UTM_E), modern depositional setting, core top elevation in meters (m) relative to mean tidal level (MTL), approximate core length (m), and compaction or length shortening (%) in representative cores.Core top in T7 is taken from a modern marsh core top below artificial fill (disturbed core top).Depositional settings are as follows: intertidal flat (ITF), intertidal accretionary bank (IAB), subtidal channel accretionary bank (STC), intertidal marsh (ITM), and tidal inlet shoal (TIS).
Quantitative grain size analyses of vibracore deposits from Willapa Bay (Table 1) are based on 30-50 cm long channel-samples (~1x1x50 cm dimensions).Large shell fragments were removed from the ~0.5 m long samples prior to homogenization and analysis for 1) sediment size fraction abundances, sand:mud:gravel by wet sieving at <62 µm, 62-2000 µm, and >2000 µm screen sizes and 2) mean size and standard deviation (± 1 σ µm) of the sand size fraction from Sedigraph TM settling in glycol.Due to bioturbation, paleo-liquefaction, and/or vibracoring compaction, primary bedding was substantially disturbed in some sandy tidal flat deposits.Small-scale sedimentary structures were not analyzed in this reconnaissance study.The average grain size distributions from the bulk samples are used to 1) confirm interpretations of general depositional sequences and 2) quantify the relative importance of littoral sand sinks in the intertidal settings in Willapa Bay.Representative 14 C samples were selected from undisturbed vibracore sections to establish local sedimentation rates and a general sediment level curve in the broad tidal flats of Willapa Bay.Conventional and AMS radiocarbon ages were adjusted for 14 C marine reservoir (shell) and were calibrated with the online version of CALIB 7.10 (Calib7.10, 2017;Stuiver et al., 2017).

Modern Depositional Settings
Vibracore and gouge core top sediments (modern) from Willapa Bay (Figure 6, Table 1) were compared for depositional setting, composition-textural relations and elevation relative to Mean Tidal Level (MTL).The modern deposits in the marginal tidal flat sites in the central reaches of Willapa Bay (Figure 2) ranged from mud to sand at elevations between + 0.3 m and -1.7 m MTL (Figure 7).Traces of gravel-sized shell fragments were observed in some core top samples, but they were not included in the quantitative analyses of siliclastic sediments performed for this study.Mud occurred at +0.3 m MTL in site KI11 behind the wind/wave protection of the small Kindred Island spit at the north end of Willapa Bay.Tidal marsh deposits (≤ 1.5 m MTL) at sites T1-T8 were dominated by muddy peat or peaty mud, at elevations between +1.5 m and +1.2 m MTL.For the purposes of this article, peat is defined as > 50% fibrous/matted/rooted organics, muddy peat is defined as 25-50 % peaty organics, and peaty mud is defined as 5-25 % peaty organics, by visual inspection of cut-core surface area.Intertidal accretionary banks, at elevations between -0.4 m and -1.9 m MTL, ranged from muddy sand to mud.Deeper subtidal accretionary banks in major tidal channels that were sampled in the central reaches of Willapa Bay (sites W21 and W23 at elevations of -3.3 m and -3.1 m MTL) were dominated by sand.Tidal inlet or bay mouth shoal deposits, which were exposed to waves and strong tidal flows in the lowest reaches of the estuary (sites W11 and W7), at elevations of -2.0 m and -0.7 m MTL respectively, were also dominated by sand.The recovered core top data from the central reaches of Willapa Bay are generally consistent with the model of intertidal deposit fining with increasing tidal elevation (Figure 3A).However, vertical textural trends in the vibracore subsurface sections generally do not reflect the expected channel fining-up sequences, as addressed below.for 36 core sites in the central reaches of Willapa Bay.See Figure 2 for core site locations.Gouge core site T7 (Figure 6) is excluded from the plot due to possible core top disturbance from artificial fill burial.Depositional setting and elevation data are from Table 1.Deposit texture types (sand, muddy sand, sandy mud and mud) were verified and/or calibrated by quantitative grain size analyses, as are presented below in Section 4.3.

Vibracore Depositional Sequences
Vertical sequences of sediment composition and texture in the Willapa Bay vibracores (2-5 m depth) are grouped by depositional setting from north to south (Figure 6, Table 1).A group of vibracores W1, W3, and W4, from the northernmost tidal flats in Willapa Bay, contained multiple buried peaty horizons in basal layers (Figure 8).The rooted peat, muddy peat, and peaty mud horizons have sharp upper contacts with mud or sandy mud (Figure 9), indicating episodic coseismic subsidence of 1-2 m.The basal peaty mud in vibracore W3 overlies a consolidated mud, interpreted to be a late-Pleistocene terrace deposit.An adjacent vibracore W2 included a basal mud layer, which was abruptly buried (sharp contact) by a thin gravel lag, some 1-2 centimeters in thickness, which graded up-section into a muddy sand layer.Anomalous thin gravel lag laminae were also present in vibracores W3 and W4.Though a peaty horizon was absent in vibracore W2, the presence of multiple buried peaty layers in adjacent core sites suggests that the abrupt burial of the mud layer in site W2 could have also occurred by coseismic subsidence.All four of the vibracores from the exposed northern tidal flats (W1-W4) grade up-section into sand or muddy sand layers (Figure 8).Those four tidal flat sites are currently exposed to winter winds and wind-waves from the southwest.One nearby vibracore KI11, from a tidal flat site that is protected from wind-wave propagation behind the Kindred Island spit (Figure 2), grades up-section from alternating sand and muddy sand layers into mud.Vibracores W5 and W6 were located at or near modern tidal channel banks in the northern part of Willapa Bay.They both included mud and sand couplets that are diagnostic of accretionary bank deposition (Smith, 1988).Both vibracores W5 and W6 graded up-section into muddy sand deposits.Vibracores W8 and W9, from a more exposed tidal flat setting on the south side of the Willapa River channel (Figure 6), contained alternating sand and muddy sand layers.Vibracore W7, from a site that is exposed to waves and currents from the tidal inlet, was dominated by cross-bedded sand.
Vibracores W13 and W16, in an eastern tidal flat setting in the central reaches of Willapa Bay (Figure 6), contained multiple buried peaty horizons (Figures 9 and 10).At least five coseismic subsidence events were recorded in vibracore W13, nearly 5 m in length.Two intervening vibracores W14 and W15, in the northeastern tidal flat setting, contained multiple abruptly buried mud horizons, of possible coseismic subsidence origins.Vibracores W13, W14, W15, and W16 grade up-section to sand or muddy sand.Those four vibracores contrast with two adjacent vibracores W12 and W17 that yielded short sections (~ 2 m) of buried accretionary bank deposits, i.e., sand and mud couplets.Vibracore W10, from a tidal flat site, which is located near both a major tidal channel and the tidal inlet, consisted of alternating muddy sand and cross-bedded sand layers.Vibracore W11, from a tidal inlet shoal that is exposed to ocean surf and strong tidal currents, was dominated by wave compacted sand.Vibracore W18, from a progradational shoreline on the west side of Willapa Bay, records one coseismic subsidence event below a modern tidal marsh.
Figure 10.Vibracore logs from the middle reaches of Willapa Bay Vibracore logs include 1) sediment compositions, 2) coseismic subsidence (CS) contacts (sharp) between underlying rooted or peaty mud and overlying muddy sand or sandy mud, and 3) radiocarbon (RC) dated samples in thousand years (ka) before present.Possible coseismic subsidence events (?) are shown at sharp contacts between underlying laminated mud and overlying muddy sand, sand or gravel lag (GL).
Vibracores W19 and W21, from subtidal accretionary bank settings in the major Nahcotta tidal channel, were dominated by sand (Figures 10 and 11).Vibracores W25, W26, W28 and W29, from intertidal accretionary bank settings (Figure 6, Table 1), included distinctive sand and mud couplets (Figure 12).Vibracore W27, from a western tidal flat, contained two possible coseismically subsided mud layers above a coseismically subsided peaty horizon at about -3 m MTL.The alternating sandy mud and mud layers in vibracore W27 grade up-section to sand at the deposit surface.Some thick mud layers in vibracore W27 are intruded by small sand dikes and sills, suggesting post-depositional coseismic paleoliquefaction.It is possible that convolute sand beds in vibracore W24 are also from coseismic paleoliquefaction.A short vibracore W30, from the southernmost tidal flat in Willapa Bay, reached one coseismic subsidence event at ~2 m depth subsurface.Mud deposits above the subsided peaty horizon graded up-section to sandy mud in the southernmost tidal flat, at a modern elevation of about -0.5 m MTL.

Deposit Grain Size Distributions
The vibracore sections recovered from Willapa Bay were analyzed for quantitative grain size distributions from bulk channel-split samples (30-50 cm in length) at 0.5 m depth intervals (Table 2).The quantitative grain size data help to verify interpretations of up-section depositional sequences and to constrain estimates of littoral sand accumulation rates in Willapa Bay.One vibracore KI11, taken from Morton et al. (2002), was not analyzed for the means or standard deviations of grain sizes in the sand fractions.Gravel was absent above trace levels (<1 %) in all but a few of the vibracore samples that are presented in Table 2. Sand abundances ranged from 0 to 100 percent, and averaged 70 percent for all analyzed samples (n=159).Sand mean grain size ranged from 86 to 249 µm and averaged 185 µm for all analyzed samples (n=144).A linear least squares analysis between percent sand abundance (sand %) and mean grain size (µm) of the sand fraction yields a correlation value of R 2 =0.7 (n=144).The positive correlation between relative sand abundance and sand mean grain size in the 144 vibracore channel-split sample indicates that sand abundance was likely influenced by energy of deposition and/or winnowing in the shallow depositional settings in Willapa Bay, rather than by relative proximity to sand source(s).The mean grain size (185 µm) of the analyzed sand size fractions, as reported in Table 2, are within the grain size range of littoral sand supplies in the nearby CRLC beaches (sample means 0.19-0.27µm) and inner-shelf (sample means 0.13-0.23 µm) (Peterson et al., 2016).12:88/ na 08:92/ na Notes: Channel-split samples (1 cm width, 1 cm depth and 30-50 cm length) taken at 0.5 m intervals from 0 to 4 m depth below surface.Large shell fragments (>0.5 cm diameter) were removed by hand from the sample prior to grain-size fraction analyses by wet sieving.Grain size fractions sand:mud:gravel are shown as percent weight.Gravel was absent above trace levels (< 1 %) in all but a few samples.The sand size fractions (62-2000 µm) were analyzed for mean and standard deviation (mean ± 1 σ µm).Some samples were not analyzed (na) due to disturbed or missing core sections.One sample from vibracore W29, at 3.5-4.0m depth, is presented separately in these Notes (W29 65:35/195±29 µm).
The grain size data for 20 selected vibracores from Willapa Bay are summarized in plots of sand abundance (sand %).Coarsening-upward trends of sand abundance, relative to silt and clay, are shown in vibracores W2, W3, and W4 (Figure 13) from the exposed northern tidal flat settings (Figure 6).These coarsening-up trends correspond to the burial of peaty mud and sandy mud deposits by winnowed sand deposits, as shown in textural facies core logs (Figure 8).Vibracore KI11, from a tidal flat that is protected from wind-wave propagation behind the Kindred Island spit (Figure 2), shows a fining-up trend of decreasing sand abundance.Vibracore W5, from a modern accretionary bank also records a fining-up trend.Vibracore W6 includes a coarsening-up trend that corresponds to the burial of accretionary bank deposits by muddy sand tidal flat deposits, as interpreted from observed textural facies sequences.The relative abundances (sand %) of sand size fractions (62-2000 µm) are plotted relative to depth (m) subsurface in selected vibracores from Willapa Bay (Figure 6).Silt and clay make up the non-sand components because gravel is generally < 1% abundance by weight.Data are from Table 2. Sand abundance plots in vibracores W13, W14, W15, and W16, from the eastern tidal flat setting in Willapa Bay, include coarsening-up trends of sand abundance (Figures 13 and 14).The coarsening-up trends correspond to up-section transitions from peaty mud and sandy mud deposits to sand deposits in wind-wave winnowed tidal flats on the eastern side of Willapa Bay (Figures 2, 6, and 10).One nearby vibracore (W12) from a modern accretionary bank shows a fining-up trend, as expected.Another nearby vibracore (W17) records a coarsening-up trend that corresponds to a tidal flat burial of underlying accretionary bank deposits.Two vibracores from modern accretionary banks (W25 and W26) demonstrate fining-up sequences.Vibracore W27 includes a coarsening-up trend from a western tidal flat setting in Willapa Bay.Vibracore W29, from a modern accretionary bank bounding the southern tidal flat setting in Willapa Bay (Figure 2), records a complex up-section vertical sequence of fining-up in basal deposits and slight coarsening-up in core top deposits.The reversing up-section grain size trends in vibracore W29 correspond to an apparent burial of accretionary bank deposits (Figure 12) by a recent onset of tidal flat deposition.The relative abundances (sand %) of sand size fractions (62-2000 µm) are plotted relative to depth (m) subsurface in selected vibracores from Willapa Bay (Figure 6).Silt and clay make up the non-sand components, because gravel is generally < 1% abundance by weight.Data are from Table 2.

Deposit Ages
Vibracore samples are dated by 1) radiocarbon analyses, 2) the last Cascadia paleo-subsidence event (AD1700), and 3) the deepest presence of introduced in-situ oyster shells (~ AD1900) as shown in Table 3. Median probabilities of radiocarbon ages from the deepest vibracore samples range from 126 yr BP (vibracore W29) to 6,992 yr BP (vibracore 3).The young age (~1.3 ka) of the basal sample (2.8 m depth subsurface) in vibracore W29 corresponds to rapid sedimentation in a channel accretionary bank (Figure 11).The old age (~7.0 ka) of the basal sample (3.8 m depth subsurface) in vibracore W3 corresponds to a wetland, which pparently developed on a late-Pleistocene bay terrace surface (Figure 8).The anomalous wetland age (7 ka) in vibracore W3 likely predates intertidal deposition at corresponding paleo-tidal elevations (-4.8 m MTL) in Willapa Bay.

Potential Contributions of Paleotsunami Deposition
What role, if any, might nearfield paleotsunami inundations associated with coseismic subsidence events in Willapa Bay (Atwater, 1997) have on tidal flat deposition?For example, it is not known whether several anomalous gravel lag layers (1-5 cm thickness) in vibracores W2, W3, W4, W15 and W24 (Figures 6, 8, and 10) represent deposition by paleotsunami, storm surge, debris flow, or other extraordinary mechanisms in Willapa Bay.Of particular-concern to this study, is the potential importance of nearfield (Cascadia) paleotsunami surges in the reworking/deposition of sand in Willapa Bay's broad intertidal flats.Paleotsunami deposition of resuspended sand in sandy tidal flats could be difficult to discriminate from more common processes of deposit reworking by combined tidal flow and storm wind-waves.Eight tidal marsh sites in the central reaches of Willapa Bay (Figure 6) were selected to serve as proxies for paleotsunami sand deposition that could have occurred in adjacent marginal tidal flats.The eight tidal marsh sites are located within several 10's of meters of the modern tidal flat shorelines, but they are protected from shoreline retreat in small embayments and/or behind narrow prograded sand berms (< 1.0 m height).The proxy tidal marsh sites are located away from narrowing tidal creek valleys and/or large tidal creeks to avoid runup amplification effects that could have been associated with constricted tidal creek morphologies.The marsh sites were gouge cored to establish paleotsunami sand thicknesses associated with the last several coseismic subsidence events in Willapa Bay (Figure 3B) (Atwater et al., 2003).
The sandy paleotsunami deposits in the eight paleotsunami gouge core sites ranged from 0 to 3.0 cm in thickness (Figure 15).These values are based on the total estimated thickness of distinct sandy laminae above each coseismic subsidence contact.The textural composition of the sandy paleotsunami deposits ranged from muddy sand to sandy mud.Non-sandy laminae of silt and detrital organics were not included in the measured thicknesses of the sandy paleotsunami deposits.Two anomalous sandy laminae, which were not associated with coseismic subsidence contacts, occurred in multiple cores taken at site T1.Those sandy laminae occurred at 18 cm and 110 cm depth, respectively, and might correspond to the historic 1964 Gulf of Alaska rupture tsunami (Landers et al., 1993) and a possible northern Cascadia rupture paleotsunami at ~800 yr BP (Clague & Bobrowsky, 1994;Peterson et al., 2015).The paleotsunami sandy deposits measured in the proxy marsh sites (Figure 15) averaged 0.5 cm in thickness above 17 coseismic subsidence contacts, including zero thickness values.The thickest sandy paleotsunami deposits (2-3 cm total thickness) were associated with the last nearfield paleotsunami inundations in tidal marsh localities T1 and T2, which are the closest in distance (10-15 km) to the tidal inlet.All eight of the shallowest marsh subsidence contacts, which are assumed to represent the AD 1700 rupture event (Satake et al., 1996), contained paleotsunami sandy deposits.Only one half of the second subsidence (down-core) event contacts (total n=6) contained sandy paleotsunami deposits, yielding an average thickness of 0.3 cm (n=3).None of the third (down-core) subsidence contacts (n=3) yielded sandy paleotsunami deposits.
Figure 15.Paleotsunami core logs from selected tidal marsh sites in Willapa Bay Subsidence contacts as numbered 1-3 might represent Cascadia rupture events at 0.3, 1.1, and 1.3 ka (Atwater & Hemphill-Haley, 1997), but are not 14 C dated in this reconnaissance study.Paleotsunami sand layer thicknesses (total laminae thickness in centimeters) are shown at corresponding coseismic subsidence contacts, as numbered.
Paleotsunami laminae textures include muddy sand (Ms) and sandy mud (Sm).Two muddy sand laminae that were not associated with coseismically buried peats occur at 18 cm depth and 110 cm depth in tsunami core site T1.

Discrimination between Tidal Flat Deposition from Channel Migration and Bay Flat Transgression
Vibracore deposit vertical trends (2-5 m subsurface depths) from the central reaches of Willapa Bay (Figures 2  and 6) fall into three general groups based on composition/texture descriptions (Figures 8,10,and 11) and quantitative grain size trends (Table 2).Group 1 vibracores are from tidal flats exposed to high-energy wave/current conditions near the tidal inlet at sites W7, W9, W10, and W11, and in subtidal accretionary banks of the major tidal channels at sites W19, W20, W21, and W23, which are dominated by sand.Little to no significant textural variation occurred up-section, as based on 0.5 m sample intervals, in those high-energy depositional settings.The group 2 vibracores from the lower-energy intertidal flats and intertidal accretionary channel banks demonstrate substantial variability in deposit composition and/or sediment grain sizes (Figures 13  and 14).Specifically, the marginal tidal flats are generally characterized by 1) peaty mud or sandy mud basal deposits that coarsen up-section to muddy sand or sand deposits in settings that are exposed to local wind-wave propagation at sites W1, W2, W3, W4, W13, W14, W15, W16, W20, W22, W24, and W27.Group 3 vibracores are from intertidal accretionary banks that accumulated sand and mud deposits.The shallow accretionary bank deposits generally fine up-section (sites W5, W6, W12, W17, W25, W26, W28, and W29), though several vibracores at sites W5, W6, W17, and W25, topped-out in slightly-coarser muddy sand, possibly reflecting transitions to lower-intertidal flat settings.With one exception (site W17) the intertidal accretionary bank deposits were obtained from sites within ~0.5 km of modern tidal channels.
The coarsening-up tidal flat sequences (Figures 8,10,and 11) are interpreted to represent 1) bay ravinement of upper-intertidal marsh and mud flat deposits and 2) subsequent transgression of lower-intertidal sand flats over the buried peaty mud and sandy mud deposits, during latest-Holocene sea level rise.These sequences are somewhat analogous 'in miniature' (Kraft, 1971) to the transgressive ravinement and upper-shoreface deposition in the inner-shelf located offshore of Willapa Bay (Peterson et al., 2016).Multiple vibracore sites in a northern tidal flat (sites W2, W3, and W4) and an eastern tidal flat (W13, W14, W15, and W16) do not contain the same number of abruptly buried peaty or mud horizons.Locally variable marsh development and/or variable erosional ravinement by estuary wind-waves and high tides could account for the apparent lack of stratigraphic continuity.
Recent ravinement of an eastern tidal flat locality near vibracore sites W16 and W17 (Figure 6) has exhumed a previously buried supratidal wetland (Figure 16A).The wind-wave erosion is producing a remnant ghost forest of standing stumps and exposed tree roots in the tidal flat.The outer-margin of a rooted stump at site G1 has a 14 C age of 0-303 cal yr BP (Table 3).Due to the 14 C calibration curve within the last few hundred years, this age range corresponds to the last Cascadia subsidence event at AD1700 (Figure 3B) (Satake et al., 1996).The ghost forest being exhumed at G1 was apparently subsided during the last Cascadia mega-thrust rupture, and then buried by mud, prior to the marsh/shoreline retreat that is now exposing the remnant tree stumps.
Figure 16.Remnant ghost forest exposed in eroding tidal flat in Willapa Bay Part A, Modern bay ravinement of a paleo-tidal marsh with exhumed tree stumps/roots, referred to here as a remnant ghost forest (Atwater, 1997) east of vibracore sites W16 and W17 (Figure 6).Peaty mud deposits are being eroded by wind-wave attack during high tides.Wind-wave rippled sand is locally infilling eroded hollows (white dashed lines) in the exhumed ghost forest substratum.An in-situ rooted stump margin (G1) was sampled for 14 C analysis (Table 3), which yielded a calibrated age range that is consistent with the last Cascadia subsidence event in AD1700.A rip-rap boulder (25 cm diameter) is shown for scale in photo foreground (lower white arrow).Photo view is to the northwest.Part B, Position of the exposed ghost forest (white dashed line) near vibracore sites W16 and W17 (white solid squares).The dated in-situ stump (G1) occurs at an elevation of +0.5 m MTL at the landward-most extent of the modern tidal flat (Table 1).Linear sand bars and runoff runnels extend from the Nahcotta channel bank (left side of image) across the marginal tidal flat (white arrow) to the receding bay cliff shorelines (white line).The shorelines are temporarily stabilized with rip-rap boulders north and south of G1.
Wind-waves and tidal currents are currently transporting sand across the adjacent tidal flat from the Nahcotta tidal channel, as evidenced by wide low-amplitude sand bars and intervening runnels that are subparallel to the Nahcotta Channel (Figure 16B).The sand bars are slightly asymmetric with lee slopes (10-30 cm height) facing east, demonstrating modern sand migration by wind-wave resuspension and flood tidal currents across the tidal flat (~2.5 km in width).The transgressing sand is infilling eroded hollows in the peaty mud substratum of the exhumed ghost forest (Figure 16A).However, eventual burial and preservation of the remnant ghost forest at G1 is uncertain.Widespread sand loss is occurring along many of the eastern bay shorelines, leading to shoreline erosion and retreating bay cliffs (Figure 2), such as at the exhumed ghost forest locality (site G1).
5.2 Late-Holocene Sedimentation Rates in Willapa Bay.
Dated samples from tidal marsh settings, tidal flat settings, and tidal channel settings (Table 3) are plotted relative to sample depth subsurface (m) and age (ka) in Figure 17A.A general sedimentation trend yields a long-term sedimentation rate of ~1 m ka -1 for the combined depositional settings in Willapa Bay.Dated samples from two of the depositional settings, including tidal marsh and tidal flats, are re-plotted relative to mean tidal elevation (m MTL) in Figure 17 B to produce, respectively, a tidal marsh curve and a tidal flat sediment level curve.A sea level curve is projected at 1.0 m below the tidal marsh level curve, based on an assumed low marsh elevation of 1.0 m MTL in Willapa Bay (Peterson et al., 2000).The tidal flat sediment curve occurs about 1.0 m below the estimated sea level curve and about 2.0 m below the tidal marsh deposit curve.terrace (Figure 8).Sample depositional setting and age data are from Table 3.Part B, Sample age (ka) is plotted versus tidal elevation (±m MTL).An averaged modern low-marsh setting (large solid square) and an averaged modern tidal flat setting (large solid circle) are shown at 0 ka, though individual sites in both settings include substantial elevation variability (Table 1).Tidal marsh and tidal flat curves (trends) are approximated.A sea level curve is plotted at 1.0 m below the tidal marsh curve, assuming a low marsh elevation of 1.0 m above MTL in Willapa Bay (Peterson et al., 2000).
Tidal flat sedimentation rates are divided by core interval ages (basal age-modern core top) in ten representative vibracore sites from lower-intertidal settings (Tables 1 and 4).The mean and standard deviation (mean ±1σ) of calculated sedimentation rates (n=10) are 1.2±0.7 m ka -1 .Alternatively, the average core length (2.2 m) divided by the average interval age (2.3 ka) for the ten samples yields a net sedimentation rate of 0.9 m ka -1 .For the purposes of this article we use an estimated average sedimentation rate of 1 m ka -1 for the Willapa Bay lower-intertidal flat settings during latest-Holocene time.1).Two tidal flat sites (W6 and W12) included short accretionary bank sections.Vibracore sites KI11 and W3 are excluded, respectively, due to their positions behind an accreted sand spit and on a possible late-Pleistocene terrace.Sedimentation rates are based on core length (m) divided by interval age (ka).The dated core intervals extend from reported basal 14 C age depths to modern core surfaces (Figures 8,10,and 11).Vibracore 14 C ages (median probability of calibrated ranges) are from Table 3.

Littoral Sand Sinks in Willapa Bay
To estimate littoral sand accumulation rates in Willapa Bay during latest-Holocene time, the intertidal surface area (m 2 ) between MHHW-MLLW is multiplied by the average sediment sand abundance (percent sand fraction) and an average sedimentation rate (m ka -1 ) from analyzed vibracore and gouge core sections (Tables 2, 3, and 4).The reported intertidal surface area (189.5 km 2 or ~1.9x10 8 m 2 ) represents about 55 percent of the Willapa Bay surface area (Andrews, 1965).The average sand abundance in the latest-Holocene intertidal deposits (0-3 m depth subsurface) in Willapa Bay is estimated from 1) the textural composition (percent sand) in 0.5 m length vibracore sections (n=159) (Table 2) and 2) apparent sand abundance (0 percent) in 0.5 m length gouge core (paleotsunami tidal marsh site) sections (n=20) (Figure 15).The average relative sand abundance for the 179 intertidal deposit intervals analyzed is 60 % sand.An average sedimentation rate of 1 m ka -1 is based on dated vibracore sections of representative tidal flat sites (Table 4).The exclusion of subtidal settings, ~ 45 % of the Willapa Bay surface area, likely underestimates littoral sand accumulation rates in Willapa Bay.Rotary drilling will be required to sample and date the deeper subtidal channel deposits in Willapa Bay to establish their role(s) as littoral sand sinks.
The annual volume of tributary sediment supply to Willapa Bay, is based on 1) an annual suspended sediment discharge (131x10 3 mt yr -1 ) (Karlin, 1980), 2) the relative proportion of bedload=25% of suspended sediment discharge, and 3) an assumed mass to deposit volume ratio of 0.54.The tributary sand supply to Willapa Bay is estimated to be 1.77x10 4 m 3 yr -1 (Table 5).The estimated annual tributary sand supply was insufficient to fill the increasing intertidal accommodation space in Willapa Bay (1.9x10 5 m 3 yr -1 ), as based on an assumed sea level rise rate of 1.0 m ka -1 in the intertidal surface area (Figure 17B).The deficit in tributary sediment supply to Willapa Bay under-filled the submerging intertidal area, thereby providing the accommodation space for littoral sand influx into the estuary during late-Holocene time.

Tidal Flat Depositional Response to Tectonic Cyclic Subsidence and Uplift
A total of 30 peaty and mud layers (15-50 cm thickness) include sharp upper contacts with overlying coarser-grained deposits in 13 vibracores from the central reaches of Willapa Bay (Figures 8,10,and 11).The interpreted subsidence contacts occur above peaty-mud or muddy peat-deposits (n=18) and laminated mud or rooted mud deposits (n=12) in the 13 vibracores.The initial burial deposits include sand (n=3), muddy sand (n=13), and sandy mud (n=14), and average 50 cm in thickness (n=30).Though the episodically buried peaty mud and mud deposits are interpreted as coseismic subsidence events, the sandy tidal flat deposits must also have experienced episodic submergence.For example, the relatively rapid infilling of the tidal flat site KI11 (Morton et al., 2002) records several up-section transitions from sand to muddy sand before finning-up to mud at the core top.It is unknown whether the alternations between sand and muddy sand at site KI11 represent coseismic subsidence events or aseismic processes.
The dated intervals between coseismic subsidence events in Willapa Bay during the last 3 ka range from as little as 200 years to as much as 1,100 years (Figure 3B) (Atwater et al., 2003).However, the duration of the first stage of the interseismic strain interval, a period of interplate recoupling, might last only a few decades after coseismic subsidence.Tectonic uplift during the middle stage of the interseismic interval could emerge the tidal flats to pre-subsidence tidal elevations within another 100-200 years (Cruikshank & Peterson, 2017).The third stage of the present interseismic interval in Willapa Bay is characterized by slight subsidence (several millimeters per year), based on a decade of averaged GPS station vertical velocities.Could a sustained period of slight submergence during historic time account for the apparent bay cliff erosion that is occurring presently along some of the eastern shorelines in Willapa Bay (Figures 2 and 16B)?The periods of maximum submergence following cosesimic subsidence in Willapa Bay, are relatively brief-certainly less than some of the shorter rupture recurrence intervals of ~ 200 years between rupture events (Atwater et al., 2003).In summary, the duration of post-subsidence maximum submergence in Willapa Bay probably only lasts about one century.Such short periods of maximum submergence and associated sandy unit deposition (average 50 cm thickness) must occur before the tectonic uplift and sedimentation lead to emergence to higher tidal flat elevations, such as occurred in Cook Inlet after the 1964 Gulf of Alaska rupture (Atwater et al., 2001).The post-subsidence influx of sand (0.6x10 8 m 3 ) into Willapa Bay (Table 5) likely occurred within a short span of time (~100 yr), yielding an average annual littoral sand sink of 6x10 5 m 3 yr -1 for the assumed maximum submergence time interval.

Sources of Post-Subsidence Littoral Sand Supply to Intertidal Flats
Given the assumed short period of post-subsidence maximum submergence in Willapa Bay (~100 yr), what mechanisms could account for the rapid influx of littoral sand into the subsided tidal flats?Paleotsunami resuspension cannot account for the post-subsidence burial events due the negligible thickness of paleotsunami sand (average 0.5 cm thickness above n=17 subsidence contacts) in the proxy marsh sites (Figure 15) relative to the thickness of the initial sandy burial units (average 50 cm thickness for n=30 vibracore sections above interpreted subsidence contacts).The immediate sand supply must have been delivered through the major tidal channels (Figure 18).Combined flow from enhanced wind-wave propagation and tidal currents over the subsided tidal flats could have distributed available sand and mud to the vertically accreting tidal flats.However, substantial channel incisions would have been necessary to sustain the influx of littoral sand in the sandy burial units over the subsided intertidal flats (Figures 8,10,and 11) if the source of tidal flat burial sand was from channel banks.It is not known whether widespread channel incisions occurred in Willapa Bay following coseismic subsidence.An alternative source of littoral sand to the submerged tidal flats, via the major tidal channels, could have been from the tidal inlet.Could a post-subsidence increase in the tidal prism have increased flood tidal current transport of littoral sand into the lower reaches of Willapa Bay and then into the major tidal channels?Post-subsidence retreat of ocean beaches adjacent to the Willapa Bay tidal inlet, as shown in Figure 18, are estimated to have been on the order of 200-400 m by Bruun's Rule (Bruun, 1962;Peterson et al., 2000).How much beach retreat would be required to feed the nearshore sand delivery to the tidal inlet, the major subtidal channels, and the intertidal flats at the century time scale?This question is relevant to potential erosion of barrier beaches located adjacent to inshore sand sinks following predicted increases in sea level rise from ongoing global warming, as discussed below.energy barrier-beach systems that adjoin large shallow mesotidal estuaries.The unique conditions of cyclic submergence and transgressive submergence in Willapa Bay address the likely consequences of a rapid and progressive sea level rise, as predicted for the next century (DeConto & Pollard, 2016;Mengel et al., 2016).Littoral sand influx into submerging estuaries could occur over short time scales, as has apparently occurred at the century time scale in Willapa Bay.Longshore sand sinks could exacerbate local beach retreat from offshore sand loss, and/or alongshore sand displacements following future accelerated sea level rise.Such potential inshore sand sinks should be included, where appropriate, in predictive models of shoreline impacts from future global sea level rise.

Conclusions
Latest-Holocene depositional sequences (2-5 m depth subsurface) in the broad marginal lower-intertidal flats of Willapa Bay are dominated by sandy tidal flat transgression(s) over muddy/peaty deposits rather than by lateral channel migrations.Though fining-up sequences are associated with intertidal accretionary bank deposits near modern tidal channels, the vertical sequences in the broad tidal flats are dominated by coarsening-up textural trends.Combined flows from tidal currents and local wind-wave resuspension are likely to have been responsible for the truncation and/or burial of muddy upper-intertidal flat deposits by sandy lower-intertidal deposits during the last several meters of marine transgression in latest-Holocene time (3-0 ka).Shorter-term cycles of fining-up sequences are superimposed on the longer-term transgressive trends, which demonstrate abrupt burials of peaty mud or mud flat deposits by coarser sandy units.The abrupt burial sequences are interpreted to reflect episodic coseismic subsidence of the intertidal flat settings, following ruptures of the Cascadia subduction zone megathrust.Both the long-term transgressive sequences and short-term cyclic burial sequences demonstrate net influxes of littoral sand into Willapa Bay from the adjacent littoral zone.Willapa Bay serves as an important analog for potential littoral sand supply to inshore sand sinks from adjacent barrier beaches following predicted future sea level rise from ongoing global warming.

Figure 1 .
Figure 1.Map of Willapa Bay and the Columbia River and Grays Harbor estuaries Figure is redrafted from Peterson et al. (2016).

Figure 2 .
Figure 2. Physiography of the Willapa Bay study area

Figure 5 .
Figure 5. Several depositional settings in Willapa Bay

Figure 7 .
Figure 7. Plot of surface deposit composition versus tidal elevation in Willapa Bay Core site compositions and grain size types are shown with elevation (± m) relative to mean tidal level (MTL)for 36 core sites in the central reaches of Willapa Bay.See Figure2for core site locations.Gouge core site T7 (Figure6) is excluded from the plot due to possible core top disturbance from artificial fill burial.Depositional setting and elevation data are from Table1.Deposit texture types (sand, muddy sand, sandy mud and mud) were verified and/or calibrated by quantitative grain size analyses, as are presented below in Section 4.3.

Figure 9 .
Figure 9. Photographs of vibracore sections from sites W1 and W13 Vibracore photographs show coseismic subsidence horizon contacts (arrows) between underlying peaty deposits and overlying sand (Part A), muddy sand/sandy mud (Parts B and C), and mud (Part D).See Figures 8 and 10 for complete core logs.

Figure 11 .
Figure 11.Vibracore logs from the southern reaches of Willapa Bay Vibracore logs include 1) sediment compositions, 2) coseismic subsidence (CS) contacts (sharp) between underlying rooted or peaty mud and overlying muddy sand, or sandy mud, and 3) radiocarbon (RC) dated samples in thousand years (ka) before present.Possible coseismic subsidence events (?) are shown at sharp contacts between underlying laminated mud and overlying muddy sand, sand or gravel lag (GL).

Figure 12 .
Figure 12.Photos of vibracore sections from sites W24, W27, and W29 Sandy layers occur above coseismic subsidence contacts (arrows) in vibracores W24 and W27 (Parts A and B).Small clastic sand sills/dikes intruded into mud layers in Part B could reflect post-depositional event(s) of coseismic paleoliquefaction.Sand and mud couplets occur in tidal channel accretionary bank deposits in vibracore W29 (Part C).See Figure 11 for complete core logs.

Figure 13 .
Figure 13.Down-core plots of sand abundance in Willapa Bay vibracores

Figure 14 .
Figure 14.Down-core plots of sand abundance in Willapa Bay vibracores

Figure 17 .
Figure 17.Plots of dated sample depth and tidal elevation Part A, Sample age (ka) versus depth (m) for paleo-tidal marsh setting (solid square), paleo-tidal flat setting (solid circle), and paleo-channel bank setting (solid hexagon).See corresponding interpreted core logs in Figures 8, 10, and 11.Two basal samples from vibracores KI11 and W3, respectively, are excluded due to recent accretion behind the Kindred Island sand spit (Figure 2) and wetland development on a possible late-Pleistocene

Table 1 .
Vibracore and gouge core sites in Willapa Bay

Table 3 .
Vibracore and ghost forest deposit sample ages

Table 4 .
Sedimentation rates in Willapa Bay lower-intertidal flat vibracores

Table 5 .
Late-Holocene sediment accumulation rates in Willapa Bay