Coastal Evolution in a Wetland Affected by Large Tsunamigenic Earthquakes in South-Central Chile: Criteria for Integrated Coastal Management

: The coastal evolution of the microtidal Tubul-Raqui wetland in south-central Chile (36 ◦ S), which historically has been affected by large earthquakes and tsunamis, particularly the 1960 (Mw = 9.5) and 2010 (Mw = 8.8) subduction earthquakes and their associated tsunamis, is analyzed. Historical aerial photographs and topographic and bathymetric surveys from the 1961–2017 period, as well as salinity, sediment, and ﬂora data obtained following the 2010 earthquake were used for comparison with data from prior to the event. A steady state of the shoreline was established, with an average erosion rate of − 0.016 m/year in the 1961–2017 period. However, erosion pre-dominated in the period between these two large earthquakes (1961–2009), with an average rate of − 0.386 m/year. The wetland dried up, partially recovered saline intrusion a year later, and recovered the salinity conditions it had before the earthquake two years later. The postearthquake effects on the ﬂoristic composition were not signiﬁcant, with the species Spartina densiﬂora , which presented a high tolerance to these types of changes, predominating. Moreover, 75 percent of the taxa in pre- and postearthquake conditions coincided, with the halophyte species Spartina densiﬂora , Sarcocornia fructicosa , and Cotula coronopifolia predominating, while the best-conserved community was Spartina-Sarcocornia association located in the saltmarsh. Seven years after the earthquake, the shoreline presented an accretion rate of 2.935 m/year; if the current tectonic conditions prevail, an erosive trend can be expected in the coming decades. The morphological variability and the changes associated with the shoreline in this wetland are strongly controlled by tectonic factors. Criteria aimed at integrated coastal management to promote its occupancy and use in accordance with its evolutionary dynamics are proposed.


Introduction
Coastal wetlands are currently the object of worldwide concern, given the heavy anthropization of the coastal zone and the need for climate change adaptation mechanisms. In addition to being the most biologically productive ecosystems and important sources of biological diversity, water reservoirs are responsible for high primary productivity [1]. Due to their fragility and delicate ecological equilibrium, they are susceptible to anthropogenic and natural disturbances that determine their functioning and variability [2].
In this respect, research on natural disturbances associated with earthquakes and tsunamis that affect coastal wetlands on a geological scale is well represented (e.g., coseismic subsidence rates for large earthquakes, their relationships with paleotsunamis, and ancient sea levels) [3][4][5][6][7][8][9][10][11][12][13][14]. However, for historical and short-term scales, research has been limited, and there has been little systematization of the function of violent disturbances such as earthquakes and tsunamis in recovery and change adaptation mechanisms in coastal wetlands [15][16][17][18][19]. An understanding of such changes and the ways in which the coast readapts is fundamental for guiding the use and forms of occupancy of the coast, given the increase in extreme events in the last decade (typhoons, hurricanes, meteotsunamis, storm surges) and the increasingly destructive damage caused by earthquakes and tsunamis [20][21][22][23][24][25][26]. Moreover, these adaptive characteristics, in addition to maintaining the essential functions of the coast, are an indicator of system resilience [27,28]. Studies that integrate the physical, ecological, biogeochemical, and geomorphic responses of estuaries to sea level rise [29] or any other type of change in water level are fundamental to bridging the remaining knowledge gap for better coastal zone management.
In Chile, the study of coastal wetlands has been approached mainly in terms of biological diversity, in keeping with the worldwide trend [29][30][31][32][33][34][35][36][37], and the effects of urbanization [38][39][40]. Dura et al. [14] recently studied the frequency of large tsunamigenic earthquakes using records of sediment deposits in ancient wetlands of the coast of central Chile (33 • S). Nonetheless, little is known about the natural change factors in these environments that would allow their integration into coastal management plans such that uses could be planned in accord with the functioning and natural dynamics of the wetlands that are continually exposed to earthquake and tsunami phenomena.
In central Chile (32 • -37 • S), the Mw = 8.8 earthquake of 27 February 2010 and the associated tsunami generated widespread destruction along 600 km of the coast, causing major morphological, social, and environmental changes [41][42][43]. Coastal uplift of up to 1.8 m occurred on the coast of Lebu (37 • S), uncovering shore platforms and resulting in changes in the intertidal zone and drying of wetlands [44][45][46][47].
One of the most affected coastal wetlands was Tubul-Raqui (36 • S), the aquatic components of which underwent significant modifications, especially in the intertidal marshes, which partially dried up [48][49][50]. Here, the magnitude of coastal uplift was 1.4 m, which interrupted the entry of the tide that allowed the saltmarsh to be maintained, causing a significant ecosystem impact, mainly due to the mortality of benthic macroinvertebrates [44][45][46]51,52]. According to [53,54], coastal uplift produced limitations on the exchange of river flow with the ocean, causing a decrease in the salinity of the water in the wetland, severely affecting the two most important benthic species for artisanal fishing, the alga pelillo (Gracilaria chilensis) and the bivalve razor clam (Ensis macha), with no evidence of any degree of recolonization.
Since the 1960 earthquake, the Chilean coast had not been affected by a large-magnitude earthquake and tsunami with such violent modifications of the coastal geomorphology, which has brought about the need to investigate the change factors of a tectonic coast to generate applications on a land-use planning and coastal management scale since the urbanization processes underway in the country have developed around coastal axes, forming new metropolitan areas and in turn new coastal axes [55,56]. In addition, the need for adaptation mechanisms requires an understanding of how to plan human settlements and urban expansion in safe areas by considering given climate change scenarios [57]. In rural communities and especially fishing coves, these mechanisms are urgent due to their dependence on natural resources and the need to conserve coastal marine ecosystems.
Thus, the purpose of this investigation is to analyze historical-scale coastal evolution in the Tubul-Raqui wetland during the period between the two large subduction earthquakes (1960 and 2010) in order to propose criteria for integrated coastal management in coastal areas affected on a recurring basis by these phenomena. The paper also includes the effects of the reconstruction process following the 2010 earthquake, 10 years on from the disaster, to illustrate how an understanding of the changing factors of this type of coast can influence improvements in decision making and contribute to planning for coastal resilience by allowing the reorganization and adaptation of the natural system and its ability to maintain its functions or regain them to be considered [58].

Geographical Context
On the coast of the Biobío Region (36 • -38 • S), especially on the coastal plain, there are numerous extensive wetland zones that make up valuable environmental heritage. Among them, the Tubul-Raqui wetland stands out for its high biodiversity [49] and as the foundation for the local economy of a small rural community [59]. This coastal wetland of south-central Chile is characterized by a large saltmarsh with the development of esparto and vegetation adapted to salinity variations [29,[60][61][62]. It is surrounded by shore platforms up to 300 m and intermittently flooded plains used for livestock activity or subsistence agriculture.
The Tubul-Raqui wetland (37 • 13 S-73 • 26 W) is formed by the confluence of the Tubul and Raqui rivers, which flow through two interconnected coastal basins of about 274 km 2 that drain directly into the Gulf of Arauco (Figure 1, ±1.0 m of tidal range) [63]. The area has a Csb Mediterranean climate (warm summer) with oceanic influence and rainfall concentrated in winter [64]. In geological terms, it is located in the sedimentary basin of the Arauco Peninsula, is from the Cretaceous age, and is composed of Quiriquina, Tubul, Pilpilco, and Ranquil formations [65].
Starting in the 1960s, an organized group of fishermen and alga (Gracilaria chilensis or "pelillo") collectors established itself here, forming a settlement of about 2000 inhabitants. Until the earthquake of 27 February 2010, this town was the center of social and economic organization of Tubul Cove [59,66]. shows the marine zone next off Tubul Cove in pre-earthquake conditions, which is now an intertidal zone of Tubul Beach due to coastal uplift produced by the earthquake. Here, a beach approximately 100 m wide was formed (B); (C) shows the marine abrasion platform that emerged due to coastal uplift measured at 1.4 m by Quezada [44] after the 2010 earthquake.
Starting in the 1960s, an organized group of fishermen and alga (Gracilaria chilensis or "pelillo") collectors established itself here, forming a settlement of about 2000 inhabitants. Until the earthquake of 27 February 2010, this town was the center of social and economic organization of Tubul Cove [59,66].

Materials and Methods
Three groups of data were collected for pre-and post-2010 earthquake conditions: (a) topography and bathymetry to determine morphological changes; (b) a salinity parameter (PSU), used to demonstrate the recovery of saline intrusion; and (c) wetland floristic and vegetation composition.
In addition, through the coupling of the process-based model Delft 3D FLOW [67] for morphological response and the wave propagation model Delft 3D WAVE [68], the changes in the morphological response of the coastline induced by coseismic shifts occurring during the 27 February 2010 earthquake were evaluated for regular and storm wave conditions.  , and (C) (view from river mouth to point C toward NW) show pre-and postearthquake morphological changes, indicating their location on the coastal area of the wetland; (A) shows the marine zone next off Tubul Cove in pre-earthquake conditions, which is now an intertidal zone of Tubul Beach due to coastal uplift produced by the earthquake. Here, a beach approximately 100 m wide was formed (B); (C) shows the marine abrasion platform that emerged due to coastal uplift measured at 1.4 m by Quezada [44] after the 2010 earthquake.

Materials and Methods
Three groups of data were collected for pre-and post-2010 earthquake conditions: (a) topography and bathymetry to determine morphological changes; (b) a salinity parameter (PSU), used to demonstrate the recovery of saline intrusion; and (c) wetland floristic and vegetation composition.
In addition, through the coupling of the process-based model Delft 3D FLOW [67] for morphological response and the wave propagation model Delft 3D WAVE [68], the changes in the morphological response of the coastline induced by coseismic shifts occurring during the 27 February 2010 earthquake were evaluated for regular and storm wave conditions.

(a) Morphological Changes
The morphological changes in the wetland were determined based on shoreline variations in the 1979-2017 period, that is, before and after the 27 February 2010 earthquake conditions. The methodology proposed by Martinez et al. [69], employing the use of aerial photographs from different years and precise topographic surveys for the generation of a historical database, was used. The aerial photographs were georeferenced and digitized, with mean quadratic errors (RMS) of less than 1 m obtained ( Table 1). The topographical surveys were carried out in winter and summer between 2010 and 2015. For this purpose, Water 2021, 13, 1467 5 of 27 topographic instruments (GPS, Trimble, model R-4) with submeter precision (mm) were used. The changes in the shoreline were determined by superimposing maps in the GIS platform (ArcGis 9.3), using three equidistant segments every 100 m. Three beach profiles (P1, P2, P3) were also located, in which sediment was collected from the surf, beach-face, and backshore zones ( Figure 2). For each zone, a surface sediment sample (120 gr) was collected manually and stored in a plastic bag labeled with an alphanumeric code. The sediment was analyzed by sieving in the Physical Geography Laboratory of the Universidad de Concepción. The sediments were sieved using a Ro-Tap mechanical stirrer (RETSCH Mod. AS 200 brand) and dried at a constant temperature of 100 • C in a universal forced convection oven (MEMMERT Mod. ULM 700 brand) [70]. For the textural size classification, the Wentworth [71] classification was used. The mean diameter (mean size) and the dispersion parameters (asymmetry) were obtained by statistical analysis using the criteria of Folk and Ward [72]. The Trask (So) selection or classification coefficient (sorting index, therefore, =D 75 /D 25 ) was also included [73]. Table 1. Aerial photographs used in the determination of the planform changes of the shoreline and the mean squared error. In general, all of the errors were less than 1 m.

Name
Year For the determination of vertical changes at the bottom of the Tubul River and Bay, detailed bathymetries were carried out using an echosounder in the summers of 2011, 2012, and 2013. A portable echosounder with built-in GPS (GARMIN brand, model FF 400C GPsMap 420s) and a maximum depth of 400 m was used, which was installed aboard a boat. The depth data were referenced to a Chart Datum (CD), estimated by applying the approach suggested by Pawlowicz et al. [74] using two HOBO pressure sensors installed in the Tubul River and the marine zone off Tubul Cove (denoted by a star in Figure 2A). These sensors recorded tidal variations for at least 15 days, with a recording interval of 5 min, recording the spring tides of the month. For pre-earthquake conditions, a precise bathymetric survey carried out in summer 2008 by the Hydraulic Works Department of the Ministry of Public Works (DOH-MOP) [75] was available. Using both surveys, vertical changes were compared through two topo-bathymetric beach profiles, denoted by PB-1 and PB-2 (see Figure 2A).

(b) Salinity and Recovery of Marine Intrusion
Salinity (practical salinity units (PSU)) measurements were taken using a multiparameter probe (Horiba U-10) at six stations in the Tubul River and seven stations in the Raqui River ( Figure 2B). These data were collected in spring tide conditions for a tidal cycle according to the methodology of Perillo and Piccolo [76]. To determine the presence of a saline wedge in the estuary, salinity was measured in both the deep and surface layers of the water. The stations that presented a depth of less than 0.5 m were considered the surface layer [76]. To classify the data, the criteria of Madden et al. [77] were used. To analyze possible changes in the shape of tidal waves inside the wetland, two HOBO pressure sensors were installed in the Tubul River and the adjacent marine zone during the periods from 18 March to 19 April 2011, and 16 December 2011 to 16 January 2012. These data were compared with tide data from the Tubul River from before the earthquake (17-31 January 2008), available from the Ministry of Public Works (MOP) [75]. (b) Salinity and Recovery of Marine Intrusion Salinity (practical salinity units (PSU)) measurements were taken using a multiparameter probe (Horiba U-10) at six stations in the Tubul River and seven stations in the Raqui River ( Figure 2B). These data were collected in spring tide conditions for a tidal cycle according to the methodology of Perillo and Piccolo [76]. To determine the presence of a saline wedge in the estuary, salinity was measured in both the deep and surface layers of the water. The stations that presented a depth of less than 0.5 m were considered the For the tide analysis, the methodology proposed by [74] was used, which employs routines in MATLAB that are a translation of the original routines of Foreman 1977 in [74] FORTRAN program. Harmonic and nonharmonic analyses were carried out to determine their regime [78]. With the harmonic constituents obtained, a forecast for the same period was made to compare these values with the observed values. To this end, a residual and  Figure 2B). The samplings were carried out at exactly the same geographical coordinates, randomly distributed throughout the wetland area.
The vegetation sampling was carried out according to the Braun-Blanquet [79] approximation, categorizing each taxon according to a total coverage scale in quadrants, with 10 one-square-meter quadrants. For measuring abundance, the coverage percentage of the emergent marsh species, which was obtained by observing the area of the quadrant occupied by each plant species, was used. To identify the species, material samples were collected for subsequent identification, using reference information [60,[80][81][82][83][84], and direct comparison with specimens deposited in the CONC herbarium of the Universidad de Concepción. For each sampling station, the community parameters of richness (S), abundance (N), Simpson diversity index (1-D) [85], Shannon diversity index (H') [86] Pielou evenness (J') [87], and dominance (D) were estimated. To this end, the Diversity module of PAST [88] was used.

Numerical Modeling
The numerical model used in this investigation was the Delft 3D model developed by Deltares, which allows 1D, 2DV, 2DH, and 3D modeling of both hydrodynamics and the morphological evolution of the beach zone. The model uses nested grids to increase the resolution in the zones of greatest interest and couples two calculation modules-SWAN [89] and FLOW [67]-to incorporate wave propagation and then morphological and hydrodynamic effects.
To quantify the hydrodynamic effects, the model solved the unsteady shallow water equations. Sediment transport was calculated using the Van Rijn [90] equation for total transport based on the integrated vertical advection-diffusion equation. Morphological changes occur on time scales between 1 and 2 orders of magnitude greater than those of hydrodynamic changes [91]. Therefore, sea bottom changes were calculated by amplifying the sediment transport effects, estimated at each hydrodynamic time step. This amplifying factor, called MORFAC [92], varied by analysis type, that is, different MORFACs were applied for regular and storm waves. This approximation allowed us to accelerate morphological changes to hydrodynamic times and thereby immediately update the bathymetry after each iteration. The numerical model used considered the entire Gulf of Arauco area, with seven sectors modeled in detail. The model details and calibration were described by Gomez et al. [93], and in this investigation, only the results associated with the Tubul-Raqui wetland mouth were used.
In Figure 3, the dimensions of the nested grids used to model wave propagation and the smaller grids associated with each beach for studying morphological effects are shown. The detail grid for the Tubul sector has a resolution of 15 m × 20 m, ensuring a resolution sufficient for detecting local morphological changes. The tide was modeled based on 15 main tide components. In addition, the different wave conditions were determined from a statistical analysis of deep-water waves extracted from the Wave Watch III [93,94] data series for node 37 • 0 S 74 • 0 O, located off the coast of the Gulf of Arauco. Wave climate details can be found in Gomez et al. [93]. The summary parameters used for the modeling, their corresponding MORFAC, and modeling times can be seen in Table 2. shown. The detail grid for the Tubul sector has a resolution of 15 m × 20 m, ensuring a resolution sufficient for detecting local morphological changes. The tide was modeled based on 15 main tide components. In addition, the different wave conditions were determined from a statistical analysis of deep-water waves extracted from the Wave Watch III [93,94] data series for node 37°0 ′S 74°0′ O, located off the coast of the Gulf of Arauco. Wave climate details can be found in Gomez et al. [93]. The summary parameters used for the modeling, their corresponding MORFAC, and modeling times can be seen in Table 2.  To consider the effect of coseismic shifts, a seafloor uplift of 1.4 m was generated, in accordance with the observations of Martinez et al. [19] for this area. The simulations considered the same spectral wave parameters and tidal forcings but with different depths in order to determine morphological behavior changes in the area before and after the earthquake.

Pre-and Postearthquake Geomorphological Changes
The shoreline of the Tubul-Raqui wetland (segment 1, or marine segment, extends from point B to point A' in Figure 4) presented an average variability of −28.27 m in the last 36 years and a mean retreat rate of −1.32 m/year. On the banks of the Tubul River (segment 2, or estuarine segment, extends from point A to point B in Figure 4), the mean advance rate was 0.57 m/year ( Figure 4). The changes in the shoreline indicate major morphological variability, especially in low-tide conditions ( Figure 5).  To consider the effect of coseismic shifts, a seafloor uplift of 1.4 m was generated, in accordance with the observations of Martinez et al. [19] for this area. The simulations considered the same spectral wave parameters and tidal forcings but with different depths in order to determine morphological behavior changes in the area before and after the earthquake.

Pre-and Postearthquake Geomorphological Changes
The shoreline of the Tubul-Raqui wetland (segment 1, or marine segment, extends from point B to point A' in   In Figure 5A, the period between two large earthquakes (post-1960 event) was analyzed; a strong erosive trend was observed, with shoreline retreat values between −1.5 and −0.2 m/year, consistent with the interseismic period. In Figure 5B, the postearthquake period (2010-2017) was analyzed, in which the coastal uplift caused a widening of the intertidal zone and accretion processes, explaining the accretion trend in most of the estuarine zone. Even so, some parts of the exposed beach present values corresponding to stability and erosion. In Figure 5A, the period between two large earthquakes (post-1960 event) was analyzed; a strong erosive trend was observed, with shoreline retreat values between −1.5 and −0.2 m/year, consistent with the interseismic period. In Figure 5B, the postearthquake period (2010-2017) was analyzed, in which the coastal uplift caused a widening of the intertidal zone and accretion processes, explaining the accretion trend in most of the estuarine zone. Even so, some parts of the exposed beach present values corresponding to stability and erosion.
The bathymetric changes in Tubul Bay established an average variability of −0.45 m in the PB-1 profile in the 2008-2013 period (Table 3 and Figure 6). In general, the topo-  (Table 3 and Figure 6). In general, the topobathymetric profiles presented wide variability in postearthquake conditions (Figures 7 and 8), which was in turn reflected in beach profiles P1 and P2. Profile P3, because it is located in an area that is more protected from wave action, did not present major geomorphological changes during the year. bathymetric profiles presented wide variability in postearthquake conditions (Figures 7  and 8), which was in turn reflected in beach profiles P1 and P2. Profile P3, because it is located in an area that is more protected from wave action, did not present major geomorphological changes during the year.
In the three analyzed profiles, the mean beach particle size fluctuates between medium and fine sand for the surf and beach-face zones, with fine sand predominating in profile P3. The surf zone sediment was leptokurtic and with moderate sorting, while that of the beach face was mesokurtic and moderately well sorted ( Figure 9 and Table 4).   The average size in all profiles corresponded to medium or find sand.    In the three analyzed profiles, the mean beach particle size fluctuates between medium and fine sand for the surf and beach-face zones, with fine sand predominating in profile P3. The surf zone sediment was leptokurtic and with moderate sorting, while that of the beach face was mesokurtic and moderately well sorted (Figure 9 and Table 4).

Figure 8.
Beach profiles, Tubul Bay. The beach profiles were surveyed in postearthquake conditions on a seasonal basis (winter-summer). They present wide variability (A,B), with the exception of profile 3 (C), which is located in a more protected area of the bay. See Figure 2A for the location of each profile.   Table 4. Statistical parameters of the beach-face and surf zone sediment. Mean values, kurtosis, and sorting (values lower than 1.6 indicate a good classification, values between 1.6 and 1.9 indicate a moderate or medium rating, and values above 1.9 indicate a poor rating) for the three beach profiles are presented in the table. In general, the surf zone sediment presented leptokurtosis and moderate sorting, while the beach-face sediment was mesokurtic and moderately well sorted. The average size in all profiles corresponded to medium or find sand.

Salinity and Recovery of Marine Intrusion
The effect of coastal uplift caused by the 27 February 2010 earthquake on the wetland was observed until 2011, the time when the salinity values were minimal at all stations. The Tubul River was only accessible up to station 2, while the low depth of the Raqui River made it impossible to reach any station by boat ( Figure 10). In postearthquake conditions, the greatest salinity values were recorded in austral summer 2012 (February and March), when tidal wave intrusion was recorded up to 6 km upriver. These conditions remained during 2012, although with a significant variation in salinity values, until decreasing notably in 2014 in greater river streamflow conditions. was observed until 2011, the time when the salinity values were minimal at all station The Tubul River was only accessible up to station 2, while the low depth of the Raqu River made it impossible to reach any station by boat ( Figure 10). In postearthquake con ditions, the greatest salinity values were recorded in austral summer 2012 (February an March), when tidal wave intrusion was recorded up to 6 km upriver. These condition remained during 2012, although with a significant variation in salinity values, until de creasing notably in 2014 in greater river streamflow conditions.  stations (a,b). Even so, this was a short-term effect because, by early 2012, salinity was instrumentally recorded 6 km upstream in the Tubul River. At stations with a water depth greater than 0.5 m (Raqui River and stations 3 to 7 in the Tubul River), a slight recovery of marine intrusion was observed over time.
Harmonic analysis of tidal signals shows Courtier's "F" coefficient values of 0.43 (2011) and 0.404 (2012), indicating a mixed semidiurnal tide regime. The nonharmon  (a,b). Even so, this was a short-term effect because, by early 2012, salinity was instrumentally recorded 6 km upstream in the Tubul River. At stations with a water depth greater than 0.5 m (Raqui River and stations 3 to 7 in the Tubul River), a slight recovery of marine intrusion was observed over time.
Harmonic analysis of tidal signals shows Courtier's "F" coefficient values of 0.435 (2011) and 0.404 (2012), indicating a mixed semidiurnal tide regime. The nonharmonic values for pre-and postearthquake conditions were different, with greater depths observed before the earthquake, meaning that the tidal cycle penetrates differently relative to pre-2010 conditions in the estuary.

(a) Floristic Composition
In 2009, we collected a total of 19 plant species, belonging to 10 families. The main species detected correspond to the Poaceae (48%, four species), Asteraceae (22%, two species), Juncaceae (12%, three species), and Fabaceae (10%, two species). Overall, 65% of the wetland's vegetation cover was composed of native plants and 35% of introduced plants. Then, in the fieldwork carried out during 2012, we found a total of 26 plant species, belonging to 10 families. The main species detected correspond to the Poaceae (31%, five species), Asteraceae (22%, four species), Fabaceae (14%, four species), Juncaceae (13%, three species), and Cyperaceae (9%, three species). Overall, 44% of the wetland's vegetation cover was composed of native plants and 42% of introduced plants. Later, during 2014, we collected a total of 25 plant species, belonging to the same 10 families. The main species detected correspond to the Poaceae (31%, six species), Asteraceae (15%, three species), Juncaceae (13%, three species), Cyperaceae (13%, two species), and Fabaceae (8%, three species). Overall, 56% of the wetland's vegetation cover was composed of native plants and 33% of introduced plants. Finally, in 2016, we found a total of 31 plant species, belonging to 11 families. The main species detected correspond to the Poaceae (44%, nine species), Juncaceae (12%, three species), Cyperaceae (9%, four species), and Fabaceae (8%, two species). Overall, 59% of the wetland's vegetation cover was composed of native plants and 36% of introduced plants (Supplementary Materials, Table S1 and Figure 11). and Sarcocornia fructicosa, which is the only association found without major human interventions due to the extreme characteristics of the habitat.  Most of the collected flora was of herbaceous habit, mostly hemicryptophytes and chamaephytes, which confers ideal attributes for resisting the salinity stress of the saltmarsh and poor drainage.

(b) Vegetation Structure
The vegetation structure of the Tubul-Raqui wetland is constituted by nonarboreal emergent communities of vascular plants, with roots that occur in areas in which salinity range from very high (26 PSU) near the saltmarsh to very low (1 PSU) in higher zones, where fresh water dominates, and other zones, where periodic pools generated by tides occur, giving rise to very different geomorphological units according to vegetation type (coastal saltmarsh, floodplain, and palustrine wetlands).
The coastal saltmarsh proved to be extremely homogenous, dominated by the tallgrass Spartina densiflora and with the occasional participation of Sarcocornia fruticosa, both of which adapted to the high salinity of the substrate and surrounding water, an association considered typical of saltmarshes by San Martín et al. [29].
The floodplains were characterized by a greater diversity of plant species. This characteristic depended directly on the lower saline influence in the floodplains and the greater extent of human activity in the area. These factors promoted the establishment of a higher number of plant species, allowing the survival of halophyte species such as Cotula coronopifolia and the introduction of forage species such as Trifolium repens. In addition, depending on the drainage conditions of the plain, it was possible to identify native hydrophilic species such as Scirpus californicus and Juncus procerus.
A general zoning diagram establishes that the peripheral areas of the floodplains and ravines that drain into the wetland are covered by a low Myrtaceae forest, giving rise to an ecotonal zone dominated by a prairie with reeds and rushes, with variable coverage, and a low prairie with Cotula coronopifolia in the salt plains, while the lower portion is composed of a medium-tall grassland of halophyte plants dominated by Spartina densiflora and Sarcocornia fructicosa, which is the only association found without major human interventions due to the extreme characteristics of the habitat.

(d) Conservation Status
The Tubul-Raqui wetland did not present species with conservation problems since many of them are widely distributed in other saltmarsh and palustrine wetland systems [29,60,62]. Human intervention such as livestock activity and farming was observed in the field. Home construction is another reason for which it is practically impossible to find remnants of original vegetation in these areas, with the best-conserved community being the Spartina-Sarcocornia association in the lower portion of the wetland since the habitat conditions impede its use by humans.

Numerical Modeling Results
The numerical modeling results show the change in the morphological response of the sector before and after the coseismic movements on 27 February 2010. Two types of responses can be observed, the beach response to regular waves and the beach response to storm waves.

(a) Regular Waves
The modeling results for regular wave conditions, distinguished between winter and summer, do not show large differences between the seasons. However, when the uplift of 1.4 m in the area [95] is incorporated, the changes in the accretion and erosion patterns of the areas can be clearly identified (see Figure 13).

(d) Conservation Status
The Tubul-Raqui wetland did not present species with conservation problems since many of them are widely distributed in other saltmarsh and palustrine wetland systems [29,60,62]. Human intervention such as livestock activity and farming was observed in the field. Home construction is another reason for which it is practically impossible to find remnants of original vegetation in these areas, with the best-conserved community being the Spartina-Sarcocornia association in the lower portion of the wetland since the habitat conditions impede its use by humans.

Numerical Modeling Results
The numerical modeling results show the change in the morphological response of the sector before and after the coseismic movements on 27 February 2010. Two types of responses can be observed, the beach response to regular waves and the beach response to storm waves.

(a) Regular Waves
The modeling results for regular wave conditions, distinguished between winter and summer, do not show large differences between the seasons. However, when the uplift of 1.4 m in the area [95] is incorporated, the changes in the accretion and erosion patterns of the areas can be clearly identified (see Figure 13). One of the main changes that can be observed is a certain accretion trend that can be directly observed at the mouth of the wetland, which tends to close and thus decrease the One of the main changes that can be observed is a certain accretion trend that can be directly observed at the mouth of the wetland, which tends to close and thus decrease the tidal exchange of the wetland.
Another observation from the obtained results is related to a clear recovery of the intermediate part of the beach from sediment originating in deeper water, as can be observed in Figure 14, which shows that the strongest erosion occurred in deep water areas. This change, if extrapolated over time, will tend to form a wider beach with a lower slope than the beach before 27 February 2010.

Discussion
Large-magnitude natural disturbances such as earthquakes and tsunamis are recurring phenomena on the Chilean coast. This historical recurrence has been well documented, with evidence of events that have been classified as gigantic due to their magnitudes (M≥9.0) [96][97][98][99]. The last major earthquake (Mw = 8.8 in 2010) heavily affected the Tubul-Raqui area and the entire coast near the southern rupture segment [46,51]. Since Tubul-Raqui is a coastal wetland, the effects were varied, affecting the morphology, tidal wave intrusion, and biota. The types of changes undergone in this wetland, according to the collected data, can be interpreted as a function of three principal aspects:

(a) Geomorphological Changes and Marine Intrusion Effects
Microtidal estuarine environments such as Tubul are characterized by a significant natural variability of forms; however, the role played by earthquakes and tsunamis in this variability on medium-and short-term scales is little understood [19]. Two towns on this coast-Lebu and Tubul-have been affected by significant coastal uplift generated by the 27 February 2010 earthquake, which has allowed an understanding of the postearthquake recovery process. In the case of the first, the uplift was 1.4 m, which increased the width of the intertidal zone of the beach by 100 m [51,100]. However, this did not change the position of the shoreline, and thus, during the last 36 years, the shoreline of the Tubul-Raqui wetland presented a mean retreat rate of −1.32 m/year. Meanwhile, the Lebu pocket beach (37° S), affected by an uplift of 1.7 m, presented a mean accretion rate of 2.80 m/year in the 1984-2014 period. Here, there was a rapid recovery of the prior morphodynamic characteristics of the beach [19].
The differences in the mean rates of change in the towns of Tubul and Lebu can be explained by the greater uplift caused by the 1960 earthquake in Lebu. On the northern shore of the Arauco Peninsula, the 1960 earthquake generated an uplift of 0.4 cm [101]. The towns of Lebu, Llico, Punta Lavapié, and Tubul, located closer to the trench, underwent greater interseismic subsidence, with which the shoreline tends to recede due to the ensuing erosion process [95]. In Tubul, an area affected by interseismic subsidence, the predominant process is erosion generated by tectonic processes that naturally induce shoreline retreat. The coseismic movements in this area have not allowed the wetland to compensate for the erosion process, in contrast to what occurred in Lebu as a result of the coastal uplift of 1960 and 2010.
This uplift caused a significant geomorphological change in the wetland upon inter- The changes observed for storm waves ( Figure 14) are very significant for both conditions (pre-and postearthquake). This is because the initial beach profile for the calculations was the equilibrium profile for regular wave conditions. However, important differences in the evolution of the beach are observed in each case, heavily increasing the eroded zone in the post-2010 condition. It is also observed that the previous condition did not generate major erosion beyond the coastline since the line is somehow in equilibrium for storm conditions, but it did generate greater sedimentation in the mouth that tended to block the tidal exchange to and from the wetland. Furthermore, the new bathymetric scenario created following 27 February 2010 with the subsidence requires a greater amount of sediment from the coast to balance the storm wave energy, which generates a reset of the wetland mouth and is perceived as large erosion areas, but the zones located on the eastern part of the mouth remain little altered.

Discussion
Large-magnitude natural disturbances such as earthquakes and tsunamis are recurring phenomena on the Chilean coast. This historical recurrence has been well documented, with evidence of events that have been classified as gigantic due to their magnitudes (M ≥ 9.0) [96][97][98][99]. The last major earthquake (Mw = 8.8 in 2010) heavily affected the Tubul-Raqui area and the entire coast near the southern rupture segment [46,51]. Since Tubul-Raqui is a coastal wetland, the effects were varied, affecting the morphology, tidal wave intrusion, and biota. The types of changes undergone in this wetland, according to the collected data, can be interpreted as a function of three principal aspects:

(a) Geomorphological Changes and Marine Intrusion Effects
Microtidal estuarine environments such as Tubul are characterized by a significant natural variability of forms; however, the role played by earthquakes and tsunamis in this variability on medium-and short-term scales is little understood [19]. Two towns on this coast-Lebu and Tubul-have been affected by significant coastal uplift generated by the 27 February 2010 earthquake, which has allowed an understanding of the postearthquake recovery process. In the case of the first, the uplift was 1.4 m, which increased the width of the intertidal zone of the beach by 100 m [51,100]. However, this did not change the position of the shoreline, and thus, during the last 36 years, the shoreline of the Tubul-Raqui wetland presented a mean retreat rate of −1.32 m/year. Meanwhile, the Lebu pocket beach (37 • S), affected by an uplift of 1.7 m, presented a mean accretion rate of 2.80 m/year in the 1984-2014 period. Here, there was a rapid recovery of the prior morphodynamic characteristics of the beach [19].
The differences in the mean rates of change in the towns of Tubul and Lebu can be explained by the greater uplift caused by the 1960 earthquake in Lebu. On the northern shore of the Arauco Peninsula, the 1960 earthquake generated an uplift of 0.4 cm [101]. The towns of Lebu, Llico, Punta Lavapié, and Tubul, located closer to the trench, underwent greater interseismic subsidence, with which the shoreline tends to recede due to the ensuing erosion process [95]. In Tubul, an area affected by interseismic subsidence, the predominant process is erosion generated by tectonic processes that naturally induce shoreline retreat. The coseismic movements in this area have not allowed the wetland to compensate for the erosion process, in contrast to what occurred in Lebu as a result of the coastal uplift of 1960 and 2010.
This uplift caused a significant geomorphological change in the wetland upon interrupting tidal wave intrusion and resulting in the drying of the saltmarsh. As indicated by the collected data, the tide was impeded for almost a year, and subsequently, its entry was weak, according to the reported salinity data. Other studies also indicate that tidal wave intrusion was partially recovered 10 months after the earthquake, especially within three kilometers of the mouth, although with values lower than those recorded before the earthquake [48]. Starting in 2012, the salinity values were greater and comparable to those recorded before the earthquake [48]. In December 2008, salinity values of up to 32.2 PSU, 4 km from the mouth of the Tubul River, and 31 PSU, 8 km from it, were recorded [48]. These values are high and expected in the summer period with low river streamflows since both rivers have an entirely pluvial regime. The tidal wave entry capacity before the earthquake was high, with direct influence up 10 km upstream in both rivers [48,61].
In this investigation, the maximum recorded value was 27.3 PSU in summer (February 2012), 8 km upstream in the Tubul River. In the Raqui River, the maximum recorded salinity value was 27.9 PSU (station 4), 6 km upstream, also in summer conditions (March 2012).

(b) Changes in Floristic Composition and Vegetation Structure
Ecosystems are generally defined based on the plant species and principal characteristics of the site. The composition and structure of plant communities are a fundamental basis for determining the current conservation status and potential changes in restoration or conservation processes in an area [102].
The coastal saltmarsh and floodplains of the Tubul-Raqui area have different floristic compositions, a phenomenon explained by changes in salinity gradients and human disturbance. In the saltmarsh perennial species dominated, mostly halophytes, although there were others with a wide distribution without this trait [103]. The low number of species able to tolerate the high salinity of these ecosystems results in rather homogenous vegetation [103,104], with the native Poaceae Spartina densiflora, responsible for the formation of prairies or esparto fields, predominating. This species is of fundamental importance in these ecosystems, acting as an engineer species because it models and stabilizes the riverbanks against the effects of the tides and is the main energy source in the wetland, supplying large quantities of organic detritus that enter aquatic and terrestrial food chains [9].
However, in some sectors, it was also possible to find the succulent halophyte Sarcocornia fruticosa, which has been observed in other areas in association with Spartina densiflora [105], giving rise to complex vegetation mosaics called Sarcocornio-Spartinentum densiflorae [29], with Spartina densiflora much more abundant than the latter [106]. These areas are also characterized by little human intervention since there is no agriculture or livestock activity.
The proportion of native taxa was greater in undisturbed areas, in zones dominated by coastal saltmarshes where high salinity levels and extreme environmental conditions impede the proliferation of invasive plants and their use by humans; however, the greatest species diversity was observed in disturbed areas where introduced taxa contribute in large measure to species richness. The plant associations present in the Tubul-Raqui wetland are similar to those of other saltmarsh systems in south-central Chile [60][61][62].
The Tubul-Raqui wetland, similar to other saltmarshes of southern Chile, presents low levels of vegetation diversity [60], and the plant species there do not present conservation problems since they are widely distributed in other coastal wetlands; however, the various human activities such as livestock, agriculture and home construction have replaced the original vegetation of the floodplains and areas with lower salinity [49].
The vegetation structure comparison between the pre-earthquake (2009) and postearthquake (2012, 2014, and 2016) periods did not show differences, since all the stations analyzed remained without relevant vegetation changes. However, it is possible to observe slight changes in the wetland's floristic composition, mainly marked by a decrease in native flora from 65% during 2009 to 44% by 2012, followed by an increase in native plant species to 57% during 2016. These results suggest that the Tubul-Raqui wetland, although severely affected by the coastal uplift resulting from the 27 February 2010 earthquake, has gradually recovered, exhibiting no evidence of damage at the analyzed stations. It bears mentioning that this is the case only from a vegetation perspective and is probably due to the high tolerance to major environmental changes presented by the halophyte Spartina densiflora, which is very important in the physical, chemical, and biological balance and structuring of this ecosystem [49] Among wetlands on other tectonic coasts, the Copper River wetland (Alaska), affected in 1964 by an M = 8.4 earthquake that generated an uplift between 1.8 m and 3.4 m above the previous sea level, stands out. According to Thilenius [15], this uplift had a major impact on the ecosystem since the subtidal zones converted to in-tertidal zones; however, five years after the event, the nature of the taxa and vegetation structure, only was affected with slight changes in the location and abundance of plants.
In the Tubul-Raqui wetland, the vegetation structure as a foundation for a community of specialized organisms adapted to salinity gradients reacted differently to the coastal uplift caused by the earthquake, since its taxa were less affected and its recovery was faster once the tidal wave intrusion was restored. The foregoing coincides with the description of Valdovinos et al. [48], who established that the broad saline meadows or esparto fields dominated by Spartina densiflora did not present signs of having been affected by coastal uplift one year after the earthquake, which is important because it is a fundamental species in the structure of the wetland. According to the authors, this species stabilizes recent sediment deposits and contributes to the maintenance of significant biological diversity.
However, Valdovinos et al. [48,49,53] established that while there was not a significant impact on the coverage of the dominant species of the wetland (Spartina densiflora), there was significant damage to benthic invertebrates due to the drying of and lack of saline intrusion in the saltmarsh, especially the bivalve Ensis macha (navajuela) and Kingiella chilenica, the latter of which is endemic. Ten months after the earthquake, an intense recolonization of amphipods and Polychaeta in undried wetland channels was observed. In addition, the newly emerged beds were recolonized by vegetation and there was a return to the pre-earthquake estuarine conditions. Furthermore, Sandoval et al. [52] established that the most tolerant taxa to environmental perturbations (i.e., Diptera, Annelid, and Polychaeta) surpassed pre-earthquake abundance records after just ten months, whereas the most sensitive taxa were not found after the earthquake. Although the ecosystem was severely affected by the earthquake, it still has properties that make it among the most important in the country [53].

(c) Changes in Morphological Response Induced by Subsidence
Despite the changing trends in the wave climate, sediment transport off of the mouth of the Tubul-Raqui wetland was observed to be rather uniform, with a very prominent southwest direction. This homogenizes the wave patterns observed in the study area, reducing the impacts of changes in the main forcings. The shape of the Gulf of Arauco and the location of Tubul Beach at the bottom of the gulf determine the morphological response, as waves, although they come from different directions, reach the beach having undergone numerous bathymetric refraction processes and diffraction at Sta. Maria Island (see Figure 3) [93].
For the analyzed case at Tubul, the morphological changes observed before and after the earthquake are considerable given the uplift of the land: The morphology of the beach clearly went from a condition of quasistatic equilibrium (very little sediment transport for the predominant conditions) to a situation of change, in which a great quantity of sediment is transported from areas farther from the coast to the beach. The new trends of sediment transport increased the width of the beach, causing an eventual blockage of the wetland mouth that had repercussions for the exchange of water within it and the artisanal fishing activities of the nearby communities.
Due to the lack of obvious changes in principal forcings and the existence of measurements that indicate erosion conditions on Tubul Beach, it is difficult to suggest reasons for this behavior. The numerical modeling of the geomorphological evolution of the beach before and after 27 February 2010 corroborates that the observed erosion response could be related to the coseismic shifts measured in the sector and a reduction in beach recovery capacity due to heavy storms.
The morphological response of the coastline does not present great seasonal variations; thus, the beach should not vary from historical observations. However, storm waves will cause significant erosion in some sectors of the mouth, the repercussions of which will depend on the actual availability of sediment transported by the longshore drift toward the W-SW. The construction of new infrastructure (after 27 February 2010) at the northeast side of the beach that took place as part of the reconstruction process to provide greater functionality to the fisherman's cove could interrupt the normal flow of sediment and alter the sedimentary dynamics of the area, possibly bringing about as yet unknown consequences for the hydrodynamics of the wetland, which must be assessed to prevent further damage to its biota. Although the ecosystem was severely affected by the earthquake, it has maintained the properties that make it one of the most important coastal wetlands in the country.
In addition, since 2015, the Chilean coast has been affected by intense, recurrent storm surges, which have caused drastic morphological changes to numerous sandy shores and severe damage to coastal infrastructure [107,108]. Recent research has established a direct relationship between warm ENSO phases and the recurrence and intensity of storm surges, and that amid the scenario of an 8.5 emissions climate change, they will continue to affect the Chilean coast [107]. This has resulted in 80% of beaches in Chile undergoing a significant erosion process, possibly triggered by storm surges but certainly also by anthropogenic factors [109].

The Chilean Coast in Crisis: Guidelines for Integrated Management
Given the tectonic context of the Chilean coast and the predicted climate change scenario, having guidelines that allow decisions to be grounded in territorial planning and integrated coastal zone management is a priority. Unfortunately, similar to many countries in other parts of the world, Chile does not have concrete instruments that allow land uses to be planned according to the evolutionary dynamics of the coast and climate change scenarios to be included. The most recent studies on vulnerability to climate change on the Chilean coast establish that the increase in sea level will be important in the coming 15 years and that extreme events associated with storm surges will continue to recur over the next 45 years (IPCC 8.5 emissions scenario) [110]. In addition, experiences with reconstruction processes in the last 10 years show a high social and economic cost, with at least 13 reconstruction processes having taken place after tsunamigenic earthquakes, mudslides, and volcanic eruptions that have affected the coast. In particular, research on the reconstruction process following the earthquake and tsunami of 27 February 2010 has established that it has promoted the development of new risk areas instead of an effective risk reduction, considering the great opportunity such activities present to increase the social and urban resistance of the affected areas [59,111,112].
In this context, and amid the dramatic reduction in marine coastal ecosystems in the country, particularly due to urbanization and urban growth [113,114], a coast law in Chile is being pushed for [100]. Among the criteria proposed in this law are disaster risk reduction and integrated coastal management focused on human welfare sustainability (HWS) and ecosystem-based management. The consideration of the precautionary principle, environmental justice, and the recognition of the public nature of the coast are also important [115]. These efforts, while medium or long term, will allow the foundations to be laid for future governance focused on sustainable development.
Given the strong tectonic control of the Chilean coast and the urgent need to adapt to climate change, the main challenges to increasing coastal resilience and sustainability are closely linked to the geomorphic adaptation of estuaries to new streamflow conditions and forcings and the dynamic response of intertidal vegetation communities. Indeed, resilient ecosystems require a rigorous, science-informed coastal planning approach, implemented at the appropriate time scale [27]. To achieve this, the decision-making process for estuaries requires extensive interdisciplinary studies such as that presented in this paper. Only by analyzing the implications for hydrodynamics and morphodynamics, as well as ecological, biogeochemical, and dependent processes, can uncertainties be reduced, enhancing the predictive ability of models to provide management advice [28] and improving the quality of the measures that are implemented. The effectiveness of these measures relies on the scientific effort to understand each system and on the political will to translate scientific language into the language of public policies and methodological decision-making processes [116].

Conclusions
The Tubul-Raqui wetland presents an erosive trend of −1.32 m/year , facilitated by tectonic factors (interseismic subsidence). This process could be increased by anthropogenic factors capable of producing alterations in the sedimentary balance of the estuary.
The recovery of the salt wedge was fast, even considering a large-magnitude coseismic uplift (1.4 m). A year after the event, the saline intrusion was recorded up to the first 3 km, and two years after, it reached 10 km, returning to pre-earthquake conditions. The foregoing statement establishes the high resilience of these environments to violent geomorphological changes. The effects of the 2010 earthquake and tsunami, despite the ecological damage, especially to benthic macrofauna, did not significantly affect the floristic composition and vegetation structure of the wetland.
The Tubul Beach, located east of the estuary mouth, currently presents erosive behavior and induces vulnerability in the system since its sediment dynamics directly affect the hydrodynamics of the Tubul-Raqui wetland. The sediment transport trends have the potential to block the flow exchange at the mouth of the estuary, with consequences for the physical and biological environment and social implications for the fishing community. Moreover, human interventions in the area such as the construction of a new dock on the western beach could further affect the littoral dynamics of the area.
The increasing urbanization that coastal areas have undergone in recent decades and current trends in climate change highlight the need for a better interdisciplinary understanding of coastal areas to achieve more sustainable and resilient natural systems [58]. It is a priority to know how the coastline changes and recovers after a high-impact natural disturbance to contribute to a paradigm shift toward more integrated and holistic coastal zone management.

Data Availability Statement:
The data presented in this study may be available in a redacted form on request from the corresponding author.