Sensitivity of tidal marshes as recorders of major megathrust earthquakes: constraints from the 25 December 2016 M w 7.6 Chiloé earthquake, Chile

: We present evidence of land ‐ level change resulting from the 2016 M w 7.6 Chiloé earthquake from tidal wetlands along the southern coastline of Isla de Chiloé, Chile, to test criteria for the detection of low ‐ level, < 0.1 m, coseismic land ‐ level change. In order to record coseismic land ‐ level change in tidal wetland sediments, both the creation and preservation thresholds must be exceeded. High ‐ resolution diatom analyses of sediment blocks at two tidal marshes reveal that the 2016 earthquake exceeded the creation threshold and a statistically significant change in diatom assemblage is recorded. In contrast, the preservation threshold was not exceeded and the record of coseismic land ‐ level motion is not preserved at any location visited. After nine months, interseismic and coseismic changes are statistically indistinguishable. The most sensitive part of the tidal wetland is not consistent between research locations, possibly as a result of changes in sedimentation after the earthquake. We compare records of change from great earthquakes in Alaska with the record from the Chiloé earthquake to explore the detection limit. We propose that coastal palaeoseismological records are highly likely to underestimate the frequency of major (M w 7 – 8) earthquakes, with important implications for recurrence intervals and assessment of future seismic hazards.

Palaeoseismological investigations of Holocene coastal sediments provide a means to assess the temporal and spatial variability of different earthquake rupture modes in many subduction zones, including Alaska (e.g. Shennan et al., 1999Shennan et al., , 2018, Cascadia (e.g. Nelson et al., 2018), Chile (e.g. Cisternas et al., 2005), Japan (e.g. Sawai et al., 2004), Indonesia (e.g. Dura et al., 2016) and New Zealand (e.g. Hayward et al., 2015). While palaeoseismological investigations based on coastal sedimentary records infer alongstrike segmentation of seismogenic megathrust interfaces (e.g. Briggs et al., 2014;Ely et al., 2014;Cisternas et al., 2017;Shennan et al., 2016Shennan et al., , 2018, instrumental records of recent earthquakes (e.g. Moreno et al., 2018) and modelling studies (e.g. Herrendörfer et al., 2015) also highlight down-dip segmentation, with more moderate earthquakes potentially occurring on deeper segments. Like their larger counterparts, these magnitude 7-class earthquakes may cause abrupt coastal uplift and/or subsidence, yet their recurrence, role in cycles of strain accumulation and release, and lasting impression on the geological record remain uncertain. The November 2018 M w 7.1 Anchorage, Alaska, earthquake produced vertical deformation in the range of +2 to -5 cm (West et al., 2020), below our current estimates for the detection limit from marshes (Shennan et al. 2016). This earthquake impacted many sectors of society (West et al.,2020) and illustrates the importance of knowing more about the recurrence of damaging earthquakes in the magnitude 7-9 range.
Changes in bio-and/or lithostratigraphy within tidal marshes, which may result from seismic or non-seismic drivers, provide records of emergence and/or submergence due to relative sea-level (RSL) changes (Shennan et al., 2016). In regions of coseismic subsidence freshwater wetlands drop instantaneously into the intertidal zone and become rapidly buried by muddy intertidal sediment. This creates a peat-mud couplet with peat overlain by fine-grained muddy sediment and a sharp stratigraphic boundary between the two. Sediment sequences with multiple peat-mud couplets suggest multiple earthquake cycles. Well-developed methodologies allow identification of the largest earthquakes (M~9), which produce metre-scale vertical displacement of the crust across hundreds of kilometres (Nelson et al., 1996;Shennan et al., 2016). These studies led to six criteria for establishing unequivocal coseismic submergence or emergence of 0.1-0.2 m or greater: a) lateral extent of peat-mud or mud-peat couplets with sharp contact, b) suddenness of submergence or emergence, replicated at multiple locations within a site, c) amount of vertical motion replicated at multiple locations within a site, d) synchroneity of submergence or emergence based on age modelling, e) spatial patterns of submergence or emergence, and f) possible additional evidence, such as tsunami inundation or liquefaction concurrent with submergence or emergence (Shennan et al., 2016).
We know little about the lower detection limit of vertical displacement in tidal marsh sequences. This places a critical, but currently poorly known constraint on the lowest magnitude of earthquake that coastal sediments record, and similarly on the crustal deformation at the spatial boundaries of larger ruptures, as the vertical displacement decreases to zero. Understanding the range of earthquake sizes and the spatial extent of ruptures that can be identified provides an important contribution to integrating palaeoseismic records into seismic hazard assessments. Furthermore, identifying and characterising smaller segment ruptures is critical for assessing their significance in releasing strain and evaluating the hazards they pose.
Detection limits of tidal marsh records are dependent on tidal range, stratigraphy and transfer function error. Recent work on late Holocene earthquakes in Alaska suggests a vertical deformation detection limit as low as 0.1-0.2 m, and this range is partly dependent on tidal range (Shennan et al., 2014(Shennan et al., , 2016(Shennan et al., , 2018. In Alaska, Shennan et al. (2016) suggest that a detection limit of 10-15% of Great Diurnal Range or better can be achieved. The land-level change required to exceed the detection limit is dependent on tidal range, as that is usually a control on the elevational distribution of microfossils within tidal marsh sequences. Detection limits of tidal wetland sediment sequences that yield recognisable palaeoseismic evidence therefore depend upon two critical thresholds that record relative sea-level change: creation thresholds and preservation thresholds (Nelson et al., 2006;McCalpin and Nelson, 2009). To exceed creation thresholds, palaeoseismic indicators must be distinct from similar indicators that may result from non-seismic processes (Shennan et al., 2016). Preservation thresholds are exceeded when indicators withstand subsequent alteration and therefore remain part of the palaeoseismic record. Importantly, palaeoseismic records from tidal marsh sequences therefore include only those earthquakes which exceed both the creation and preservation thresholds (McCalpin and Nelson, 2009;. In addition, the error terms of typical transfer function models are~10-15% of the tidal range (Barlow et al., 2013). For locations with low-level land-level changes, the error terms of reconstructed elevations can be larger than the modelled subsidence. Assessment of the lower limit for detecting pre-instrumental earthquakes using established palaeoseismological methods and the potential for geological records to underestimate the frequency of major earthquakes is critical. The 25 December 2016 Chiloé magnitude M w 7.6 earthquake offers a valuable opportunity to assess creation and preservation thresholds within tidal marsh records for an earthquake of known magnitude and geodetically measured surface deformation.
The 25 December 2016 Chiloé M w 7.6 earthquake On 25 December 2016, a M w 7.6 subduction zone earthquake struck south-central Chile. The earthquake occurred within the rupture zone of the AD 1960 Valdivia M w 9.5 earthquake, the largest earthquake ever recorded during the instrumental era (Cifuentes and Silver, 1989; Fig. 1). Palaeoseismic studies provide evidence for eight earthquakes over the millennium preceding 1960, including historically documented earthquakes in 1575, 1737 and 1837 (Cisternas et al., 2005, 2017, Garrett et al., 2015Lomnitz, 2004).
Several studies have explored the characteristics of the 2016 Chiloé earthquake, based predominantly on seismic, interferometric synthetic aperture radar (InSAR) and global positioning system (GPS) datasets (e.g. Lange et al., 2017;Moreno et al., 2018). Initially, it was posited that the maximum coseismic slip, up to 5 m, exceeded the post-1960 slip deficit (Lange et al., 2017;Melgar et al., 2017;Ruiz et al., 2017;Xu, 2017), although this was subsequently revised to~2.9 m (Moreno et al., 2018). Modelled land-level change along the southern coastline of Isla de Chiloé suggests uplift in western locations of <0.5 m (Moreno et al., 2018;Fig. 1). Maximum coseismic subsidence for the 2016 Chiloé earthquake is modelled as <0.1 m to <0.15 m (Garrett et al., 2019;Moreno et al., 2018) in eastern locations along the southern coastline of Isla de Chiloé (Fig. 1). Informal interviews provide no evidence for unusual wave activity following the earthquake, but the tide gauge at Castro, Isla de Chiloé, records a 44 cm oscillation on the falling limb (NGDC/WDS Global Tsunami Database, 2020). Garrett et al. (2019) provided the first field geological evidence of land-level change based on the vertical extent of mortality of coralline algae, recording uplift of 25.8 ± 14.3 cm at Isla Quilan, based on 103 measurements across seven sites (Fig. 1). In this paper, we aim to further develop new field constraints of land-level change additional to those produced by Garrett et al. (2019) and therefore determine the detection limits for reconstructing land-level change in eastern locations.
Tidal wetlands within the zone of subsidence provide an opportunity to assess vertical land-level change and to test methods for determining the detection limit of submergence resulting from coseismic subsidence. Therefore, our aims are, first, to determine whether millimetre-to centimetre-scale coastal subsidence, as occurs in M < 8 megathrust earthquakes, creates a record in tidal wetlands, and second, to identify whether the record is likely to be preserved in the marsh stratigraphy.

Location
We focus on three tidal marsh sites on the southern coastline of Isla de Chiloé across the transition from uplift to subsidence -Quilanlar, Asasao and Ayentema (Fig. 1). Southern Isla de Chiloé lies within the coseismic subsidence zone of the 1960 Valdivia earthquake, with 2.1 + /-0.2 m of subsidence recorded at Quilanlar (Plafker and Savage, 1970). Along the southern coastline, there are extensive ghost forests in tidal marsh and beach locations at Inio (Fig. 1c), Quilanlar (Fig. 1d), Asasao (Fig. 1e) and Ayentema (Fig. 1f), which we conclude result from rapid RSL rise caused by coseismic subsidence in 1960. In comparison, during the 2016 earthquake vertical coseismic deformation occurred on a centimetre to decimetre scale along the southern coastline of Isla de Chiloé (Lange et al., 2017;Melgar et al., 2017;Ruiz et al., 2017;Xu, 2017). Our study sites lie to the east of the sites investigated by Garrett et al. (2019), in a transect across the region of modelled coseismic subsidence. Based on deformation estimates derived from Moreno et al. (2018), we anticipate that Quilanlar subsided by 19 mm, Asasao subsided by 43 mm and Ayentema subsided by 72 mm (Fig. 1).

Methods
To establish signatures of the 2016 earthquake, we completed two field surveys, seven and nine months after the earthquake occurred. The surveys occurred in winter, with the highest tides. We collected multiple blocks of the uppermost 50 mm of sediment at three estuaries, with sampling locations ranging from 0.3 to 1.5 m above mean sea level. Sediment blocks were collected from Asasao and Ayentema in August 2017 and from Quilanlar and Ayentema in October 2017. We collected sediment blocks during both field surveys at Ayentema to establish temporal constraints on creation and preservation potential. Tidal observations were recorded at each location using a portable ultrasound tide gauge (Wesson et al., 2015), with tidal levels derived through comparison with the TPXO 8 Atlas tidal model (Egbert and Erofeeva, 2002). Fig. S1 in the supporting information compares tide gauge data and tidal model predictions for each site. We recorded evidence for shaking and/or liquefaction at each marsh location using GPS-located photography.
We photographed sediment blocks in the field and laboratory and scanned individual sediment blocks using X-ray computed tomography to allow the identification of changes in sedimentation. Example imagery is presented for two sediment blocks (Fig. 2). We counted diatom assemblages for five samples in each block from the uppermost 20 mm of sediment: at 0-1 mm from the surface, 1-5 mm, 5-10 mm, 10-15 mm and 15-20 mm. Typically, 5 mm sediment samples are used in palaeoseismological research. Depending on sedimentation rates this may include some effects of post-seismic elevation change. We adopt the strategy outlined whereby 1 mm minimises the effect of post-seismic changes and allows comparison with the typical sampling interval of fossil sequences. Diatom samples were prepared using the standard techniques (Palmer and Abbott, 1986) and a minimum of 250 diatoms were counted per sample. Fig. 2 includes a summary diatom figure for two sediment blocks (see Results section), with full diatom assemblages presented in the Supplementary Information (Appendix S2 and S3). Diatoms are chosen due to the availability of a large training set (Hocking et al., 2017) and known preservation potential in marsh sequences in Chile (Garrett et al., 2015).
To reconstruct marsh surface elevations from diatom assemblages we follow an approach outlined by Kucera et al. (2005). We take six quantitative methods frequently adopted in sea-level reconstructions (Kemp and Telford, 2015), rather than depending on one single method. We use the software program C2 (Juggins, 2004) to run the following transfer function models: maximum likelihood (ML), weighted averaging-partial least squared (WA-PLS), locally weighted-weighted averaging (LW-WA) (classical and inversed), modern analogue technique (MAT) and weighted modern analogue (WMA). These transfer function methods are mathematically distinct and model, either implicitly or explicitly, species responses in different ways, and it is useful to compare reconstructions from a range of techniques. If downcore biological changes are primarily driven by changes in the reconstructed environmental variable, in this case, elevation, we would expect reconstructions from different models to follow similar trajectories (Juggins and Birks, 2012). We use the 20th percentile of the dissimilarity coefficients calculated between all modern samples as the cut-off between 'close' and 'poor' modern analogues for fossil samples.
Within the transfer function models, we convert sediment block elevation to standardised water level index (SWLI) values to account for differences in tidal range between locations in the fossil and modern datasets (Hamilton and Shennan, 2005). A SWLI value of 100 equates to mean sea level and 200 represents mean higher high water. The value of one SWLI unit for each site is dependent on tidal range: one SWLI unit at Ayentema is 12.4 mm, at Asasao is 11 mm and at Quilanlar is 10.2 mm. Modern datasets are based on Hocking et al. (2017) with updates to nomenclature. We first used the fossil data for each counted sample (0-1, 1-5, 5-10, 10-15, 15-20 mm), and then also investigated merging the uppermost two samples (i.e. into one 5 mm sample) to replicate a commonly employed sampling strategy in palaeoseismic studies.
With the error term of typical transfer function models 10-15% of the tidal range (Barlow et al., 2013), the error terms of reconstructed elevations are larger than the modelled subsidence at each of our field sites. At 24 to 32 SWLI units, the 1σ uncertainties of the transfer function models are approximately an order of magnitude larger than the modelled subsidence in southern Isla de Chiloé in 2016. One part of our approach, therefore, is to test the hypothesis that the change between the samples immediately above and below the earthquake horizon is consistent in sign across the six transfer function-model estimates, rather than a random pattern of submergence and emergence, and different to the changes between adjacent samples in the sediment below the earthquake horizon.

Ayentema (FS1 and FS2)
Ayentema is the easternmost location studied along the southern coastline of Isla de Chiloé ( Fig. 1), with the highest modelled subsidence (72 mm, see Fig. 1; Moreno et al., 2018). We collected sediment blocks along the same transect from the tidal flat to the high marsh in both Field Survey 1 (August 2017, FS1) and Field Survey 2 (October 2017, FS2). Sediment blocks range from blue-grey silty sand on the tidal flat to rootlet-rich herbaceous peat on the high marsh. We found no laterally extensive peat-mud or mud-peat couplet with a sharp boundary during field investigations, within laboratory photographic analyses or X-ray imagery of sediment blocks from either FS1 or FS2. Examples of laboratory photographs and X-ray imagery for two sediment blocks from FS1, collected seven months after the earthquake, are presented here to demonstrate this lack of signal within sediment lithologies (see Fig. 2).
FS1 diatom assemblages at Ayentema are dominated by epipsammic and epiphytic taxa (Fig. 2). Assemblages demonstrate a range of trajectories (positive change, negative change and change close to zero) within the uppermost analysed samples (0-1 mm and 1-5 mm), which are not consistent across all sampled locations on the marsh transect. We present two example assemblages, from sediment blocks AY21 (tidal flat) and AY6 (high marsh; Fig. 2) to demonstrate the changes recorded. Full diatom data for each sediment block are available in the Supplementary Information. Transfer function outputs from WA-PLS, LW-WA, ML, MAT and WMA model runs all produce a modelled decrease in SWLI values between the uppermost samples at AY21 (Fig. 2). The sign of change is therefore consistent across all models used, indicating submergence at AY21. Five of the six models also demonstrate a decrease in SWLI value between the lowermost two samples, although this is unlikely to relate to the 2016 earthquake. At AY6, there is also a notable visible change in diatom assemblage between the 1-5 mm and 0-1 mm samples (Fig. 2). In contrast to AY21, reconstructed elevation change between In FS1 at Ayentema, two of seven locations across the site record submergence between the top two samples (0-1 and 1-5 mm in each box) using all six transfer function models. The other five locations record a mixed record of submergence or emergence that is model-dependent. At the site scale (seven locations across the transect), the change between the top two samples differs from the changes between the other adjacent samples (1-5, 5-10, 10-15, 15-20 mm in each box), summarised by the distribution plot and median change (Fig. 3). The median change between the top samples (0-1 mm and 1-5 mm) is -0.66 SWLI units, compared with +3.81 SWLI units between the lower samples. A two-sample t test assuming unequal variances indicates a statistically significant difference between the two distributions (Table 1). Dissimilarity analysis (>20th and >10th percentiles) shows that many fossil samples at Ayentema from FS1 have poor modern analogues (filled sectors in pie charts in Fig. 4). In FS1, two sediment blocks have good modern analogues across the contact, two sediment blocks show 50% analogy, and three sediment blocks show poor modern analogues (>20% percentile threshold).
Samples collected nine months after the earthquake at Ayentema are dominated by epipsammic and epiphytic taxa (see Supplementary Information). In contrast to FS1, the visible changes in percentages of individual diatom taxa are no longer evident in FS2. Transfer function outputs demonstrate a range of modelled land-level motion across the site with no clear difference between changes at the top and those lower down ( Fig. 3 and Table 1). In FS2, two sediment blocks have good modern analogues across the contact, one sediment block shows 50% analogy, and four sediment blocks have poor modern analogues (>20th percentile threshold, Fig. 4). In addition to sedimentological and diatom evidence, we assessed the marsh at Ayentema for evidence of ground shaking and liquefaction. During the two field surveys, we observed a range of evidence, including detached and toppled sediment blocks and exposure and detachment of tree roots (Fig. 5). Submergence of the marsh would lead to erosion through both increased water depth and greater lateral incursion of individual tides onto the marsh. These additional lines of evidence may therefore provide support for marsh submergence at the location but no quantitative measure of land-level change (cf. Nelson et al., 1996).
Alternatively, erosion may have occurred as the result of a tsunami post-earthquake. There are no observations of a tsunami at the location, but a tide gauge at Castro, eastern Isla de Chiloé, measures wave amplitudes of up to 44 cm (NGDC/ WDS Global Tsunami Database, 2020) on 25 December 2016. Similar observations of the impact of larger waves have been noted in sheltered inlets in Canada following a M w 7.8 earthquake and tsunami in 2012 (see Bednarski, 2014, 2015). There is no sedimentological evidence of a tsunami at Ayentema.

Asasao (FS1)
Asasao is the central research location along our transect, with modelled deformation of 43 mm (Moreno et al., 2018; Fig. 1). We analysed three sediment blocks collected seven months after  the earthquake from the lower (AS9) and high (AS7, AS5) marsh. All three sediment blocks comprise brown silty herbaceous peat with abundant herbaceous rootlets. There is no evidence for a peat-mud or mud-peat couplet with sharp contact across the site. Laboratory photographic and X-ray imagery also fail to identify a signal within the sediment lithology. Diatom assemblages are characterised by epipelic and epipsammic taxa. Upon visual inspection, there are a number of changes to the proportions of individual taxa in the uppermost diatom samples. However, this change is not recorded at a consistent depth, ranging from 0-5 mm to 0-10 mm across the individual sediment blocks.
Transfer functions produce a range of modelled land-level changes at Asasaotwo sediment blocks demonstrate a mix of emergence, submergence and close to zero change, depending on the model (AS9, AS5), with the third mid-elevation sediment block (AS7) showing submergence in all model outputs. Dissimilarity analysis shows that two sediment blocks at Asasao have good modern analogues across the contact, and one sediment block has a 50% analogy between modern and fossil  Fig. 4). If we use the >10th percentile threshold, all samples have poor modern analogues. There is no consistent signal of emergence or submergence which can be traced between samples across the site.
In a similar fashion to Ayentema FS1, at the site scale, despite variable land-level change reconstructions across the marsh, the median change between the top two samples (0-1 and 1-5 mm) at Asasao is -1.53 SWLI units compared with +6.43 SWLI units in the lower samples and a t test indicates a statistically significant difference in the distributions (Table 1). This suggests that the uppermost samples demonstrate a reversal of the sign of relative sea-level change at Asasao and represent a submergence signal when high-resolution sediment analyses are conducted.
Evidence for ground shaking at Asasao was limited to toppled trees close to the tidal marsh. We observed no other lines of evidence to suggest land-level change.

Quilanlar (FS2)
Quilanlar is the westernmost site of our transect (Fig. 1). Analysis of site stratigraphy produces no evidence of change to sediment lithology in the uppermost sediments at Quilanlar. We analysed five sediment blocks for diatoms from the low marsh to the high marsh, which were collected during FS2. Sediment blocks range from blue-grey silty sand (lower marsh) to dark brown sandy peat with abundant rootlets (high marsh). Diatom assemblages include epipsammic, epipelic and epiphytic taxa.
Transfer function outputs demonstrate a range of signals modelled from the uppermost diatom samples. The lowest elevation sediment block (QL11) demonstrates subsidence in all model outputs, with all other sediment blocks showing a combination of modelled emergence, submergence or zero change. Dissimilarity analysis shows that two sediment blocks have good modern analogues across the contact, two sediment blocks have 50% analogy, and one sediment block has poor modern analogues (>20th percentile threshold, Fig. 4). When we use the >10th percentile threshold, all samples have poor modern analogues. There is therefore no consistent signal of land-level change within the model outputs from Quilanlar.
In contrast, the medians show a reversal in sign of change, from positive to negative, at the top (Table 1); the uppermost diatom samples show a median change of -0.53 SWLI units and the lower samples demonstrate a median change of +2.34 SWLI units. Despite this change in sign, the t test result does not indicate a statistically significant difference between the distributions.
In addition to the diatom evidence, we observed cracking within the beach sands which could be traced over a distance of 680 m (Fig. 5). In contrast to Ayentema and Asasao, we observed no evidence of toppled trees or detached sediment blocks at Quilanlar.

Evidence of land-level change from the 25 December 2016 Chiloé earthquake
Model results predict Quilanlar subsided by 19 mm (Fig. 1) and our diatom-based estimates, nine months after the earthquake, suggest a median change of -8 mm between the upper two diatom samples. Model results predict that Asasao subsided by 43 mm (Moreno et al., 2018) and our diatom evidence suggests a median change, seven months after the earthquake, of -19 mm between the upper two diatom samples. At Ayentema, subsidence of 72 mm is predicted (Moreno et al., 2018); in contrast, diatom-based estimates from samples collected seven months after the earthquake at Ayentema suggest median change of -11 mm between the upper two samples. The standard deviation, a measure of the variation in the different transfer function models and the multiple locations within each estuary, of diatom-based estimates of change between uppermost samples are 190 mm at Quilanlar, 45 mm at Asasao and 137 mm at Ayentema (FS1).
In addition to evidence of vertical deformation, evidence of ground shaking resulting from the 2016 earthquake is found along the southern coastline of Isla de Chiloé (Fig. 5). Residents reported minor damage to infrastructure, including wooden buildings (Fig. 5b), alongside extensive cracks in supratidal beach environments (Fig. 5c, d). Tree roots exposed by marsh front erosion (Fig. 5e), toppled sediment blocks (Fig. 5f, g), erosion of marsh surfaces (Fig. 5h) and deposition of mud clasts on marshes (Fig. 5i) may also result from tidal marsh submergence during the earthquake. This additional evidence may suggest coseismic subsidence in some locations, but does not provide a quantitative estimate of land-level motion. Alternatively, these data may provide evidence for a tsunami post-earthquake, with decimetre-scale oscillations noted in tide gauge data (NGDC/WDS Global Tsunami Database, 2020).

Detection limits of tidal wetlands to record coseismic land-level motion
The variable diatom-based records of coseismic submergence described above provide valuable constraints on the limits of vertical surface deformation that we may expect to observe in palaeoseismic records from coastal sediment sequences at varied plate-boundary locations, not just south-central Chile. In all such settings, both creation thresholds and preservation thresholds, and hence the detection limit, vary across tidal wetlands as the sediments accumulating become more or less sensitive recorders of elevation with respect to sea level. The subtidal zone represents the end member at the seaward end. With increasing water depth in the subtidal zone, sediment lithology and biostratigraphy (e.g. diatom species) vary little with changes in sediment surface elevation and the detection threshold increases rapidly. We soon reach the point where coseismic deformation of even 2 m would not be recorded. From the level of the lowest tides, across the tidal mud flat, to the lower limit of pioneer marsh vegetation, microfossil species vary with elevation although tidal currents and sediment movement mean that microfossil-based transfer function models, used to quantify elevation, have relatively large uncertainty terms. Within vegetated marsh environments, from the pioneer marsh to the upper limit of tides, sediment lithology and microfossil assemblages provide increasingly more precise model estimates of elevation. At some point above the upper limit of tides we reach the landward end member, a point at which freshwater terrestrial environments are independent of sea level and would only record coseismic deformation if subsidence was sufficient to lower the site into a marsh or tidal flat elevation.
For submergence >1 m the uncertainty term of the transfer function models (usually given as the root mean square error of prediction, RMSEP), used to reconstruct the amount of submergence, is smaller than the reconstructed submergence, with error bars for samples either side of the stratigraphic contact not overlapping, even when multiplying the RMSEP by 1.96 to estimate 95% confidence limits, and allowing for large tidal ranges, such as Upper Cook Inlet, Alaska (Shennan and Hamilton, 2006). For progressively smaller amounts of submergence, even for areas with smaller tidal ranges, the error bars start to overlap, as we see for our sites in Chiloé. The absolute values for model uncertainty terms, and hence the detection threshold, is partly dependent on tidal range, arising from the standardisation required in using SWLI in the transfer function models. Recent work on late Holocene earthquakes in Alaska suggests a vertical deformation detection limit as low as 0.1-0.2 m equivalent tõ 10-15% of the tidal range,~20-30 SWLI units (Shennan et al., 2014(Shennan et al., , 2016(Shennan et al., , 2018. The estimates of coseismic subsidence for our sites in Chiloé are mostly much smaller, many by an order of magnitude ( Fig. 3 and Table 1).
With submergence >1 m we see variations in the estimated amount of submergence between different sampling points for the same event across a marsh, even while using just a single transfer function model (Shennan and Hamilton, 2006). We still, however, expect elevation change across the contact greater than changes between adjacent samples within the sedimentary units below the contact, i.e. during the interseismic phase. We also expect, for multiple locations for the same event at the site, quantitative estimates of elevation change across the contact all indicating the same sign (emergence or submergence) rather than a random distribution either side of zero (Shennan et al., 2014(Shennan et al., , 2016. Here we compare these expectations using reconstructions from Alaska of variable, but larger, amounts of coseismic subsidence with those from Chiloé (Fig. 6).
Marshes at Girdwood, Alaska, underwent~1.5 to 2.4 m subsidence during the AD 1964 M 9.2 earthquake (Plafker, 1969), the range reflecting variations in sediment compaction of unconsolidated marsh sediments in addition to~1.5 m of regional crustal subsidence. Fig. 6 illustrates the change across a peat-silt contact that marks submergence in 1964, from freshwater marsh to tidal flat, and the change between the three samples at 1 cm intervals below the contact, recorded at three sites across the marsh at Girdwood, and estimated using the same six different transfer function approaches described above. The frequency distributions (coseismic vs. pre-earthquake) are very different, with submergence across the contact much greater than the small amount of submergence indicated by the samples below the contact. This is not surprising, given the large amount of coseismic submergence, greater than the error terms of the transfer function models (e.g. Shennan and Hamilton 2006).
The second example from Kalsin Bay on Kodiak Island had a smaller amount of submergence,~0.4 m in AD 1788, recorded by a peat-silt couplet along with evidence of a tsunami (Shennan et al., 2014), sampled at six locations estimated using the same six transfer function approaches. With smaller amounts of submergence, the frequency distributions overlap more, though they are still clearly different.
These two examples illustrate a trend that we would expect to see continue as the amount of land-level change decreases and the detection limit is reached: the coseismic and nonseismic land-level changes become indistinguishable. Fig. 6 shows the changes recorded from Ayentema FS1, Ayentema FS2, Asasao FS1 and Quilanlar FS2, alongside those from Alaska. As would be expected, given the small amount of submergence at the Chilean sites, there is much overlap between the distributions for change at the top of the sediment blocks and those below. The samples collected seven months after the earthquake at both Ayentema and Asasao, however, show a small difference in the frequency distributions (coseismic vs. pre-earthquake) and a difference in their medians. The samples collected nine months after the earthquake show virtually no difference in either. From these data we infer that multiple samples across a site reveal changes in diatom assemblages that support small-scale subsidence, but the signal is temporary and sedimentation processes between our two field surveys mean that the signal is not preserved.
Combining the data from multiple locations across each marsh, as in Fig. 6 and Table 1, may obscure a signal that is only recorded in the part of the marsh most sensitive to recording change. Fig. 4 shows the records of change at the top of each sediment block, plotted against elevation and sedimentary environment, as well as the percentage of poor (>20th percentile) modern analogues for the top two samples and for all samples at each sampling location. From this evidence we conclude that there is no preferred elevation zone that is more likely to record the small amount of submergence. We also note the preponderance of poor modern analogues for fossil samples analysed (Fig. 6). This is consistent with previous studies in Chile, where many sites demonstrate poor modern analogues for fossil records (e.g. Garrettet al., 2013). The lack of modern analogues across our study sites in southern Isla de Chiloé remains a limitation of this study. Despite this, the reconstructions are consistent with the observed distributions of key taxa in Chilean marshes, even if differences in relative abundances between modern and fossil samples contribute to dissimilarity coefficients exceeding chosen thresholds. We note the need for further modern samples from a wider range of locations in Chile and the importance of this endeavour for the further development of the existing modern training set. These observations, along with the conclusion above regarding no preserved record of submergence, have broader implications with respect to coastal palaeoseismology in Chile.

Implications for palaeoseismological research in tidal wetland environments
Our evidence from the southern coast of Isla de Chiloé contributes to the assessment of the lower limit of detection of pre-instrumental earthquake activity, building on established palaeoseismological methods and assessment criteria (Nelson et al., 1996;Shennan et al., 2014Shennan et al., , 2016Shennan et al., , 2018. Not all of the criteria (d, e and f described in Context and aims) are relevant to contemporary observations from a recent earthquake, but we can assess the remaining three, which were found useful in identifying 0.1-0.2 m of land-level change in Alaska (Shennan et al., 2014(Shennan et al., , 2016(Shennan et al., , 2018: a) laterally extensive peat-mud or mud-peat couplet with a sharp contact at each outcrop or core, b) for multiple locations along the couplet, sudden elevation change determined by assessing the change across the contact in comparison with the changes below and above the contact, c) for multiple locations along the couplet, quantitative estimates of elevation change across the contact all indicating the same sign (emergence or submergence) rather than a random distribution of mean values either side of zero. Our results demonstrate that records of low-level landlevel change can be created, for a limited time period in limited parts of the tidal wetland, following a major earthquake, but even when preserved, these changes would not be detected based on criteria a) to c) due to the lack of a peat-mud or mud-peat couplet. These changes are only detected using high-resolution sampling techniques (sample thicknesses of 1 to 5 mm) and are not shown as a laterally extensive peat-mud or mud-peat couplet, but instead are detected through analyses of the uppermost sediments from the tidal wetland. We also identify that different sections of the marsh are more sensitive to change in different locations, possibly as a result of increased sedimentation following the earthquake. Assessment of lithological changes over large parts of a marsh does not provide a key identifier of coseismic low-level, <0.1 m, land-level motion. Despite this, for palaeoseismological studies, it would be difficult to establish possible earthquakes within the sediment profile without this change in sediment lithology and a sedimentary couplet with a sharp boundary therefore remains an important criterion for establishing low-level land-level change (Shennan et al., 2014(Shennan et al., , 2016(Shennan et al., , 2018. In addition, tsunami deposits may also provide a useful guide for diatom sampling strategies in instances of low-level coseismic deformation, although we identify no such deposits in the sediment blocks analysed in this study. It is important to note that the detection of submergence from the first field survey was only possible when adopting a millimetre-scale sampling strategy. When we merge the data for the 0-1 mm samples with those from 1-5 mm, to reflect a commonly used sampling strategy of 5 mm thick samples (e.g. Garrett et al., 2015;Shennan et al., 2018), and compare the change with respect to the 5-10 mm depth samples, the signal of submergence disappears from the samples from the first field survey. For example, all of the medians shown in Fig. 6, which are negative, indicating submergence, become positive when using 5 mm thick samples. We do not propose the adoption of millimetre-scale sampling strategies for palaeoseismological research but instead highlight the implications for records of earthquake activity from tidal wetland environments. We note that adoption of the standard sampling strategies further reduces the potential to detect low levels of land-level change in tidal wetland environments, with implications for the establishment of complete, longer-term records of megathrust earthquake activity.
The results from Isla de Chiloé therefore highlight the potential for the underestimation of magnitude 7-class earthquakes in palaeoseismological records from tidal wetland environments. Within the 1960 Valdivia M w 9.5 earthquake rupture zone, palaeoseismological studies have identified three predecessors over the past millennium that probably resembled the 20th century's largest earthquake (Cisternas et al., 2005Garrett et al., 2015;Moernaut et al., 2014). A further five earthquakes, including historically documented ruptures in 1737 and 1837, were probably also great (>M w 8) earthquakes . Our evidence shows that the M w 7.6 Chiloé earthquake is not encapsulated within the tidal wetland sediment or diatom record, meaning that previous earthquakes of similar magnitude may also not be preserved due to the creation and/or preservation threshold not being exceeded. It is highly probable that magnitude 7-class earthquakes are underrepresented in palaeoseismic records from tidal wetland environments in the region, as these earthquakes are either close to or below the detection limit. Furthermore, as damage caused by such earthquakes may not be dramatic when they occur away from population centres, it is likely that they are also underrepresented in the historical record. Consequently, the recurrence interval of magnitude-7 class earthquakes may be difficult to assess.
Tidal marsh records from other coastal regions around the Pacific may similarly underestimate or not record the occurrence of major earthquakes. Consequently, if recurrence intervals are calculated or seismic hazard assessments are developed based on tidal marsh palaeoseismology; this limitation must be recognised and accounted for. Furthermore, the underrepresentation of major earthquakes limits our ability to use tidal marshes to assess the role of earthquakes of this magnitude in releasing accumulated strain over millennial timescales.
We suggest that other major and great earthquakes of known magnitude are investigated to determine the point at which records of earthquakes are created, preserved and distinguishable from non-seismic processes (Fig. 5). Our results add support to recommendations from investigations in Alaska that identification of decimetre-scale coseismic changes requires quantitative estimates of elevation change for multiple locations along the couplet or proposed earthquake horizon. We also recommend using more than one transfer function model, proposed coseismic submergence or emergence should be replicated both spatially and with different quantitative models.

Conclusions
Diatom evidence from tidal wetlands indicates vertical landlevel motion during the 2016 M w 7.6 Chiloé earthquake, consistent with model estimates of low-level subsidence of <0.1 m. Quantitative analyses of diatom data show that records of the earthquake were created but not preserved within sediment sequences. Transfer function models provide a range of trajectories of change (positive, negative and close to zero) for individual sediment blocks and while there was a consistent pattern of statistically significant change between the uppermost samples which we could identify seven months after the earthquake, we could not identify a consistent signal after nine months. The deformation magnitude was at or below the detection limit for sites of this type, as coseismic, postseismic and interseismic changes cannot be distinguished in the sediment record. The most sensitive part of the tidal marshes on the southern coast of Isla de Chiloé is variable and may be dependent on sedimentation after the earthquake. In addition, the transient records of submergence seen after seven months are not visible when the standard sampling strategy for palaeoseismological research is adopted due to the low-level change identified, with important implications for the detection of earthquakes in tidal wetland records.
Our results support the third criterion for the detection of low-level coseismic land-level changes developed in Alaska for samples collected seven months after the earthquake (Shennan et al., 2014(Shennan et al., , 2016(Shennan et al., , 2018: quantitative estimates of elevation change across the contact all indicating the same sign (emergence or submergence) rather than a random distribution of mean values either side of zero.
Sampling strategies can be guided by changes in sediment stratigraphy and tsunami deposits, although neither are evident in the sediment blocks analysed as part of this study (Criterion A). We note that statistical analysis of change between marsh surface elevation reconstructions can be a powerful tool for the assessment of land-level motion (Criterion B). In turn, these analyses provide an important insight into the consequences of sampling strategy and timing of sample collection. There are important implications for records of palaeoseismicity from tidal wetlands, which are highly likely to underrepresent the occurrence of major earthquakes due to the creation and/or preservation threshold not being exceeded. Consequently, attempts to assess future seismic hazards in south-central Chile and elsewhere must recognise that tidal marshes provide an incomplete record of major earthquakes. of the Monte Horeb. The authors would also like to thank Neil Tunstall and Chris Longley for their assistance with X-ray imaging. This is a contribution to International Geoscience Programme (IGCP) Project 639. The authors declare no conflicts of interest.

Data availability statement
The data that support the findings of this study are available in the supplementary material for this article.

Supporting information
Additional supporting information can be found in the online version of this article. Figure S1. Comparison of tide gauge data collected using a portable ultrasound tide gauge (blue), tidal observations collected using a level and staff (black outlined circles) and tidal predictions (red) from the TPXO 8 Atlas tidal model (Egbert and Erofeeva, 2002). Tidal model amplitude scaling to reduce the misfit between model and observations: Ayentema = 1.1, Asasao = 1.1, Quilanlar = 1.05.
Appendix S2. Diatom assemblage diagrams for each field sample Appendix S3. Raw diatom counts for each field sample