A Localized Particle Filtering Approach to Advance Flood Frequency Estimation at Large Scale Using Satellite Synthetic Aperture Radar Image Collection and Hydrodynamic Modelling

Abstract

This study describes a method that combines synthetic aperture radar (SAR) data with shallow-water modeling to estimate flood hazards at a local level. The method uses particle filtering to integrate flood probability maps derived from SAR imagery with simulated flood maps for various flood return periods within specific river sub-catchments. We tested this method in a section of the Severn River basin in the UK. Our research involves 11 SAR flood observations from ENVISAT ASAR images, an ensemble of 15 particles representing various pre-computed flood scenarios, and 4 masks of spatial units corresponding to different river segmentations. Empirical results yield maps of maximum flood extent with associated return periods, reflecting the local characteristics of the river. The results are validated through a quantitative comparison approach, demonstrating that our method improves the accuracy of flood extent and scenario estimation. This provides spatially distributed return periods in sub-catchments, making flood hazard monitoring effective at a local scale.

1. Introduction

One way to reduce flood risk is to develop tools for flood monitoring and prediction, as these tools help raise awareness about where and when inundation may occur and are crucial components of early warning systems [1,2,3]. Satellite remote sensing in addition to more traditional shallow-water modelling are widely used to generate flood hazard maps (i.e., here defined as flood extent with attributed probability of occurrence). These maps are essential for understanding potential flood risks in an area [4,5,6,7,8,9]. SAR sensors stand out for flood monitoring because they can capture Earth’s surface images anytime, regardless of day, night, or cloud cover. This unique ability makes them highly valuable for observing flood extent and develo** new flood monitoring methods [10,11,12,13,14,15]. The value of SAR data in map** flood extent is well-established, with numerous studies proposing various approaches [16,17,18,19,20].

Moreover, hydrodynamic modelling is a powerful tool to simulate flood events, by predicting how floodwaters will spread over time, providing valuable insights like flood extent and depth across an area [21,22,23,24,25,26]. Such models also allow for the computation of flood extent and water depth maps using synthetically generated input data for hypothetical, previously unobserved flood scenarios with varying probabilities of occurrence. In recent decades, incorporating remotely sensed flood data has significantly improved how we assess water conditions, calibrate flood models, and assimilate data for real-time flood monitoring and forecasting [8,27,28,29,30,31,32,33]. Recent research shows a surge in new methods for flood hazard assessment [22,34,35,36,37]. These methods combine statistical analysis with physical modeling to address uncertainties in data from various sources. Despite their strengths, both SAR data and shallow-water models have limitations for near-real-time (NRT) flood hazard assessment. While SAR images provide valuable snapshots of flood extent, they lack the ability to show how the flood is evolving or its severity in a historical context. Shallow-water models, while powerful for simulating floods, become computationally expensive to run in NRT for large areas. Furthermore, the alluvial phenomenon in a catchment takes on a different evolution depending on the local character of the hydro-geomorphologic dynamics [38,39,40,41,42]. To accurately assess flood risks, we need to consider how different parts of a river respond to floods. This means studying the variations in landforms and the natural processes happening along the river [43,44,45,46,47]. A promising way to circumvent model and satellite data limitations is to optimally combine them through data assimilation [48]. When applied with localization, data assimilation allows flood hazard maps to account for a river’s unique behavior during extreme events. Here, we focus on the need to estimate flood hazards by integrating model and satellite data at a local scale to account for the dynamism of flooding during extreme events.

This study aims to develop and test a new method for the NRT assessment of local flood hazards. Our approach uses satellite data to estimate the return period of a flood event as it occurs. This method is designed for NRT because it quickly incorporates data from a SAR image into a set of pre-calculated flood maps. This eliminates the need for slow, real-time simulations with shallow-water models. To evaluate our method, we applied it to a series of images captured by the Envisat satellite. These images show floods of varying severity (low to high magnitude) that occurred over the Severn River floodplains in the United Kingdom. We chose this case study because ground truth data, which allow us to verify the accuracy of our results, are available for this area. The remainder of the paper is structured as follows: Section 2 presents the methods; Section 3 describes the study area, available data, and experimental set up; Section 4 presents the validation approach; Section 5 presents the results and discussion; and Section 6 is the conclusion.

2. Methods

The proposed approach uses a particle filter [49,50] in order to optimally combine satellite and model data. It is a further development of the data assimilation algorithm first proposed in [48] where SAR-derived flood map are probabilistic maps [51] and simulated flood extent maps are binary maps. As explained and shown in [48,51], the SAR-derived probabilistic maps account for the flood extent classification uncertainty, providing a probability of flooding given the SAR backscatter value. Here, we therefore propose to follow the method proposed in [51] for SAR image processing and adapt the data assimilation approaches proposed in [48] to enable a particle filter localization and a local flood return period estimation (Figure 1). In particular, each SAR-derived flood map (SAR-derived PFM in Figure 1; see Section 2.1) is integrated into a series of simulated flood extent maps (DFM in Figure 1; see Section 2.2) using a particle filtering approach (blue circle in Figure 1; see Section 2.3). This process assigns a return period (estimated from the model) to the flood extent (observed from the satellite) in each hydrographic sub catchment (subregion map in Figure 1). The final map displays the flood extent (red box in Figure 1) along with a range of possible return periods, reflecting the dynamic nature of the flood event and considering the local fluvial behaviour.

2.1. Map** Flood Inundation Extent from SAR Images

Most SAR-based flood map** algorithms primarily focus on detecting water in areas with little to no vegetation or bare soil. Typically, algorithms for identifying flooding in bare soil areas rely on a single SAR image, identifying pixels with low backscatter values indicative of standing water. The primary factor creating a significant contrast between flooded and non-flooded areas in the radar image is the strong specular reflection from the water’s surface. The smooth water surface reflects the radar signal mostly in the specular direction, resulting in lower backscatter values. In contrast, rough, dry landscapes scatter the signal in multiple directions, leading to higher backscatter values. These principles frame flood map** as a classification problem with two categories: water and non-water [52].

SAR-based flood extent observations are prone to uncertainty that needs to be accurately estimated in the framework of data assimilation. In this context, we take advantage of a previously developed algorithm [51] that enables us to estimate, for each pixel of a SAR image, its probability corresponding to water given its backscatter. This unsupervised and automatic method returns probabilistic flood maps from large SAR images. For any pixel of the SAR image, its probability of being covered by water given its backscatter is derived using Bayes’ theorem, as follows:

$$\left(w|{\sigma }^0\right)=\frac{p\left({\sigma }^0|w\right)p\left(w\right)}{p\left({\sigma }^0\right)}$$

(1)

where $p\left(w|{\sigma }^0\right)$ is the probability of a pixel being covered by water given its backscatter value, $p\left({\sigma }^0|w\right)$ is the probability distribution of backscatter values of open water pixels, $p\left({\sigma }^0\right)$ is the marginal probability distribution of backscatter values of a pixel, and $p\left(w\right)$ is the prior probability of a pixel being water. Setting the latter equal to 0.5 by default (as proposed in [51]), one only needs to estimate $p\left({\sigma }^0|w\right)$ and $p\left({\sigma }^0\right)$ to derive the probabilistic flood map from Equation (1). To do so, the hierarchical split based approach of [52] is applied. First, the image is iteratively split into tiles of various dimensions to identify those exhibiting a bimodal histogram. Next, the parameters of a bi-Gaussian mixture probability density function are fitted to the histogram of the previously selected image sub-tiles. The first Gaussian distribution (i.e., with the lowest backscatter values) is used to estimate $p\left({\sigma }^0|w\right)$ and the sum of the two gaussian distributions is used to estimate $p\left({\sigma }^0\right)$. The output of this approach is therefore a probabilistic flood inundation map, called here a probabilistic flood map (PFM), in which each pixel exhibits a probability of being water and has a value ranging between 0 and 1.

Water is potentially not detectable by the proposed flood map** algorithm in some areas (i.e., in forested and urban areas); we considered these as exclusion areas [53] and we masked out the corresponding pixels during the assimilation using freely available land cover maps.

2.2. Generating an Archive of Design Flood Inundation Maps

We define design flood inundation maps as flood extent maps simulated by a shallow water model—namely Lisflood-FP [54,55]—using as boundary conditions design streamflow hydrographs. These correspond to synthetic hydrographs representing flood events of different return periods derived from observed streamflow time series. To do so, we first estimate the return period of each annual largest flood from the available discharge time series. This process involves deriving the annual maximum streamflow values from this time-series. Subsequently, these annual maximum streamflow values (x_max) are arranged in ascending order, and a generalized extreme value (GEV) distribution with three parameters—namely the location µ, scale σ, and shape κ [56]—is fitted to the estimated annual maxima using the maximum likelihood estimation (MLE) method:

$$F\left(x_{max};\mu ,\sigma ,k\right)=\left\{\begin{array}{c}\exp \left\{-{\left[1+k\left(\frac{x_{max}-\mu }{\sigma }\right)\right]}^{-1/k}\right\} if k\ne 0 \\ exp\left\{-exp\left[-\frac{x_{max}-\mu }{\sigma }\right]\right\} if k=0\end{array}\right.$$

(2)

where F is the GEV cumulative distribution function of x_max, the maximum annual streamflow at the considered gauging station. The GEV distribution includes Gumbel (k = 0), Frechet (k > 0), and Weibull (k < 0) families. According to the extreme value theory [57], this allows for the estimation of streamflow values associated with non-exceedance probabilities of flood events. The fitted GEV cumulative distribution function (CDF) provides estimates of the streamflow value Q associated to any non-exceedance probability P and, consequently, a return period T:

$$T=\frac{1}{1-P} where P=F\left(Q,\mu ,\sigma ,k\right)$$

(3)

For 15 preselected flood return periods (1, 2, 5, 10, 20, 30, 50, 75, 100, 200, 300, 400, 500, 750, 1000 years), design hydrographs are drawn using the corresponding design peak streamflow and a synthetic hydrograph shape. The decision to consider 15 flood scenarios is made to achieve a meaningful and representative set of designed flood extents, considering as many flood extents as possible and including very extreme events. Here, we propose to draw the synthetic hydrographs using the shape of an observed flood event. Based on this shape, design hydrographs are computed by rescaling the synthetic hydrograph so that its peak equals the design peak streamflow whereas its minimum remains unchanged (Figure 2).

To simulate the design flood inundation maps, we use the Lisflood-FP hydrodynamic model, a raster-based two-dimensional flood inundation model [48]. Topographic information and channel geometry, boundary conditions, and friction parameters are required for setting up and running the model. Water depth and water surface elevation are the main outputs of the model. The model is first used to simulate each design flood event utilizing the previously drawn design hydrographs as boundary conditions, assuming all river inflows experience the same flood return period at the same time. From each simulation, a maximum depth map is generated. It is subsequently converted into a maximum flood extent map assuming a grid cell to be flooded if the water depth is larger than 10 cm and dry otherwise. Finally, the generated maximum flood extent maps are recorded as the design flood inundation maps for the considered return periods. This enables us to create an archive of design flood inundation maps (DFMs) associated to the preselected return periods.

2.3. Assessing Flood Return Period Using a Localized Particle Filter

The particle filtering approach adopted to optimally combine the satellite PFM and the archive of DFMs aims to estimate the return period of the imaged event. The PF approach used here is based on that of [48]. The two main differences here are (1) PF is exploited to estimate a flood return period rather than to update a flood forecasting system and (2) PF is here localized to enable the estimation of flood return periods locally. Our implementation of PF is first presented on a general basis with localization carried out over subareas. Next, we propose various localization strategies based on hydrological decompositions of the model domain (sub-catchments). The latter represents the main contribution of the proposed approach, which enables the assimilation of flood observations in each hydrographic catchment and provides a range of return periods for a flood event.

Particle filters represent prior (before the assimilation) and posterior (after the assimilation) properties via a set of particles. They enable us to estimate the posterior probability of each model based on Bayes theorem:

$$p\left(x|y\right)=\frac{p\left(\mathrm{y}|x\right)p\left(\mathrm{x}\right)}{p\left(y\right)}\approx {\sum }_{i=1}^N{W}^i\delta (x-x_i)$$

(4)

where $p\left(x|y\right)$ is the probability of the model knowing the observation (i.e., the posterior), $p\left(y|x\right)$ is the probability of the observation knowing the model realization (i.e., the likelihood), p(x) is the probability of the model realization (i.e., the prior) and p(y) becomes a scaling factor so that the integral of the posterior is equal to 1, $W^i$ is the weight of the i-th particle, and $x_i$ is the i-th particle value. Further information and details on the theoretical basis of particle filters can be found in [29,37].

In a first step, as our observation is a probability to be water knowing the SAR backscatter, a pixel-based DFM likelihood $w_k^i$ is first assigned to each individual k-th pixel of the i-th DFM using the following formula (see [48] for further details):

$$w_k^i=p_k^i\left(w|{\sigma }^0\right)\times {\theta }_k^i+(1-p_k^i\left(w|{\sigma }^0\right))\times (1-{\theta }_k^i)$$

(5)

where $p_k^i\left(w|{\sigma }^0\right)$ is the probability of the k-th pixel being water in the observation; ${\theta }_k^i$ is the value of the k-th pixel in the i-th DFM being flooded ${\theta }_k^i=1$ or not flooded ${\theta }_k^i=0$. This likelihood represents the local probability of the observation given the considered simulated DFM.

In a second step, assuming that there is no spatial correlation in observational errors, the weights are multiplied within each predefined subarea, enabling us to estimate the likelihood ${\omega }_j^i$ over the j-th subarea:

$${\omega }_j^i={\prod }_{k=1}^{N_j}{w}_k^i$$

(6)

where $N_j$ indicates the number of pixels within the j-th subarea. Finally, the likelihoods are normalized in each subarea so that they sum up to 1, thus providing for each j-th subarea the weight $W_j^i$ of each i-th DFM, as follows:

$$W_j^i=\frac{{\omega }_j^i}{{\sum }_{i=1}^{Np}{\omega }_j^i}$$

(7)

The set of the weights $W_j^i$ evaluates the posterior probability of the DFMs given the observation. The estimated return period ${\widehat{T}}_j$ of the observed flood event in each subarea is finally computed as the expectation (weighted mean) over the j-th subarea:

$${\widehat{T}}_j={\sum }_{j=1}^{Np}{W}_j^i\times T^i$$

(8)

where Np is the number of DFM in the ensemble, and $T^i$ is the return period associated with the i-th DFM. An estimated flood inundation map is also derived from the assimilation of the satellite observation into the ensemble of simulated DFMs. To do so, the expectation of the simulated water depth map is first computed as the weighted mean (using $W_j^i$) of the ensemble of maximum water depth maps generated in Section 2.2. Next, it is converted into an estimated flood inundation map assuming a grid cell to be flooded if the water depth is larger than 10 cm and dry otherwise. As a consequence, the output of the whole procedure for any satellite observation is an estimated local return period ${\widehat{T}}_j$ and the associated estimated flood inundation map.

We define subareas as sub-catchments (i.e., drainage areas) as these enable us to distinguish various fluvial segments potentially affected by floods of various return periods within the catchment (Figure 3).

Topographic and morphological characteristics determine the limits of a catchment relative to a location along a river channel [58,59]. Here, the sub-catchment drainage areas are mapped through a standard hydrological analysis via a GIS platform [60,61,62]. Initially, a flow accumulation raster is generated from a digital elevation model [63,64] using the conventional method (including sink filling, flow direction, and accumulation map computation).

Subsequently, the hydrographic network is delineated by applying a threshold to the flow accumulation map, whereby drains and river channels are identified as commencing when a sufficiently extensive area is drained (Figure 3a). Following this, the channels are automatically segmented into sections through a river segmentation process that depends on detecting intersections of streams (Figure 3b). Lastly, the upstream regions that contribute flow to the sections determine the sub-catchments (Figure 3c). Depending on the threshold set on the flow accumulation map, the resulting sub-catchments vary in size, thereby altering the number of spatial subdivisions within the catchment.

3. Study Area and Experimental Set-Up

3.1. Study Area

The proposed methodology is applied to a downstream part of the River Severn basin around the city of Tewkesbury (UK). The River Severn has its source on the Welsh northeastern mountains and flows southeastward to the Vale of Gloucester and the Bristol Channel. Figure 4a shows the Severn catchment and river network, and Figure 4b shows the study area (hydraulic model domain) and the hydrometric gauging stations used as upstream boundary conditions of the shallow water model (black squares). The domain area is 30.5 km by 52.4 km and is located on the lower Severn. The River Teme and River Avon are tributaries of the River Severn within the model domain. The study area was rather frequently inundated over the past decades [32].

3.2. SAR Data Set and Processing

A satellite dataset composed of eleven SAR observations is available from ENVISAT ASAR imagery acquired in wide swath mode (WSM; 75 m pixel spacing). To obtain the σ⁰, the SAR images were preprocessed applying the following steps: thermal noise removal, radiometric calibration, and terrain correction. The images were acquired during flood events that occurred in March and July 2007, January 2008, and January 2010. The dates and times of the satellite acquisitions are reported in

Table 1. Satellite data acquisition date and time and relative ID. All SAR images are from the Envisat mission, acquired in wide swath mode (400 km swath width) with an incidence angle ranging from 17° to 43°.

SAR Acquisition	Date	Time	SAR ID
March 2007	5 March	10:27	1
March 2007	5 March	21:53	2
March 2007	8 March	10:34	3
March 2007	8 March	21:58	4
July 2007	23 July	10:27	5
July 2007	23 July	21:53	6
January 2008	17 January	21:55	7
January 2008	24 January	10:12	8
January 2008	24 January	21:38	9
January 2010	18 January	10:30	10
January 2010	18 January	21:53	11

. A PFM was derived from each SAR image using the method described in Section 2.1. The resulting 11 probabilistic inundation maps (Figure 5) exhibit flood extents of markedly variable magnitude, consequently corresponding to distinct flood scenarios, making these PFMs well-suited for testing the proposed methodology. In this study, we derived the forest and urban areas from the Corine Land Cover CLC 2018, with 100 m resolution for the exclusion layer; the corresponding pixels are masked out during the assimilation. The urban and forest classes are masked out as the flood map** algorithm used here is not able to detect water in these two landcover classes. Moreover, since the dates of acquisition of the SAR images and the dates of generation of CLC are not very close in time, we verified that land cover in these areas remained unchanged during the time interval between the year the CLC dataset was compiled (2018) and the year of the previous CLC version (2012). As a consequence, using either CLC 2012 or 2018 would not substantially change our results.

In addition to the satellite imagery, aerial photography was used to derive a ground truth flood extent map of the event that occurred on 24 July 2007 via photo-interpretation [51]. This aerial photograph campaign was organized for the day following two satellite acquisitions (SAR ID 5 and 6). This flood map was used to evaluate the flood extent (see Section 4) obtained from the assimilation procedure.

3.3. Generating an Archive of Design Flood Extent Maps Using a Shallow Water Model

Daily flow time series observed at hydrometric gauging stations were used to draw design flow hydrographs following the approach proposed in Section 2.2. Long streamflow time series, spanning more than 50 to 60 years on average, at Bewdley (River Severn), Knighsford Bridge (River Teme), Evesham (River Avon), Besford Bridge (Brook Bow), Kidderminster (River Stour), Harford Hill (River Salwarpe), and Hinton (River Isbourne) were considered. In this study, 15 scenarios were defined (see Section 2.2) to draw a meaningful number of scenarios and generate a representative ensemble of design flood extents.

The derived design hydrographs were used as input inflow for all contributing rivers (Severn, Avon, and Teme) and tributaries (Brook Bow, River Stour, River Salwarpe, and River Isbourne) of the hydraulic model at daily time steps. Here, the setup of Lisflood-FP is identical to the one of [32,48]. The model grid cell size is 75 × 75 m and the river channel is represented via the model sub-grid approach [55]. Surveyed cross-sections of the main rivers were used to define channel width and depth and the Manning’s roughness parameter was fixed to 0.0426 in the channel and to 0.06 in the floodplain for the entire domain. A “free” downstream boundary condition to fix the slope of the main channel was set to 0.00007, determined as the overall valley slope.

As proposed in Section 2.3, we derived maximum flood extent binary maps for every scenario from simulated water depth maps (Figure 6). These binary maps comprise the ensemble of design flood extent binary maps that are used as particles in the assimilation procedure.

3.4. Subcatchments as Subareas in the Data Assimilation

As introduced in Section 2.3, various subdivisions of the model domain are derived through hydrological analysis, which is applied in GIS by processing Shuttle Radar Topography Mission 1 DEM (SRTM1) data with a spatial resolution of 30 m and a height accuracy exceeding 9 m. The model domain was tested and evaluated (Figure 7): (a) no subdivision, (b) 3 subareas drawn from a thresholding of the flow accumulation at 2 × 10⁵ draining pixels, (c) 9 subareas drawn from a thresholding of the flow accumulation at 5 × 10⁴ draining pixels, and (d) 26 subareas drawn from a thresholding of the flow accumulation at 1 × 10⁴ draining pixels. The definition of these model domains is associated with the necessity of differentiating the various contributing areas of the respective river channels to take into account the dynamism of hydrological and geomorphological processes. Here, (a) considers the entire test area as part of a single basin (Severn catchment), (b) separates the Severn and Avon River catchments, (c) also separates the River Teme catchment, and (d) distinguishes other sub-catchments along minor segments of the three main rivers. Figure 7 shows the generated subarea maps with 1 spatial unit, 3 spatial units, 9 spatial units, and 26 spatial units.

4. Validation Approach

Based on the available data over the study area, the proposed data assimilation method is evaluated both in terms of the estimated flood inundation map and in terms of the estimated flood return period.

For flood extent evaluation, the estimated flood inundation map (one for each localization scenario) and flood inundation map from the SAR image (binary map from PFM; see Section 3.2) corresponding to the fifth SAR acquisition on 23 July 2007 (ID 5, see Table 1) are compared with the inundation map derived from the aerial photographs acquired on 24 July 2007 (see Section 3.2). This comparison is carried out through a visual interpretation and the computation of performance metrics, namely the critical success (CSI), the Cohen’s kappa [65], and the overall accuracy indexes. These are all based on the confusion matrix that is composed by the number of correctly estimated pixels (flooded/not flooded) and the number of incorrectly estimated pixels (associated with overestimation and underestimation). This comparative approach helps to evaluate the accuracy of an estimated flood extent map with respect to an assumed ground truth (here the flood map from aerial photographs).

For the return period evaluation, two comparisons are carried out. First, we compare the ranges of estimated return period ($\widehat{T}$) and the ranges of observed return period (T_g). Daily streamflow long-term time series (about 90–70 years, from 1920/40 to 2012) at the gauging stations of Bewdley (River Severn) and Evesham (River Avon) are considered for the estimation of observed return periods. Gauged discharge on satellite acquisition day and gauged discharge on peak day of flood event were used to extract the P probability from the GEV CDF curve (see Section 2.2) and determine the range of observed T_g at the gauging station (i.e., T_g in acquisition day—maximum past T_g during the event). This approach is based on the assumption that estimated $\widehat{T}$ relates to the flood extent observed by the satellite that does not necessarily correspond to synchronously gauged T_g at the upstream boundary condition [37]. Therefore, the range of observed T_g is compared to the range of estimated $\widehat{T}$ at the river scale (i.e., a single value of $\widehat{T}$ if the river is included in a unique sub-catchment or a range of $\widehat{T}$ if the river is running over several sub-catchments depending on the segmentation used). Second, we compare the ranges of estimated $\widehat{T}$ and the ranges of predicted return period (T_m) derived from a long-term simulation of the Lisflood-FP model (about 40 years) used to predict the water volume time series in sub-catchments. This simulation uses gauged daily streamflow time series as boundary conditions and provides a 40-year daily time series of water depth maps. These maps were next converted into the volume of water at the river scale. Based on this volume time series, GEV distributions were fitted (see Section 2.2) and return periods of this volume T_m were estimated for each SAR image acquisition time. Then, a comparison at river scale in the Severn and Avon Rivers was applied between predicted T_m and estimated $\widehat{T}$ (a single value of $\widehat{T}$/T_m if the river is included in a unique sub-catchment or a range of $\widehat{T}$/T_m if the river is running over several sub-catchments, depending on the segmentation used). The Severn and Avon Rivers are selected because their floodplains are large enough to enable water detection on SAR images like ENVISAT, while other rivers with very small floodplains are not.

5. Results and Discussion

5.1. Localization

The test case over image ID 5 acquired on 23 July 2007 was first selected to analyze the results of the proposed data assimilation method using various localization sub-catchments (Figure 8), as this image exhibits the largest observed flood inundation extent.

As shown in Figure 8a, the estimated return period is 5 years all over the domain when no localization is carried out. Using a localization based on three sub-catchments enables us to distinguish the Severn and Avon Rivers and to estimate different return periods for these two rivers, namely 5 years for the Severn and 20 years for the Avon (Figure 8b). Going further in the subdivision of study area (i.e., 9 and 26 sub-catchments) enables us to estimate the return period at a smaller scale. Figure 8c and d show a continuous decrease of the estimated return period along the River Severn. This is consistent with the fact that the image was acquired after the flood peak, when the flow recession was starting. However, in Figure 8d, the estimated return periods over the River Avon are decreasing toward the downstream. This variation in the estimated return periods demonstrates the relative sensitivity of the proposed method to the domain segmentation: fairly small reaches considered along the River Avon and fairly narrow river valley in corresponding sub-catchments make these areas more challenging for the method as the sensitivity of the flood extent to the return period is limited.

The results of the comparison of the estimated flood maps and SAR flood map corresponding to the fifth SAR acquisition on 23 July 2007, along with the observed flood map (derived from aerial photography), are shown in Figure 9 and

Table 2. Confusion matrix and performance coefficients for flood extent validation. OF = observed flood; ON = observed not flooded; AF = assimilated flood; AN = assimilated not flood; SOF = SAR observed flood; SON = SAR observed not flood; CSI = critical success index.

		Confusion	Matrix		Indices
Experiment		OF	ON	Accuracy	CSI	Kappa
Assimilation 1 spatial unit	AF	7153	1796	0.912	0.681	0.752
	AN	1562	27,402
Assimilation 3 spatial units	AF	7572	1863	0.921	0.716	0.782
	AN	1143	27,335
Assimilation 9 spatial units	AF	7408	1373	0.929	0.734	0.801
	AN	1307	27,825
Assimilation 26 spatial units	AF	6983	1416	0.917	0.689	0.763
	AN	1732	27,782
SAR image	SOF	5845	804	0.903	0.614	0.702
	SON	2870	28,394

. In Figure 9, correctly estimated flood/not flood pixels are shown in dark and light blue, while overestimation and underestimation are respectively represented in orange and yellow. Table 2 reports the metrics derived from this comparison for all types of experiment, i.e., the confusion matrix and the accuracy, CSI, and kappa indices.

Figure 9 shows that flood maps derived from the proposed data assimilation method are generally in good agreement with the ground truth flood extent. Moreover, the flood extent is better estimated as a result of the data assimilation compared to the one derived from the SAR image. Green and red boxes in Figure 9e highlight underestimation in the SAR-derived flood map along the River Avon and in the urban area, respectively. Furthermore, the estimated flood extent obtained using localization is more accurate than that obtained without localization (Figure 9a–d). The quantitative analysis confirms this evaluation as all metrics are higher with localization.

The indices in Table 2 indicate that the best results in the comparison between estimated flood maps and validation flood maps belong to the experiments with a domain segmentation of three and nine spatial units. Considerations from this evaluation approach are (1) the proposed approach enables a better estimation of the flood extent, with an improved accuracy compared to the flood maps derived from the SAR image; (2) in the extreme event analysis and in the estimation of inundation scenario, flood phenomena should be observed distinctly on the River Severn, River Avon, and the River Teme. The segmentation should thus separate the three corresponding sub-catchments, without the need for further splitting the three main channels.

5.2. Assessment of Estimated Return Periods

The estimated flood return periods were assessed by applying a comparative approach at river scale. The ranges of $\widehat{T}$ were compared to T_g and T_m (see Section 4) for six SAR acquisitions (ID 1, 3, 5, 7, 8, 10; see Table 1) by selecting only one of those images acquired the same day. The results are shown in

Table 3. Ranges of observed return period and ranges of estimated return period along the Rivers Avon (grey cell) and Severn (green cell). T_g = return period derived from long observation time series at Evesham and Bewdley gauge stations; $\widehat{\mathrm{T}}$ = estimated return period using the proposed data assimilation method.

SAR ID	SAR Observation (Date, Time)	Tg	$\widehat{\mathit{T}}$
			No Loc	Loc on 3	Loc on 9	Loc on 26
1	5 March 2007 10:27	6.35	1.01	2.00	2.00	1.01–20.00
3	8 March 2007 10:34	1.00–4.35	1.01	2.00	2.00	1.01–10.00
5	23 March 2007 10:27	4.02–308.55	5.00	20.00	20.00	1.01–500.00
7	17 January 2008 21:58	2.54–5.10	1.01	2.00	2.00	1.01–5.00
8	24 January 2008 10:12	1.00–1.13	1.01	1.01	1.01	1.01–2.00
10	18 January 2010 10:10	1.10–3.19	1.01	2.00	2.00	1.01–30.00
1	5 March 2007 10:27	1.00–1.02	1.01	1.01	1.01	1.01
3	8 March 2007 10:34	1.01–1.02	1.01	1.01–5.00	1.01–5.00	1.01–5.00
5	23 July 2007 10:27	1.52	5.00	5.00	1.01–10.00	1.01–10.00
7	17 January 2008 21:58	1.10	1.01	1.01–2.00	1.01–2.00	1.01–2.00
8	24 January 2008 10:12	1.14–2.42	1.01	1.01	1.01	1.01
10	18 January 2010 10:10	1.06	1.01	1.01	1.01	1.01–5.00

and Figure 10.

In Table 3, the first columns contain the SAR image ID and the satellite image acquisition date. The following columns contain ranges of gauged return periods and estimated return periods. The T_g are relative to intervals of return periods observed at Evesham and Bewdley gauging stations (grey and green cells respectively); $\widehat{T}$ are relative to intervals of return periods estimated for the Rivers Avon and Severn. Then, a single value of T_g means that the observed return period at the time of SAR acquisition is equal to maximum return period during the flood event and a single value of $\widehat{T}$ means that the river is included in a single area or that estimated return periods are equal in the different sub-catchments of the river. Sub-columns in $\widehat{T}$ distinguish values of estimated return periods for different localizations of the assimilation filter.

Figure 10 graphically shows the comparison between ranges of estimated return period and ranges of predicted return period (described by error bars) for each subdivision of the model domain (different color in the figure) along the Rivers Avon and Severn (Figure 10a and b respectively). The results show that the estimated return periods T_m/$\widehat{T}$ vary depending on the analysis area. In some cases, the return period is specific to a single sub-catchment (like the Avon River in “loc 3” and “loc 9” of Figure 10a).

The results confirm two key points. First, there is a good agreement between the estimated flood scenarios, the observed floods, and the predicted floods. This demonstrates the effectiveness of the data assimilation method for assessing flood frequency. Second, segmenting the river by location (see “loc 3” and “loc 9” in Table 3 and Figure 10) improves the estimation of flood scenarios. This highlights the importance of considering the river scale as a spatial unit in flood hazard assessment. A key advantage of our method is that it can provide average or distributed return periods within sub-catchments, accounting for the variations in river and floodplain processes across the area. However, the results also indicate a limitation (see the range of values for “loc 26” in Table 3 and Figure 10). While location-based filtering is beneficial, it is important to maintain a distinction between the main river courses and avoid excessive segmentation within a sub-catchment, as this can reduce the method’s effectiveness.

When evaluating the method’s ability to estimate flood scenarios, it is important to consider how return periods are typically determined. Traditionally, return period refers to the peak flood discharge measured at a specific gauging station. However, our method estimates the return period based on the flood extent captured in a SAR image, which covers a larger area and reflects upstream conditions that may differ from those at the gauging station. This approach provides a range of possible return periods instead of a single value, leading to a more realistic but less precise hazard assessment.

In simpler terms, while NRT data assimilation from SAR images allows us to estimate return periods at the river scale, these estimates are expressed as ranges because flooding is dynamic. A single SAR image cannot capture the peak flood extent across the entire area of interest at a specific moment. Additionally, by the time floodwaters reach the floodplain, the upstream flow might already be receding. This is why our validation process compares the estimated return periods to two benchmarks: the observed range of return periods at gauging stations and the return periods predicted through long-term simulations (explained in Section 4). The positive outcomes of this comparison demonstrate the method’s effectiveness in integrating data from both SAR and models, while also considering variations in river hydrology and floodplain characteristics throughout the catchment.

Our study demonstrates the value of our NRT flood frequency assessment for monitoring flood events. Indeed, the enhanced accuracy in the estimation of flood extent and flood scenarios and the availability of a spatially distributed return period represent the main advantages of our data assimilation method. The estimated flood maps help visualize how floodwaters behave in different parts of the river, highlighting areas potentially affected by the flood’s impact. In the perspective of flood monitoring, a map of the return period ranges for different sub-catchments directly addresses the needs of flood management across various zones within the river catchment.

Our method does have some limitations to consider. First, it relies on the spatial and temporal resolution of satellite data. The flood event needs to be visible in SAR images for the data assimilation process to work effectively. Ideally, the method would be applied to a large collection of SAR observations to capture the full extent of potential floods in the study area. Additionally, the accuracy of flood probability predictions is limited by the number of DFMs used. A larger archive of simulated scenarios would increase the likelihood of matching the observed flood extent to the most accurate DFM.

These limitations can potentially be addressed in future advancements. As the availability of satellite data and faster hydrodynamic modelling techniques continues to improve, the method’s effectiveness should increase. Furthermore, applying the method with global models could enable large-scale flood analysis around the world, overcoming restrictions related to data and tool availability

6. Conclusions

This study introduces a novel method for near real-time assessment of a flood event’s return period. The method integrates flood extent maps derived from SAR data with pre-computed flood maps using a localized particle filter. We evaluated the approach in the Severn River basin and found that it achieves the following:

• Confirms the significant value of SAR data for flood monitoring.
• Improves the accuracy of flood extent maps and enables real-time return period estimation.
• Provides spatially distributed return periods, accounting for the variations in flood processes along the river.

This method has the potential to be a valuable tool for flood monitoring, especially considering the growing need from the insurance industry to estimate flood severity based on spatially distributed return periods. For example, this information is crucial for improving the sector’s parametric insurance models, which rely on objective data to determine payouts.

7. Patents

R. Hostache, M. Chini, P. Matgen, M. Zingaro. Method for quantifying the extent and recurrence of natural events. Filed: 18/12/2020. EP20210840907, CN116670728(A), LU102319(B1), US2024046639(A1), WO2022129611(A1).