Probabilistic characterisation of coastal storm-induced risks using Bayesian Networks

Coastal areas are often affected by inundation and erosion storm-induced risks. Detailed local risk assessments usually propagate a source (storm) through a pathway (coastal morphology) to characterise hazards (i.e. erosion and inundation) at the receptors and assess corresponding consequences. A probabilistic estimation of hazards based on the coastal response requires assessing large amounts of source characteristics. In addition, the coast is a dynamic environment, 10 and factors such as climate change projections or existing background erosion trends require performing risk analyses under different scenarios. This work applies Bayesian Networks (BNs) following the source-pathway-receptor-consequences scheme aiming to perform a probabilistic risk characterisation at the Tordera Delta (NE Spain). The BNs allow an efficient assessment of results from a large number of storms (179) and their simulated consequences at the receptor scale (~4000 receptors). Presented results highlight the storm characteristics with higher probabilities to induce given risk levels for 15 inundation and erosion, and how these are expected to change under given scenarios of shoreline retreat due to background erosion. The BNs also output probabilistic distributions of the different risk levels conditioned to given distances to the beach inner limit, allowing for the definition of probabilistic setbacks.

To characterise the forcing, the present work used hindcast waves from the Downscaled Ocean Waves dataset (Camus et al.,80 2013) derived from the Global Ocean Waves (Reguero et al., 2012). Hindcast surge from the Global Ocean Surge dataset (Cid et al. 2014), obtained at 4 locations close to the Tordera Delta at ~20 m depth, covering the period from 1954-2014 ( Figure 1), was also used. The simultaneous astronomical tide was added to the Global Ocean Sampling (GOS) dataset to obtain the total water level. The astronomical tidal range in the study area was about 0.25 m.

General framework
The methodology used in this work adapts the general approach of Jäger et al. (2018) where BNs were applied to implement the SPRC framework to assess storm-induced coastal risks. This approach has been previously implemented by Sanuy et al. 95 (2018) at the Tordera Delta to compare, in a deterministic manner, different risk reduction measures. In this work, the scheme was upgraded to a fully probabilistic risk characterisation and consisted of the following steps: (i) Storm characterisation. This step consisted of defining the local storm climate from long-term wave time-series. This stage corresponded to the (probabilistic) characterisation of the source. In practice, the result of this step was a storm dataset containing the hourly evolution of wave parameters during each event for a long period (multiple decades). 100 (ii) Hazard assessment. Once the forcing was characterised, the next step was the assessment of the storm-induced hazards, i.e. erosion and inundation, which were simulated using a process-based morphodynamic model, XBeach.
This stage corresponded to the characterisation of the pathway. To ensure the probabilistic representation of the https://doi.org/10.5194/nhess-2020-99 Preprint. Discussion started: 11 May 2020 c Author(s) 2020. CC BY 4.0 License. hazards, this step was performed for all the events of the storm dataset (first step) or for a subset of events that ensures an equivalent representation of the multivariant population representing the source. 105 (iii) Risk characterisation. In this step, simulated storm-induced hazards across the study area were transformed into risk values at the scale of individual receptors (existing buildings and infrastructure). To this end, vulnerability rules were defined as a function of the receptor typology and analysed hazard. In this stage, the receptor and consequence phases of the SPRC framework were tackled.
(iv) Scenario definition. This step consisted of defining the conditions for the assessment in terms of climate and 110 geomorphological scenarios of interest. This might require repeating steps (ii) and (iii), for all identified storms in (i).
Here, the entire storm dataset was used to characterise the baseline scenario (current conditions), while the additional scenarios were assessed with a representative subset to reduce the computational effort. The subset was also used to assess the baseline scenario to later verify that it statistically represents the same population as the original dataset from the perspective of the obtained results (for the validation of the method). 115 (v) BN integration. The obtained results for each event at the receptor scale were related to the variables characterising the storms (e.g. bulk features) and receptor properties (e.g. location) and integrated within the BN. Therefore, the BN outputs risk probability distributions accounted for the variability in the forcing conditions as well as the spatial distribution of receptors.

120
In the following sections, specific methods used in each step to analyse the Tordera delta case study are presented. Although some specificities are included, adopted methods are general enough to be applicable at nearly any site.

Storm characterisation
Coastal storms have been identified from wave time-series by employing the peak-over-threshold (POT) method using a 125 double threshold criterion as in Sanuy et al. (2020a). The first threshold, the 0.98 quantile (Hs = 2 m, in agreement with Mendoza et al. 2011 for NW Mediterranean conditions), is used to identify storm start and end times, and thus, controls the event duration and inter-event fair-weather periods. The second threshold, the 0.995 quantile (Hs = 2.6 m), is used to filter events that do not reach this value at the peak and would not be significant in terms of induced impacts.
The obtained dataset is composed of 179 storms (~3 storms per year), each being characterised by the hourly evolution of 130 wave conditions (significant wave height, Hs; peak period, Tp; storm surge; wave direction; and directional spreading). Of the 179 events, 43 correspond to multi-peak storms. These events occur when fair-weather conditions (Hs below the first threshold) between consecutive peaks last less than 72 hours ( Figure 2); they are relatively frequent in this part of the NW Mediterranean . In 12 cases, storms are formed by 3 or more peak sequences, leading to a total number of 237 individual storm peaks. 135 https://doi.org/10.5194/nhess-2020-99 Preprint. Discussion started: 11 May 2020 c Author(s) 2020. CC BY 4.0 License. To reduce the computational effort when assessing multiple scenarios, a storm subset is built aiming to maintain the statistical representativeness while avoiding the repetition of simulations of strongly similar storm conditions. The procedure consists of grouping the main variables defining the storm (Hs, Tp, duration, and direction) in homogeneous intervals 140 covering the entire range of local conditions (see Table 1). Each storm from the dataset falls into one of the resulting 4 × 4 × 3 × 3 = 144 combinations of bulk characteristics. Some combinations are populated with several storms, while others are empty groups, i.e. storm characteristics that have not been recorded and, therefore, not present in the storm dataset. To produce the subset, one storm is selected for all combinations populated with at least one event. To ensure a probabilistic representation of the source, the number of storms belonging to each combination is counted for later use as a weight 145 (multiplicity factor) when feeding the BN with results from that event. As was previously mentioned, one of the local characteristics of the storm climate in the study area is the presence of multipeak storms. As the impact of successive storms separated by relatively short fair-weather periods may be different to that of single events depending on storm characteristics and initial beach configuration (e.g. Dissanayake et al. 2015;Eichentopf et al. 2020), we retained these storms in the analysis. Thus, to properly account for their potential effects, all existing identified https://doi.org/10.5194/nhess-2020-99 Preprint. Discussion started: 11 May 2020 c Author(s) 2020. CC BY 4.0 License.
multi-peak storms in the original time-series were included in the subset. To this end, their impact was simulated with the 155 XBeach model saving the cumulative output after each peak. This allowed treating each simulation of a multi-peak event as a representative of two storms: (i) the multi-peak event itself and (ii) a single-peak storm with properties matching the first recorded peak. Thus, the storm subset comprised of 69 storms, including 43 multi-peak storm events (see Table 1). Note that the intervals used to classify storm variables are more refined than in the BN bins (Section 6), to later ensure intra-bin variability during the training. 160 The statistical representativeness of the subset with respect to the full storm dataset was tested using the methodology to compare histograms proposed by Bityukov et al. (2013). This method assumes that values at each bin of the histogram follow a normal distribution with expected value ni,k and variance σ 2 i,k (with "i" representing the bin and "k" the histogram). Thus, the significance is defined as: where ̂, is an observed value at bin "i" of histogram "k" and ̂, = ̂, . Therefore, we consider the root mean square (RMS) of the distribution of significances as: where ̅ is the mean value of ̂ and M is the number of bins of the histogram. The RMS represents a distance measure with the following interpretation: If RMS = 0, both histograms are identical; if RMS =0~1 both histograms are obtained from the same parent population; if RMS >> 1, histograms are obtained from different parent distributions. The method is applied to compare the output distributions resulting from training the BN with the whole dataset vs. training it with the subset. 175

Hazards assessment
Storm-induced hazards (erosion and flooding) have been modelled using the XBeach model (Roelvink et al. 2009), which has been previously calibrated for the Tordera Delta (see Sanuy et al. 2019b). The calibration of the model achieved a Brier Skill Score (BSS) (Sutherland et al. 2014) of 0.68. The model was implemented using a curvilinear grid with a variable cell size around the Tordera River mouth (Figure 1). The extension of the mesh is approximately 1.5 km in the cross-shore 180 direction, with a cell size ranging from 5-6 m at the offshore boundary (20 m depth) to 0.7-0.8 m at the swash zone. In the alongshore direction, the model has an extension of 4.5 km with cell size ranging from 25 m at the lateral boundaries 2-3 m around the river mouth. Storm input consists of time-series of wave conditions characterising each storm obtained from the DOW dataset at the 4 nodes at the offshore boundary (Figure 1), with a time-step of 1 hour, which is the time resolution of the original data. The model was used to simulate storm-induced hazards under 455 different events, which correspond to 185 179 original storms, plus a subset composed of 69 storms under 4 different scenarios.

Risks
To assess the induced risk, first, receptors in the study area are individually considered by their footprint polygons (~4000) and delineated using a Geographic Information Systems (GIS)-based tool to account for their exact position and dimension.
Once they are defined, a direct correspondence between each receptor with the underlying XBeach model mesh is available 190 in such a way that each receptor is associated with the model nodes directly affecting it (see Figures 3 and 4). To spatially characterise the risk of the area as a function of the variability of the local geomorphology and coastline orientation at both sides of the river mouth, five different sectors along the coast were defined (Figure 1). Two of them, SBN and SBM, are https://doi.org/10.5194/nhess-2020-99 Preprint. Discussion started: 11 May 2020 c Author(s) 2020. CC BY 4.0 License. located northwards of the river mouth ( Figure 1), with SBM being limited to the south by the river mouth. The main distinctive feature of SBN is the existence of a promenade limiting the inner part of the beach. Southwards of the river, there 195 are three sectors (Figure 1): MSM being closest to the mouth; MS1, which is located southwards of a coastal revetment; and MS2 located furthest to the south, with wider beaches, and sheltered against Eastern storm waves, which are dominant in the area ).
The vulnerability of each receptor is individually characterised as a function of their structural properties. Receptors in the study area comprise hard constructions, such as houses and infrastructures, and softer elements such as campsite elements 200 (e.g. bungalows) (Sanuy et al., 2018). To assess the flooding-induced risk, the relative damage to receptors is calculated using flood-damage curves (Table 2) using the maximum-modelled water depth within the receptor polygon. No specific damage curves are available for the Catalan coast, and due to this, we used the curves recommended and used by the Catalan Water Agency (ACA, 2014) for the development of inundation management plans. Risk to life was also been included in the assessment by using the water-depth-velocity product as input (Table 3, Priest et al., 2007) within the receptor's boundaries. 205 For the erosion hazard, the magnitude of the associated risk is based on the distance from the significantly eroded XBeach nodes to the receptors. Significant erosion was set to 0.25 m of the vertical bed level change and assumed as the common minimum depth for light structure foundations. The closest distance from the receptor corners to that erosion level was compared with the erosion risk thresholds indicated by Jiménez et al. (2018) (Table 4).
Therefore, the result of each simulation (hazard maps) was transformed into a risk value at the individual receptor. Figure 3  210 shows an example of simulated inundation water depth for a long return period event and its transformation into relative inundation damages to receptors: None (0%), Low (0-30%), Moderate (30-60%), High (60-90%), and Extreme (>90%). Figure 4 shows, for the same event, results corresponding to the erosion hazard. Individual results were stored at each of the ~4000 receptors for each of the simulated events, leading to a total number of 716,000 and 276,000 cases to feed the BN with the entire dataset and with the subset, respectively. 215

Scenario definition 230
When assessing risks in coastal areas under changing conditions, it is necessary to consider these potential variations in the assessment, otherwise, its utility for medium-long term risk management will be limited. As previously mentioned, the study area is a highly dynamic sedimentary environment subjected to a background coastal retreat (Jiménez et al. 2018). Thus, in this step, different scenarios characterising future configurations were built based on the expected future coastal changes. This was accomplished by using decadal-scale background erosion rates estimated for the different beach sectors by Jiménez 235 and Valdemoro (2019) by analysing shoreline changes from aerial photographs. The estimated average shoreline retreat at each sector is 1.1, 4.0, and 1.9 m/y at SBN and SBM, MSM and MS1, and MS2, respectively (see Figure 1 for locations). It is assumed that current evolution trends remain constant during the timeframe of the analysis, which is limited to 20 years.
However, this could be substituted by time-varying evolution rates provided this should be the case.

Figure 5: Changes in the bed level grid for the future scenarios. Difference between baseline bed level and scenario bed level (upper). Profile retreat at both sides of the river mouth at the different time horizons (lower).
Thus, to account for this background response, each scenario was defined based on a given coastal morphology at a given time horizon. The baseline morphology, which corresponds to the current scenario, is the one described in Section 2 ( Figure  245 1) that was directly measured. Future coastal morphology for each scenario corresponding to different time horizons (+5 years; +10 years; and +20 years) were built by retreating the active part of the shoreface, from a -10 m-depth to the subaerial beach, according to erosion rates at the different areas. To ensure alongshore smoothness after retreating, linear transitions between sectors affected by different retreat rates were applied. Resulting configurations for two scenarios are shown in Figure 5, along with example profiles at locations under different levels of background retreat. Local constraints due to the 250 lack of accommodation space due to the existence of hard structures at the hinterland were also considered. As an example, Figure 5 shows the beach profile retreat at two locations with different hinterland characteristics: P1 has no hard limit, whereas P2 is limited at the back by a promenade. This results in a continuous retreat of P1 for all scenarios, whereas the retreat of P2 is limited at the beach after 10 years.

Bayesian Network integration 255
The BNs are probabilistic models based on acyclic graph theory and Bayes theorem (Pearl, 1988;Jensen, 1996). They have demonstrated their versatility and utility in efficiently combining multiple variables to predict system behaviour. Within the context of this work, they can be used to represent the SPRC scheme through the dependency relations between the different steps (see e.g. Straub 2005;Jäger et al. 2018). In this sense, they can easily be adapted to assess different natural hazards and their impacts on many kinds of receptors, for both descriptive as well as predictive applications (see e.g. Beuzen et al. 260 2018b).
In this work, two BN configurations were used to characterise the system response to the impact of coastal storm events.
This was done to optimise the BN structure by limiting the number of variables per network while solving the different parts of the SPRC framework. In practice, one BN solved the source-consequences relationships (BN-A), while the other https://doi.org/10.5194/nhess-2020-99 Preprint. Discussion started: 11 May 2020 c Author(s) 2020. CC BY 4.0 License. characterised the receptor-consequence spatial distribution (BN-B), providing complementary information on the local risk 265 profile. BN-A ( Figure 6) links storm-defining variables (Hs, Tp, duration, direction, and water level) and impacts to the receptors (erosion impact, risk to life, and structural relative damage). The central variable of the network (indicated by * in Figure 6) was used to perform conditioned assessments. Depending on the objective of the analysis, it can be (i) the total number of affected receptors by inundation within a storm event; (ii) total number of affected receptors by erosion within a storm event, 275 or (iii) receptor area (SBN, SBM, MSM, MS1, and MS2), as shown in Figure 6. To account for the spatial extension of the impacts, we included the total number of affected receptors as an output variable. To characterise the impact of inundation, for each storm, all receptors presenting a relative damage other than 0%, or a risk to life other than "None" were counted.
Similarly, to characterise the impact of erosion, all receptors presenting an impact level different than "None" were counted.
In practical terms, this means that, in general, the number of affected receptors by erosion was larger than by inundation. 280 This is because, with the used criteria, it is quite probable to have receptors affected by "Very Low" to "Moderate" erosion risks representing the loss of protection provided by the beach, although this does not imply that they will be directly exposed to wave impact. However, inundation-related impacts are always associated with the presence of water at the receptors. This has to be taken into consideration when interpreting the obtained results.
https://doi.org/10.5194/nhess-2020-99 Preprint. Discussion started: 11 May 2020 c Author(s) 2020. CC BY 4.0 License. BN-B (Figure 7) links the simulated impacts on the receptors to their position, characterised by their location along the coast (area) and the distance to the public domain (DPMT), which is the limit between the beach and the hinterland. Additionally, two variables accounting for the number of impacts per receptor were included. These variables provided additional insight 290 into the system response, as the obtained distributions with the BN merge storm-climate variability and the spatial distribution of receptors. For the inundation risk, the number of impacts with damage different from 0%, and/or with risk to life different from "None" was counted at each receptor. For the erosion risk, the number of impacts different from "None" was counted per receptor. This has the same consequence as that described in the previous case (BN-A) for interpreting the obtained results. It must be noted that from all receptors displayed in Figures 1d, 3, and 4, only those presenting at least one 295 impact for the entire storm dataset, by either inundation or erosion, were used for the BN training. Otherwise, the choice of receptor population to include in the assessment would be arbitrary, affecting the obtained distributions.
The presented BN-model was designed to assess storm-induced risks in a coastal hotspot where the storm climate and coastal response are well known (e.g. Jiménez et al. 2018;Sanuy et al. 2020). Due to this, the discretisation of variables (Figures 6 300 and 7) was done manually, enabling better accuracy than automatic unsupervised methods and closer accuracy to supervised discretisation with less associated variability on model performance (Beuzen et al. 2018a). original dataset using the baseline morphology to ensure that it properly represents the local storm climate. This is followed by the presentation of the risk characterisation of the Tordera Delta, starting with risk probabilities integrating all storms and receptors (global risk probabilities), and then, with conditioned probabilities between forcing-area risk (BN-A, Figure 6) and 310 area-distance risk (BN-B, Figure 7). Table 5 shows the obtained statistics using Eq. 1 and 2 to compare the discrete probability distributions obtained with the BN using the 179-storm dataset against those from the 69-storm subset. This is done for the different BN-outputs, i.e. global risk probabilities, which are the impact distributions in Figures 6 and 7; probabilities of storm characteristics (distributions of Hs, 315 duration, direction, and water level) conditioned to different risk levels and areas; and risk probabilities conditioned to receptors locations (area and distance to the beach limit). This involves the comparison of more than one variable output (e.g. impact results are always three variables), and therefore, results are given as a mean and standard deviation.

Subset validation
All obtained values of the mean significance ̅ and its root mean square (RMS) are close to 0; therefore, from the perspective of obtained results, it can be assumed that the obtained distributions by feeding the BNs with the subset almost identically 320 represent the same source population as that of the complete dataset. This is true both for global distributions and for conditioned discrete probability density functions (PDFs).  Storm characteristics conditioned to risk levels Impact probabilities at different areas, and conditioned to Hs, duration, water level, and direction (e.g. Fig 8 and 9) 0.0006±0.02 0.05±0.03

Risk probabilities conditioned to receptors locations
Impact probabilities at the different areas and distance to the beach (Figures 10 to 12) 0.0041±0.02 0.04±0.08 325 Table 6 shows the obtained probability levels for different tested scenarios in the study area. These so-called prior (unconstrained) probabilities represent the expected frequency of the different risk levels in the study area and account for the variability of the source (storm climate), spatial distribution, and extent of the impacts on the receptors. In general, under current conditions, the probability of receptors being affected by significant (high and extreme) risks is low (1-2%). 330

Risk characterisation
However, the existence of background erosion in the study area results in a significant increase in future risks. Under the baseline scenario, the computed probability of moderate-high risks associated with erosion is larger than the ones for inundation. However, when we only consider those cases where erosion results in exposing receptors to direct impact (high and extreme risk), the obtained probability values are of the same order of magnitude as those obtained for moderate damages and risks associated with inundation. Additionally, results of number of affecter receptors from BN-A (not shown 335 in the table) show an increase in the % of storm conditions affecting a large number of receptors along the study area. As an example, storm conditions with the potential to affect more than 200 receptors with any level of inundation risk increases https://doi.org/10.5194/nhess-2020-99 Preprint. Discussion started: 11 May 2020 c Author(s) 2020. CC BY 4.0 License. from 4% under current conditions to 20% and 40% after 10 and 20 years, respectively. Simultaneously, storm conditions affecting more than 450 receptors with any level of erosion risk will rocket from the current 4% to 100% in 10 years. Here, it is important to remember that erosion risk is not only related to direct impact but also the loss of protection function 340 (decrease of beach width in front a of a given receptor), while inundation risk implies the direct effect of water on the receptor. In general, estimated probabilities associated with erosion-induced risks are larger than those due to inundation when comparing similar risk levels. are located northwards of the river mouth, whereas areas southwards of the river mouth are more affected by inundation (higher probability values). The time evolution of the affected receptors is also different, reflecting existing spatial variations in shoreline retreat rates. Thus, the largest relative increase in the number of impacted receptors under future scenarios 350 occurs southwards of the river mouth. Notably, the MS2 sector is the most sensitive to future risks, as currently, although it is well protected by a relatively wide beach, this protection will fade after 10 to 20 years.
BN-A was also used to characterise the conditioned probabilities of storm characteristics associated with the highest risks and assess whether these probabilities vary along the study area. As seen in Figure 9, under current conditions, the main 355 storms driving the highest inundation-induced risks are characterised by Hs higher than 4 m and from the E direction. This is valid for the entire area, although their relevance slightly varies along the coast. Thus, the only exception is found in the SBN sector, where the promenade is so close to the shoreline that lower Hs can induce inundation damages. For future conditions (20 years scenario), the relative importance of storms with smaller Hs increases, and the relative importance of present secondary wave directions, S and SE, also increases in relative terms. 360

375
The spatial distribution of the expected impacts across the study area was analysed using the BN-B. The objective of the analysis was to assess the probability damage occurring at the receptors located at a given distance from the beach (i.e. limit between beach and hinterland). Figures 10 and 11 show obtained results in terms of % of inundation-induced damage and risk to life, respectively, for different time horizons. Consistent with the results shown in Table 6, under current conditions 380 (baseline), storms cannot induce extreme structural damage (>90%) (Figure 10) nor extreme risk to life (Figure 11). High damages (> 60%) are mainly concentrated at the outer fringe of the hinterland of the two locations (MSM and MS1) with associated conditioned probabilities of 21% and 5%, respectively. These two areas also show the highest probabilities of risk penetration into the hinterland. Northwards of the river mouth, the SBN sector presents a large probability of moderate damages, but it is limited to the external fringe. Regarding risk to life, a similar spatial pattern is observed, with MSM 385 showing the largest probability of high risk (20%) at the external fringe, SBN at the north with 12%, and MS1 only showing a residual 3%. The obtained results reflect the role played by the current coastal morphology, where the southern area is characterised by narrow and low elevation beaches (MSM and MS1), whereas the SBN sector in the north is composed of a narrow beach backed by a promenade. Notably, SBM with a narrow beach but higher topography without a promenade and MS2 with low topography but wider beaches are the areas presenting the lowest risks. 390 Under future conditions (+20-year scenario), significant changes are observed in the intensity of risks and extension across the territory (Figures 10 and 11). The spatial modulation on induced risks as a consequence of the beach narrowing due to background erosion is especially evident in the southernmost area, MS2. Whereas this sector does not experience any risk under current conditions, significant probabilities of moderate and high damage and risk to life is expected to occur in 20 395 years, not only at the outer fringe but also in inner positions of the hinterland. The other sectors along the coast also show significant increases in the probability of occurrence of any type of risk and extension of the impacts landward (Figures 10   and 11). https://doi.org/10.5194/nhess-2020-99 Preprint. Discussion started: 11 May 2020 c Author(s) 2020. CC BY 4.0 License.  The spatial distribution of erosion-induced risk under current conditions (Figure 12) reflects the existence of hard elements and varying beach widths along the study area. Thus, SBN presents the largest probability of extreme risks at the promenade (54%), followed by MSM (9%), and MS1 (3%). In SBN, the promenade acts as a physical boundary for erosion; the 410 distribution of risk levels into the hinterland shows a linear pattern reflecting its position. At the southernmost end, MS2 is https://doi.org/10.5194/nhess-2020-99 Preprint. Discussion started: 11 May 2020 c Author(s) 2020. CC BY 4.0 License. currently well protected by a wide beach and no risk is predicted under the current conditions. Under the +20-year scenario, the effect of the promenade in SBN is reflected through the unaltered spatial pattern of affected locations and computed probabilities. MSM and SBM show the largest relative increase of extreme risks (probabilities of 45% and 47%, respectively) at receptors located closest to the beach, along with the largest spatial propagation of risks into the hinterland, 415 as no hard elements are present to limit the retreat of the shoreline. At MS1, the probability of extreme risks increases to 18% at the beach limit with small changes at larger distances, while MS2 starts presenting significant probabilities of low risks indicating that the beach will begin to decrease its protective function against storm impacts after 20 years.

Discussion
In contrast to previous applications of the BN-SPRC concept presented in Jagër et al. 2018(e.g. Van Verseveld et al., 2015420 Plomaritis et al., 2018;Ferreira et al. 2019;Sanuy et al., 2018), this paper presents a fully probabilistic characterisation of the source using all available storms in a 60-year long wave time-series hindcast, modelling their induced erosion and inundation risks over all the identified receptors at the study site.
The methodology was successful in identifying storm characteristics with higher probabilities to induce given risk levels for different coastal hazards (inundation and erosion). It was efficient in assessing the expected changes in storm characteristics 425 and probabilities under different scenarios, which were developed based on the background mid-term coastal evolution. In this sense, the obtained relation under current conditions of erosion and inundation risks with storm direction and Hs depicts the general characteristics of storm-induced hazards in the study area ). The BN output showed a lack of correlation between high risks and water levels, consistent with the previous findings of Mendoza and Jiménez (2008) on the non-relevance of storm surges. Under future conditions, the background shoreline erosion changes the sensitivity of the 430 area to storms. Thus, for the tested scenarios, the population of storms with potential to significantly impact the area increases and higher risks will be associated with storms characterised by lower Hs with currently secondary wave directions ( Figure 8). If we combine this larger exposure to southern storms with the large sensitivity of the area to the impact of such S storms (Sanuy and Jimenez, 2019), this may have serious implications for the future risk management of the area.
The method has been designed to provide a detailed spatial assessment to assess the sensitivity of the area, which permits the 435 association of the local risk profile with different morphological characteristics such as beach orientations, height, and the presence of hard structures. In this sense, the local response affected by the presence of the promenade at S'Abanell (SBN), and revetment in Malgrat North (MS1) were adequately characterised by the BN. This spatial analysis also permitted the assessment of a differentiated variation in future risks along the study area. Thus, whereas some areas being currently exposed linearly increased the probabilities of higher risks, other areas currently well protected will be subjected to higher 440 future risks without any variation in storminess.
The method can also be used for testing risk management measures such as the performance assessment of different setbacks. While this measure is effective in reducing coastal damages in eroding coastlines, especially in the context of climate change (Sanó et al. 2011), it has to be defined for given time horizons and driving conditions (e.g. Wainwright et al. 2014). To this end, the framework presented herein permits the definition of probabilistic setbacks at the study site. 445 Moreover, as this definition is based on the probabilistic distributions of the different risk levels and impacts per receptor at different locations across the coastal domain, it differs from existing approaches that are essentially based on the probabilistic definition of the shoreline position (e.g. Jongejan et al. 2016). As an example, Table 7 shows the calculated minimum distances landward of the inner limit of the beach according to different risk levels for different time horizons (scenarios). As the BN output combines the natural variability of the storm climate with the spatial variability of the 450 impacted receptors, setbacks can be defined from these (total probability, as in Figures 10 to 12) or by assuming that the presence of a given risk level must be completely tackled, focusing then only on the spatial distribution of receptors under https://doi.org/10.5194/nhess-2020-99 Preprint. Discussion started: 11 May 2020 c Author(s) 2020. CC BY 4.0 License. such levels. The second approach will result in more conservative (wider) buffers. Table 7 shows the calculated buffer distances using both perspectives. The setbacks obtained that account for the total probabilities can be used as proposals for managed retreats, as they reflect the areas with a high number of impacts per receptor; the setbacks defined by the presence 455 of a given risk level can be used to inform self-preparedness against risk, as they highlight zones where the existence of risk is possible but highly infrequent. It must be noted that all scenarios have been simulated without any assumption of receptor re-allocation, and therefore, hard limits for erosion remain homogeneous across scenarios. Therefore, the distances presented in Table 7 must be interpreted as the evolution of the baseline setbacks at different horizons in a business-as-usual situation. The presented method is based on the response approach (Garrity et al., 2006;Sanuy et al., 2020a) as it produces 465 probabilities based on how hazards (erosion and inundation) affect the receptors in each of the storm events derived from a long dataset of 60 years; it does not allow the extrapolation of the storm conditions out of the range of the ones registered in such datasets. This has relatively less impacts on the results when compared to the impacts from other sources of uncertainty, such as morphological variability or model error (Sanuy et al., 2020b). Nonetheless, it allows the simulation of all storm events with their real shapes (time evolution of storm characteristics) without introducing large uncertainty in hazard 470 estimation associated with the use of synthetic storms that are commonly used to define the shape of statistically extrapolated storm events (see e.g. Duo et al., n.d.).

460
In this study, hazards were computed using a robust model to simulate the storm-induced coastal response, XBeach, calibrated for an event representative of extreme conditions (see . They were converted to risk by using damage curves recommended for use in the study area. However, the BN methodology is flexible for any kind of 475 model, as well as to include model uncertainties (using different models or setups) and measurements (e.g. Sanuy et al., https://doi.org/10.5194/nhess-2020-99 Preprint. Discussion started: 11 May 2020 c Author(s) 2020. CC BY 4.0 License. 2020b for cross-shore parametric models) to extend the data training and improve the results while testing its predictive capacity.
With regard to building future scenarios to assess future risks, we have limited the present study to mid-term scenarios, i.e. at the decadal scale (20 years). They were built based on decadal-scale shoreline rates of displacement retreat measured by 480 Jiménez and Valdemoro (2019), which were used to build future coastal configuration assuming that no changes in evolutive conditions will occur. Even in this case where no changes in forcing conditions were applied (no changes in storm conditions nor sea level rise), this approach permitted the identification of significant changes in the storm-induced risk profile. This could be extended or adapted to changing future conditions using midterm morphological simulations with varying climatic forcing or under different adaptation scenarios, and then, used as an efficient way to test risk management strategies. 485

Summary and Conclusions
Bayesian networks have proven to be an efficient tool to develop an SPRC-based framework for probabilistic storm-induced risk assessment and risk mapping at a local scale (few kilometres). In this work, BN training has been carried out using storm events identified in a 60-year long wave time-series, and simulated hazards and corresponding risks were evaluated at the receptor scale (few metres). This resulted in a fully probabilistic characterisation of risks that accounted for climate 490 (storms) and geographic (receptor location) related variabilities. The framework is also able to predict how risks will evolve in the near future, both in intensity and spatial distribution, provided that climate and/or geomorphology scenarios are built.
One of the advantages of the framework is that it permits the identification of conditional probabilities, and thus, the identification of which are the storm characteristics that induce risks of a given magnitude. This is a very useful property in designing disaster risk reduction (DRR) strategies and measures including the design of early warning systems. 495 Concerning the analysed case study, the Tordera Delta (NW Mediterranean coast) presents, under current conditions, a larger susceptibility to storm-induced erosion than to inundation, which was identified through computed probabilities of high-risk levels associated to both hazards along the coast. Storms inducing the largest impacts are characterised by high Hs (>4 m) for inundation and long duration (>60 hours) for erosion. In both cases, these correspond to Eastern events, which are the most energetic in the area. 500 The application of the framework for future scenarios predicted an increase in the local risk as a larger number of storms will be able to induce higher risk levels. As these scenarios were built by projecting the coastal configuration up to two decades from now (based on background erosion), the framework reflected the morphodynamic feedback resulting from the loss of protection provided by progressively narrowing beaches. In addition to the increase in risk levels, it also identified a change in storm threshold conditions affecting the area in a significant manner, characterised by lower Hs values and with an 505 increasing importance of southern events.
Finally, the obtained spatial distribution of risks permitted the identification of the most sensitive areas and their evolution over time. This can be used to make decisions on the required DRR measures both along the coast and across the hinterland.
The use of the BN to obtain probability distributions of the different risk levels across the hinterland allowed for a probabilistic definition of setbacks.