Determination of Current and Future Extreme Sea Levels at the Local Scale in Port-Bouët Bay (Côte d’Ivoire)

Abstract

The Port-Bouët Bay shoreline is threatened by extreme sea level (ESL) events, which result from the combination of storm tide, wave run-up, and sea level rise (SLR). This study provides comprehensive scenarios of current and future ESLs at the local scale along the bay to understand the evolution of the phenomenon and promote local adaptation. The methodological steps involve first reconstructing historical storm tide and wave run-up data using a hydrodynamic model (D-flow FM) and the empirical model of Stockdon et al. Second, the Generalized Pareto Distribution (GPD) model fitted to the Peaks-Over-Thresholds (POT) method is applied to the data to calculate extreme return levels. Third, we combine the extreme storm tide and wave run-up using the joint probability method to obtain the current ESLs. Finally, the current ESLs are integrated with recent SLR projections to provide future ESL estimates. The results show that the current ESLs are relatively high, with 100-year return levels of 4.37 m ± 0.51, 4.97 m ± 0.57, and 4.48 m ± 0.5 at Vridi, Petit-Bassam, and Sogefiha respectively. By end-century, under the SSP1-2.6, SSP2-4.5, and SSP5-8.5 scenarios, the future SLR is expected to increase the current ESLs by 0.49 m, 0.62 m, and 0.84 m, respectively. This could lead to a more frequent occurrence of the current 100-year return period, happening once every 2 years by 2100, especially under SSP5-8.5. The developed SLR scenarios can be used to assess the potential coastal flood risk in the study area for sustainable and effective coastal management and planning.

1. Introduction

Extreme sea levels (ESLs) pose a serious threat to many coastal areas around the world [1,2,3,4,5,6,7,8,9,10,11]. They result from the combination of different ocean dynamic factors at the coast [3,9,10,11,12,13,14]. These factors or ESL components mainly include astronomical tides, storm surges, ocean waves and sea level rise. Each of these factors constitutes a separate process but together can result in extreme sea levels. In general, when at least one among these factors is at an uncommonly high level, the resulting total coastal water level can exceed a certain threshold [15], leading to devastating consequences, especially along densely urbanized coastal areas. The potential impacts of ESLs on the open coast mainly include coastal erosion and flooding [16].

Over the past few decades, ESL events have become more frequent and intense on a global scale, as underscored by recent studies [3,12,17,18,19]. The Ivoirian coastal areas are not spared from this situation. The city of Port-Bouët is one of the sites most affected by this phenomenon, where ESLs causing significant damage have been reported since the 1940s [20]. More recently, on August 2007, dozens of houses were destroyed by an ESL event in which the offshore wave contribution was estimated to be 2.87 m high [21]. In August 2011, three consecutive days of ESLs were observed, which caused the destruction of hundreds of homes and left more than a thousand families without shelters [22]. Today, coastal populations are experiencing ESL events almost every year. Global projections also showed future exacerbation of the phenomenon and its impacts due to rising sea levels induced by climate change [16,18,23]. According to [24,25], the once-a-century ESL event that is expected to occur in the present day will probably become an annual event by 2100 in many locations on the planet, even with the lowest sea level rise scenario. Additionally, the tropical region, to which Côte d’Ivoire belongs, is reported to experience the most dramatic future ESL increases.

These observations and research outcomes have led to extensive mitigation and adaptation planning efforts in order to reduce the potential physical and societal impacts associated with ESL events [26,27,28]. However, effective decision-making and planning require a good understanding of present-day ESLs and accurate projections at appropriate spatial and temporal scales. ESL determination can be approached at different scales, including local, regional, and global scales, but the local-scale approach appears more suitable for identifying extremes for planning and adaptation purposes [29,30]. This is because local climate, geological framework, and coastal morphological features strongly influence factors contributing to ESLs [29], causing ESLs and their consequences to vary by location and become more site-specific. Therefore, the local-scale approach would allow for better quantification of spatial changes in ESL frequency and intensity, as well as their impacts within a given coastal area, compared to larger-scale approaches [30].

There is a limited number of studies that have considered the aforementioned factors in determining ESLs along the Ivorian coast. This shows that questions regarding the variations in ESL elevations along the coastal areas of Côte d’Ivoire and, more particularly, along the Port-Bouët Bay are yet to be answered. The present study aims to fill these gaps by providing a better understanding of present-day ESLs and their potential evolution under climate change by the end of the century at a local scale. The objective is to develop robust ESL scenarios that can help coastal flood management and development planning for the different sectors of the study.

This paper consists of five sections, including this first introductory section. Section 2 provides a description of the study area, the collected data, and the methodology used. The presentation of study results can be found in Section 3, while Section 4 contains the discussion. Finally, the study’s conclusions are presented in Section 5.

2. Materials and Methods

2.1. Study Area

Port-Bouët Bay is a small extension of the Atlantic Ocean in the coastal area in the central part of Port-Bouët, one of the 10 towns, which constitute the Ivorian economic capital Abidjan. Its presence is related to higher rates of shoreline recession in this location [31]. This strong erosive tendency is mainly due to the narrowness of the continental shelf, positioned on the extension of a submarine canyon, and the total cut-off of the littoral drift (from West to East) generated by the opening of the Vridi canal [31,32]. The bay is located between the latitudes 5.247° N–5.255° N and the longitudes 3.951° W–4.002° W. It stretches for 5.5 km from the Vridi canal to the city’s lighthouse. It is connected to the waters of the Ebrié lagoon by the Vridi canal. The coastline of this bay is the only one in Port-Bouët that, being exposed to ESL events, has not been subjected to relocation yet; i.e., key infrastructures and an important number of people still live very close to the shoreline. This explains the interest in this site within the framework of the present study. The coastline of Port-Bouët Bay has a W–WSW and E–ESE oriented concavity and can be divided into three sectors based on the bordering neighbourhoods: from west to east, Vridi, Petit-Bassam and Sogefiha (see Figure 1). To be in line with the local-scale approach advocated in this study, these different neighbourhoods will be used as the local unit for the analysis.

2.2. Approach Overview

An extreme sea level event is caused by a combination of factors, including tide, storm surge, wave and sea level rise [3,9,10,11,12,13,14]. Tide refers to the daily rise and fall observed in coastal water levels every day due to the gravitational pull of the moon and sun. Storm surge is an abnormal rise of water generated by a storm on top of the expected tide. The combination of storm surge and the tide is known as storm tide. Waves contribute to extreme sea levels through the wave run-up process, which is the increase in the mean still-water sea level as ocean waves interact with the coastline [33]. Sea level rise refers to the long-term increase in mean sea level (MSL) due to global warming and land subsidence. In this study, we define extreme sea level (ESL) as the sum of storm tide (ST), wave run-up (R_2%) and sea-level rise (SLR) as shown in Equation (1) with all estimates based on the local MSL datum.

$$\mathrm{ESL}=\mathrm{ST}+{\mathrm{R}}_{2\%}+\mathrm{SLR},$$

(1)

At any given time, these three components can be combined to calculate the extreme elevation(s) that the sea may reach. However, due to their different time-scale variations, accurate evaluations of contributions from parameters such as storm tide and wave run-up can be made better for the present time. Therefore, assuming that the present SLR is negligible, we define the present-day and the future ESLs as expressed in (2) and (3), respectively.

$${\mathrm{ESL}}_{\mathrm{present}-\mathrm{day}}=\mathrm{ST}+{\mathrm{R}}_{2\%},$$

(2)

$${\mathrm{ESL}}_{\mathrm{future}}={\mathrm{ESL}}_{\mathrm{present}-\mathrm{day}}+\mathrm{SLR},$$

(3)

Solving these equations is equivalent to determining each element that composes them. However, the process is not straightforward because extreme sea level (ESL) events are rare and can only be accurately understood by observing sea level over a long period of time [33]. In other words, the calculation of ESL requires historical sea level data that is both long and of high quality. For accurate ESL calculations, historical sea level data should ideally be measured at a reference point (datum), sampled at hourly intervals, and free of any data gaps. Meeting these criteria is generally challenging, particularly for developing coastal countries that may lack the necessary resources and equipment. To overcome the lack of quality data, techniques have been developed to calculate extremes from shorter records [34]. Other techniques involve reconstructing historical sea level time series to increase the necessary data length and/or to fill any gaps. The latter techniques have the advantage of enabling efficient extreme analysis and thus provide accurate estimates. In this study, we opted to reconstruct sea level records over the past 40 years for both storm tide and wave run-up using various modeling methods. Figure 2 provides an overview of the methodological approach. Below, we provide detailed explanations of the modeling process, statistical methods, and data used.

2.3. Data

To perform this study, some data were necessary. They fall into two categories. On one side, the metocean data, which include offshore wave, wind, atmospheric pressure, tide gauge data and sea level rise projections. On the other side, the morphological data include beach profiles and bathymetry data.

For the study area, water level data are available through the tide gauge stations of the Port of Abidjan. We collected the 33 years (from 1983 to 2016) of available measurements at the former tidal station “Quai-Nord” located on latitude 5°18.230′ N and longitude 4°1.570′ W. This hourly dataset was mostly digitized from archived records and has undergone extensive quality checks as part of the study of [35], and has been made available for use in this study. Even though the collected dataset has an acceptable overall length, it is limited to the port domain and contains several gaps ranging from a minimum of a single hour to a maximum of 2.6 years of consecutive missing values. The total gaps represent about 46% of the overall length, making the dataset fail to meet the quality requirements for direct use in extreme value analysis. Alternatively, it was used as a base dataset for storm-tide modeling.

Offshore wave, wind and atmospheric pressure data were all obtained from the ERA5 global reanalysis database [36]. In Côte d’Ivoire, there is a huge lack of in situ ocean climate data. In such a situation, the ERA5 reanalysis data, with their acceptable performance in several regions of the globe [37,38], constitute a good alternative for providing an overview of the offshore condition and overcoming the lack of in situ measurement. For the analysis, we extracted 40 years (1981–2020) of historical hourly data for each pre-cited process. The wave data include the significant height, period and direction of the wave. The collected wind data are the zonal and meridional components of the surface wind field.

The information on the sea-level rise for Port-Bouët was obtained from the 6th IPCC AR6 Sea-Level Projection Tool [39,40,41], which provides estimates of future sea-level changes relative to a baseline of 1995–2014 at the regional level. In their latest sea level projections, the IPCC recommends the use of the Shared Socioeconomic Pathways (SSPs) scenarios, which encompass a wider range of greenhouse gas and CO₂ emissions than previous scenarios [42]. The core set of five SSP scenarios is designed from combinations of processes for which projections can be made with medium confidence, while two additional low-confidence scenarios account for deeply uncertain processes, resulting in a total of seven SSP scenarios. For this study, we focus on three medium-confidence scenarios, namely SSP1-2.6, SSP2-4.5, and SSP5-8.5, which correspond to low, intermediate, and very high greenhouse gas and CO₂ emissions, respectively. Figure 3 illustrates the SLR projections for these three SSP scenarios at the coordinate point (4°0.000′ N; 4°0.000′ W) off the coast of Port-Bouët.

The morphological data used in this study included bathymetry and beach profiles, which were collected from a variety of secondary sources. Bathymetry data for the areas of interest and surrounding regions were obtained from multiple sources, including the General Bathymetric Chart of the Oceans (GEBCO) at 15 arc-second intervals [43], a lagoon bathymetry survey by Shore Monitoring and Research [44], and a 2019 bathymetric survey conducted by the port of Abidjan. Historical beach profiles from previous studies [31,45,46,47], covering the period from 2007 to 2016, were also used to provide information on beach slope variations across the study area. It should be noted that although other morphological factors, such as vertical land movements, could significantly influence ESLs in certain regions [48,49], this study did not take these factors into account directly. This was mainly because of the limited information available on past uplift or subsidence in the study area. However, projections for future land motion have been accounted for in sea-level rise projections, with an estimated subsidence of about 5 mm by 2100 and were considered as such in this assessment.

Figure 3. Sea level rise projections for the SSP1-2.6 (in blue), SSP2-4.5 (in yellow), and SSP5-8.5 (in red) scenarios off the coast of Port-Bouët, based on data from the NASA Sea Level Projection Tool. The solid lines indicate the median range, while the shaded areas represent the 17th to 83rd percentile ranges (the likely ranges).

Figure 3. Sea level rise projections for the SSP1-2.6 (in blue), SSP2-4.5 (in yellow), and SSP5-8.5 (in red) scenarios off the coast of Port-Bouët, based on data from the NASA Sea Level Projection Tool. The solid lines indicate the median range, while the shaded areas represent the 17th to 83rd percentile ranges (the likely ranges).

2.4. Methods

2.4.1. Storm Tide Modelling

Hydrodynamic modelling of storm tide allows spatial extrapolation of water level from the tide-gauge locations (generally sheltered areas) to other surrounding areas, e.g., the open coasts. It also helps in filling small to large gaps in the water level records. In this work, storm tide simulations were undertaken using the D-Flow Flexible Mesh software developed by Deltares. This tool solves the depth-averaged shallow water equations. Detailed information on the model formulation can be found in [50]. The model domain includes the marine environment and the entire Ebrié lagoon system, both being connected through the Vridi canal (see Figure 4). This was necessary to account for the interactions between the sea and the lagoon system, which could somehow influence water levels in the area of interest. The computational grid was updated from the one of Brière et al. [51]. It is an unstructured grid made up of curvilinear and triangular meshes. The mesh resolution gradually increases as we move towards the study area and the port domain where the tidal station is located. Orthogonality and smoothness are two important grid properties for model accuracy. They respectively denote perpendicularity between flow links and net links and the ratio of neighboring grid cell dimensions. The grid is in good agreement with the orthogonality and smoothness thresholds defined in the modeling guidelines [52] and therefore satisfies the grid quality requirements. The depth schematization of the computational grids has been carried out by interpolating the different sets of bathymetric data described earlier.

The model was forced with tidal and meteorological forcing conditions. The tidal forcing at offshores boundaries as water levels were generated using ten astronomical constituents (M2, S2, N2, K2, K1, O1, P1, Q1 and M4) obtained from the DTU10 global ocean tide model [53]. The meteorological forcings, composed of hourly data of near-surface wind and atmospheric pressure, were uniformly inputted over the model domain.

The model calibration and validation were performed through comparisons against the water level measurements collected at the tidal station “Quai Nord” of the Port of Abidjan. The calibration process involved parameter changes to achieve the best possible match between the model and observation. For the present study, model calibration mainly focused on the uniform bottom roughness parameter represented by manning’s coefficient $\left(n\right)$. As suggested by Brière et al. [51], a range between $n=0.015 \mathrm{s}/{\mathrm{m}}^{1/3}$ and $n=0.030 \mathrm{s}/{\mathrm{m}}^{1/3}$ was first preselected. This range actually reflects the complexity of the area, characterized by a deep sea and shallow lagoon connected together by a relatively narrow channel (the Vridi canal). Then a series of sensitivity runs were carried out to find the best uniform roughness value. Calibration results showed that the best possible agreement between modelled and observed storm tides is attained when the Manning coefficient $\left(n\right)$ is equal to $0.029 \mathrm{s}/{\mathrm{m}}^{1/3}$.

The validation plot is presented in Figure 5, in which the modelled and observed storm tide levels at the tide-gauge station are compared for the period of 10 December 2014 to 2 January 2015. A critical look at Figure 5 shows that both the modelled and observed storm tides are in phase, but there is a minimal difference between the amplitudes. However, the calculated model performance through the Nash–Sutcliffe efficiency (NSE) coefficient [54] was estimated as 0.79, which is between 0.75 ≤ NSE ≤ 1.00. This indicates a “very good” model performance in simulating the storm-tide conditions for the study areas. Following this validation process, the storm-tide model was computed from 1981 to 2020 (plus an additional one-month “December 1980” for spin-up time) at the different locations within the research area. For each of the selected neighborhood, two to three observation points were set in the model, and their outputs were averaged to get the corresponding storm-tide time series that will be used in the extreme value analysis.

2.4.2. Wave Run-up Calculation

Several different models exist to determine the wave run-up, each with its derivation methodology and intended use [55,56,57]. For this calculation, the empirical model of Stockdon et al. [58] will be used. It was developed using field measurements collected from several sandy beach conditions [59,60]. We, therefore, expect this model to be appropriate for the sandy coast of Port-Bouët. The model’s formulation is presented as (4), in which ${\mathrm{R}}_{2\%}$ representing run-up 2% exceedance elevation is parameterized as a function of deep-water significant wave height (${\mathrm{H}}_0$), wavelength (${\mathrm{L}}_0$), and the foreshore beach slope (${\mathsf{\beta }}_{\mathrm{f}}$).

$${\mathrm{R}}_{2\%}=1.1{\left.\left(0.35\tan {\mathsf{\beta }}_{\mathrm{f}}{\left({\mathrm{H}}_0{\mathrm{L}}_0\right)}^{0.5}+0.5{\left[{\mathrm{H}}_0{\mathrm{L}}_0(0.563\tan {\mathsf{\beta }}_{\mathrm{f}}^2+0.0004)\right]}^{0.5}\right.\right)}^{},$$

(4)

The wavelength (${\mathrm{L}}_0$) can be calculated using wave period via the relationship in (5).

$${\mathrm{L}}_0=\frac{\mathrm{g}{\mathrm{T}}_0^2}{2\mathsf{\pi }},$$

(5)

The beach slope parameter ${\mathsf{\beta }}_{\mathrm{f}}$ in this study is defined, according to Ruggiero et al. [61] and Stockdon et al. [58], as the linear average slope over the beach face, which is the region of observed run-up $\pm 2\mathsf{\sigma }$, where $\mathsf{\sigma }$ is the standard deviation of the existing water-level record. ${\mathsf{\beta }}_{\mathrm{f}}$ was derived from historical beach profiles surveyed at several points across the study area. Slopes were first calculated for individual profiles, then averaged for each survey point. The representative beach slopes of each neighborhood (see

Table 1. Representative beach slopes (β_f) calculated from historical beach profiles.

Location	Vridi	Petit-Bassam	Sogefiha
Representative ${\mathsf{\beta }}_{\mathrm{f}}$	0.134	0.159	0.139

) were obtained by averaging the slopes of surveyed points along each neighborhood.

Equation (4) was used to reconstruct historical wave run-up time series of the same length as the extracted offshore wave climate, i.e., 40 years of hourly wave run-up data. However, the validation of these estimations is crucial before their use in the various analyses. Checking the validity of run-up estimates in a wave data scare location like in this study area is not easy work. Here, the calculated wave run-up was simply compared to the observed run-up elevation derived from satellite imagery [62]. The run-up boundaries can be seen on satellite images as the demarcation lines between the wet and dry beaches (Figure 6). Wet beaches generally appear darker in color than dry beaches on satellite images.

The demarcation lines have been digitized in order to extract their average altitude from the available topographical data. This obtained average altitude actually constitutes the observed total sea level. So, the prevailing storm tide value was removed from it to get the observed wave run-up that was then compared to the calculated wave run-up for the date of acquisition of the satellite image. Two statistical parameters, namely root mean square error (RMSE) and scatter index (SI), were used to assess the performance of the applied empirical model in each neighborhood (

Table 2. Statistical error parameters were obtained from a comparison between observed and calculated sea levels.

Location	Vridi	Petit-Bassam	Sogefiha
RMSE	0.376	0.479	0.279
Scatter Index (SI)	0.231	0.201	0.182

). They both evaluate the goodness of fit between observation and model. Low values (closer to zero (0)) of the RMSE and SI denote a better model performance. From Table 2, it can be seen that all the RMSE and SI values are all less than 0.5, indicating an acceptable performance and, therefore, the validity of the wave run-up estimates.

2.4.3. Extreme Value Analysis

Extreme Value Analysis (EVA) is generally used to determine the magnitude (height) and the return period of large and unusual events, e.g., extreme storm tides and waves. Here, the term ‘return period’ describes the likelihood of occurrence and, more particularly, the average recurrence interval of extreme events. It actually indicates the duration of time (in years) corresponding to the chance that an extreme event of a given magnitude will be reached or exceeded once in any one year.

The EVA theory is typically based on the extrapolation of past sea level records. Two main statistical methods emerge from the literature when it comes to carrying out an EVA. These include the annual block maxima (BM) and the peaks-over-threshold (POT) methods, which are, applied using, respectively, the generalized extreme value distribution (GEVD) and the Generalized Pareto Distribution (GPD) models [63]. However, each of these methods has a specific data requirement, i.e., the choice of which method to use depends on the characteristics of the available data (e.g., sampling frequency, erroneous/missing values and length). In general, the methods that make use of more of the available data are more accurate than those that use less data, which are, in contrast, easier to apply [64]. For this case, looking at the 40 years of reconstructed hourly sea level, we find it efficient to apply the GPD-POT method. Equation (6), where $\mathsf{\sigma }$ and $\mathsf{\xi }$ are, respectively, the scale and shape parameters of the distribution and $\mathrm{u}$ is the reference threshold, gives the GPD function.

$${\mathrm{F}}_{\mathrm{u}}\left(\mathrm{y}\right)=1-{(1+\mathsf{\xi }\frac{\mathrm{y}}{\mathsf{\sigma }})}^{-1/\mathsf{\xi }}$$

(6)

The GPD-POT method is designed to calculate the likelihood of exceeding a pre-determined threshold ($\mathrm{u}$). In this approach, it is assumed that the sea levels above the threshold in the time series have a Poisson distribution and are independent and identically distributed from the GPD [65,66,67]. The method involves the following steps:

• Defining a high threshold.
• Extracting the exceedances time-series, i.e., selection of all values above the defined threshold.
• Definition of clusters of exceedances based on the maximum storm duration (in hours). This separates the meteorological conditions that create extremes. Here, we took 72 h for the wave parameter and 25 h for storm tide as recommended by [14,68,69], respectively.
• Declustering: This step consists of selecting the peak value within each cluster in order to ensure the independency of excesses over the defined threshold and therefore justify the use of the GPD model.
• Finally, the calculation of return level: Values that are exceeded in any given year.

The determination of a return level ${\mathrm{Z}}_{\mathrm{m}}$ associated with a m-year return period using the GPD-POT method is given by (7) in which ${\mathsf{\lambda }}_{\mathrm{u}}$ is the rate of the Poisson process.

$${\mathrm{Z}}_{\mathrm{m}}=\left\{\begin{array}{cc}\mathrm{u}+\frac{\mathsf{\sigma }}{\mathsf{\xi }}\{{\left({\mathsf{\lambda }}_{\mathrm{u}}\mathrm{m}\right)}^{-\mathsf{\xi }}-1,\hfill & \mathrm{f}\mathrm{o}\mathrm{r} \mathsf{\xi }\ne 0\hfill \\ \mathrm{u}+\mathsf{\sigma }\ln \left({\mathsf{\lambda }}_{\mathrm{u}}\mathrm{m}\right),\hfill & \mathrm{f}\mathrm{o}\mathrm{r} \mathsf{\xi }=0\hfill \end{array}\right.$$

(7)

For this work, the whole process of EVA, including the threshold selection, was carried out using the software ORCA (metOcean data tRansformation, Classification and Analysis), developed by Deltares [70]. This MATLAB-based tool is composed of four main modules, with each of the modules covering a specific aspect of the metocean data analysis. To perform the analysis, two of these modules were particularly useful: the “data validation” module to inspect and structure the datasets and the “extreme conditions” module to apply the GPD-POT method to the data.

2.4.4. Joint Probability Analysis (JPM)

Joint probability refers to the chance of two or more separate processes occurring together that their combination may lead to a critical event [68,71]. This simply means that the probability of an ESL event is related to the combined probability of occurrence of all the components involved in that event, e.g., extreme storm tide and wave run-up in the case of present-day ESL. However, a m-year return period of storm tide will not always occur at the same time as a m-year return period of wave run-up, and if they do, so the resulting event would certainly be different from an ESL of the m-year return period. This is because, very often, ESL components depend on each other to a degree. The dependency between components indicates the way they should be combined to get the total ESL event. The JPM then utilizes the statistical dependence $\left(\mathsf{\chi }\right)$ between storm tide and wave conditions to calculate the likelihood of compound ESL through the relationship given in (8) [71], in which $\mathrm{T}$ stands for the return period and $\mathsf{\chi }$ values being between 0 and 1.

$${\mathrm{T}}_{{\mathrm{ESL}}_{\mathrm{Present}-\mathrm{day}}}=\sqrt{\raisebox{1ex}{${\mathrm{T}}_{\mathrm{ST}} \times { \mathrm{T}}_{{\mathrm{R}}_{2\%}}$}\!\left/ \!\raisebox{-1ex}{${\mathsf{\chi }}^2$}\right.},$$

(8)

However, the determination of the dependency factor $\left(\mathsf{\chi }\right)$ can be an exhausting process, which needs special statistical tools in some cases. For simplicity in this work, we assumed independence between the storm tide and wave run-up (i.e.,$\mathsf{\chi }=0$), then their joint return period would be the product of individual return period as given by (9).

$${\mathrm{T}}_{{\mathrm{ESL}}_{\mathrm{present}-\mathrm{day}}}={\mathrm{T}}_{\mathrm{ST}} \times {\mathrm{T}}_{{\mathrm{R}}_{2\%}}$$

(9)

Using (9), the joint return periods of various storm tide and wave run-up combinations were calculated. For the resultant joint return level, the individual return levels were summed $\left(\mathrm{S}\mathrm{T}+{\mathrm{R}}_{2\%}\right)$.

2.5. ESL Scenario Design

As the present work objective is to determine the present and future ESLs and knowing that ESL is governed by different dynamic factors, it is necessary to define scenarios that depict each case that can be faced. Scenarios are also very important for further exposure and risk assessment. ESL scenarios within the framework of this study are established as follows:

• For the present-day ESLs, six scenarios are defined based on the probability of occurrence of extreme events. These correspond to the best-fit return values for each of the 1, 5, 10, 25, 50- and 100-year return periods.
• For future ESLs, scenarios are designed by combining the present-day ESLs scenarios with the SLR scenarios. Therefore, the present-day ESL return levels are combined with three medium confidence SLR scenarios (precisely SSP1-2.6, SSP2-4.5 and SSP5-8.5) for the years 2030, 2050 and 2100. This combination follows (3) as defined in the approach overview section.

3. Results

3.1. Extreme Storm Tide and Wave Run-Up

Applying the GPD-POT method to the modelled time series of storm tide and wave run-up helps to produce extreme return levels for different return periods. Figure 7 shows the return value plots for the order of 1- to 100-year return period for each of the selected neighborhoods.

Figure 7a displays the return periods and values of an extreme storm tide. The analysis of these results indicates that storm tide hardly varies from one neighborhood to the other. The 1-year to 100-year return levels at Petit-Bassam and Sogefiha are similar, while those of Vridi differs from the two others by only 2 cm. For each location, as the return periods increase, the best-fit lines converge to the asymptotic values of 0.695 m (Vridi) and 0.716 m (Petit-Bassam and Sogefiha). The differences between the 1-year and the 100-year return levels are, therefore, all very small (1 cm). These slight variations in the extreme storm tide magnitudes could be explained firstly by the proximity of the selected neighborhoods and secondly by the fact that the meteorological variable in the storm tide (i.e., storm surge) does not have a significant influence on the total storm tide level. The storm surge contributes about 15% while the tide component is approximately 85%, making extreme storm tide events occur during high spring tide. This is particularly understandable given that the study area is recognized as a storm-free [72] and micro-tidal environment [32,73].

Figure 7b presents the EVA results for the wave run-up component. It shows that, unlike the storm tide, wave run-up variations between locations are well noticeable. The return values of extreme wave run-up obtained for Petit-Bassam are 35 to 70 cm higher than those of Vridi and Sogefiha. The difference between the two later neighborhoods is estimated to be around 9–10 cm. The absolute differences between the 1-year and 100-year return levels at each location are 0.97 m for Vridi and Sogefiha and 1.11 m for Petit-Bassam. These variations in the extreme wave run-up conditions are mainly due to the heterogeneity of beach morphological features observed on this coast [47].

From Figure 7, it is clear that the wave run-up contributes more to the total ESL as compared to the storm tide levels. For example, at Vridi, the 100-year return level of wave run-up is estimated as 3.69 m against 0.69 m for the storm tide (Figure 7a,b). Additionally, it can be seen that over the past 40 years, storm tide was quite stable in terms of temporal variability, whereas extreme wave run-up underwent both seasonal and inter-annual variations. The maximum run-up heights occurred during the meteorological event of 3 and 4 September 2016. At this particular event, the calculated run-up gave 3.61 m, 4.19 m and 3.72 m, respectively, for Vridi, Petit-Bassam and Sogefiha. About 37 occurrences of extreme wave run-up ≥ 1-year return level were identified within Port-Bouët Bay from 1981 to 2020. Among these, six were ≥ 10-year, and three were ≥ 25-year return levels. The seasonal distribution of these extreme wave events (Figure 8) shows that extreme run-up ≥ 1-year return levels occur from April to November and are most frequent during the month of August. During the assessment, no extreme wave run-up ≥ 1-year return level was observed in January and December.

In terms of the inter-annual distribution of extreme wave run-up occurrences, it has been that there has been an intensification of extreme wave conditions in Port-Bouët Bay. Both the magnitude and the number of wave run-up occurrences have increased in these recent decades. As an illustration of this, Figure 9 displays the clustering of peaks exceeding the 1-year and 10-year return levels in the neighborhood of Sogefiha. The 1-year and 10-year return levels correspond to 2.81 m and 3.35 m, respectively, as shown in Figure 7b. It can be observed in Figure 9 that approximately 75% of the extreme run-up events ≥ 1-year return level and 83% of those ≥ 10-year return levels occurred during the past two decades.

3.2. Present-Day ESLs

This section presents the outcome of the joint probability analysis. To get the present-day ESLs, extreme wave run-up was combined with extreme storm tide, assuming independence between them. Multiple combinations of these two components could be envisaged for the same ESL joint return period. For example, under the assumption of component independence, a 100-year joint return period of ESL can be obtained from the following combinations of individual return periods: 100 years and 1 years, 50 years and 2 years, 20 years and 5 years, 10 years and 10 years, with each combination having the same chance of occurring. However, only the worst-case results (the combinations that offer the maximum joint return values) are significant to the present study. They are provided with their 95% confidence intervals (in brackets) in

Table 3. Present-day ESLs deduced from the Joint Probability Method (in meters).

Joint Return Period	Vridi	Petit-Bassam	Sogefiha
1 year	3.40 [3.34–3.46]	3.86 [3.78–3.92]	3.51 [3.44–3.56]
5 year	3.78 [3.65–3.91]	4.30 [4.15–3.92]	3.90 [3.76–4.02]
10 year	3.93 [3.74–4.12]	4.47 [4.27–4.70]	4.05 [3.87–4.25]
25 year	4.12 [3.84–4.43]	4.68 [4.39–5.07]	4.23 [3.97–4.57]
50 year	4.25 [3.89–4.68]	4.83 [4.45–5.36]	4.36 [4.02–4.83]
100 year	4.37 [3.95–4.96]	4.97 [4.51–5.66]	4.48 [4.08–5.09]

. These maximum joint return values all correspond to the combination of the highest contribution of wave run-up with the smallest contribution of Storm tide in each case. The current ESLs ≥ 4 m are reached from the joint return periods of 12 years, 2 years and 8 years, respectively, in Vridi, Petit-Bassam and Sogefiha. Given the local topographic conditions, these ESL ≥ 4 m, are susceptible to cause flooding in some areas of the coast. It is important to mention that the minimal fluctuations in storm tide imply that the changes in present-day ESLs between the different neighborhoods are similar to the results obtained from the analysis of extreme wave run-up.

3.3. Future ESLs

The future ESLs were determined by combining the projected regional sea-level rise (SLR) with the present-day ESLs as stated in the approach overview. However, due to the use of this regional SLR projection, differentiation of SLR at the local scale for Vridi, Petit-Bassam and Sogefiha was not possible. Figure 10 shows calculation results for future 100-year return levels under different climate change scenarios.

Additional outcomes regarding other future ESL return levels are provided in Table S1 as Supplementary Materials. All these projections illustrate how ESLs will be impacted by climate change in the future. The results show that the ESL elevations will increase in the future as the climate keeps changing. The average future ESL increase will be 0.11 m in the next decade. This increase could reach 0.22 to 0.26 m by 2050. By the end of the century, it is predicted that the future ESL elevations will be different from those of 2020 by 0.49, 0.62 and 0.84 m on average, respectively, under SSP1-2.6, SSP2-4.5 and SSP4-8.5 scenarios. The end-century 100-year ESL at Vridi, Petit-Bassam and Sogefiha under the high emissions scenario would be 5.21 m, 5.81 m and 5.32 m, respectively.

While it is interesting to demonstrate the impact of SLR on future ESL elevations, relying solely on such data for coastal management and planning without considering the potential changes in ESL frequencies could be limiting. Therefore, we analyzed changes over time and under SSPs in the current 100-year return period, as illustrated in Figure 11. This figure shows that in Port-Bouët Bay, the current 100-year return period will decrease in all neighborhoods. Over the next decade, today’s 100-year ESLs will have return periods ranging from 52 to 63 years. This could decrease to between 24 and 37 years, i.e., more than doubling the current frequency by mid-century. By 2100, the current once-a-century ESLs in Port-Bouët Bay will occur nearly every 2 years under SSP5-8.5, and between 8 to 12 years under SSP1-2.6 and 5 to 7 years under SSP2-4.5. The neighborhood of Vridi is likely to be the most vulnerable in terms of future ESL frequency, followed by Sogefiha and Petit-Bassam in that order. This analysis allows understanding that even minor, short-term changes in SLR could significantly affect ESL frequencies. By the end of the century, ESLs that are currently infrequent in Port-Bouët Bay will become relatively common, regardless of the climate change scenario.

4. Discussion

This study used a step-by-step methodological approach that includes different modeling techniques and statistical methods to assess the spatial and temporal evolution of extreme sea levels (ESLs) along the coastline of Port-Bouët Bay. The findings demonstrate that present-day ESLs are already high and will worsen in frequency and intensity due to the impacts of climate change in the future. This implies a greater likelihood of increasing coastal hazards in the area. The ESL scenarios generated by this study provide crucial data for evaluating the exposure and risks associated with coastal flooding and devising strategies for adaptation.

The current research represents an initial attempt to study local-scale ESLs along the Ivorian coastline. Compared to previous global-scale studies by Vousdoukas et al. [74] and Kerezci et al. [75], our findings estimate much greater ESLs. The main difference in the determination of ESLs between these studies is the wave contribution, which was accounted for through the wave setup parameter in [74,75], and is about 2.5 times smaller than the wave run-up on reflective coasts [33]. This study’s estimates of extreme storm tide elevation are in agreement with those of the global-scaled study by Muis et al. [76], with present-day extreme storm tide being between 0.5 and 1 m along the coast of Port-Bouët. However, unlike the global-scale approach, the adopted local-scale approach enables the examination of spatial variations in ESLs within the bay, with a resolution of approximately 2 km.

In the West African context, similar studies by Cissé et al. [9,10] in Dakar and Saint Louis in Senegal have shown minimal contributions of storm surge components to total ESLs. This weak storm surge contribution to ESLs is due to easterly winds carrying storms away from land [53], which makes Atlantic tropical storms rarely make landfalls in this region of the globe. Instead, offshore meteorological conditions usually generate a strong swell climate [53,77], whose contributions to ESL and devastating impacts have been extensively highlighted in several studies across the region [9,10,20,78,79]. Our study also highlights an increase in ESL frequencies in Port-Bouët Bay, which is consistent with global-scale findings [24,25] but provides greater insight into the local variability for better adaptation.

In general, the accuracy of results in ESL studies may depend on data choices and assumptions in the methodology, which need to be discussed to raise uncertainties and/or find directions for further research. In this work, data choices and assumptions are mainly related to the selection of the dataset length for extreme value analysis (EVA), the determination of the representative beach slopes, and assumptions of independence between ESL components.

• Selection of dataset length

Indeed, the choice of the length of the dataset to be considered in the EVA is decisive for the consistency of results. In general, analyses conducted using datasets that are too short yield less accurate return level estimates than those using longer datasets. However, there is no practical standard for the exact length of data to be used. According to Roberts et al. [33], using at least ten years of data could provide a reliable basis for estimating up to 100-year return values via the GPD-POT method. Arns et al. [80] are more in favor of using 30 years of data for estimating 100-year return values. Nevertheless, based on his data length sensitivity analysis, Watson [14] recommended the use of a much longer dataset if available, arguing that a longer dataset is more likely to contain key information on the major extreme events that have affected the study area, thus allowing for more optimal estimates. In light of the above, the 40 years of modeled datasets used in this work can be considered sufficient for the range of estimations performed, although it is expected that accuracy will improve as more quality in situ measurements become available.

• Determination of representative beach slopes

In this study, it was observed that ESLs (Extreme Sea Levels) are very sensitive to beach morphology (beach slope). The fact of considering the beach face for the determination of the slopes, as suggested by Ruggiero et al. [61] and Stockdon et al. [58], allows for obtaining relatively high slopes that were stable and realistic for extreme conditions. Additionally, the use of historical cross-shore beach profiles helped to account for the seasonal variations of slopes. However, the irregularity of the collected profiles and their availability over a relatively short period could have an impact on the selection of representative slopes and, thus, on the quality of the estimates. It would therefore be necessary to have regular follow-ups of the morphology of the beaches in Port-Bouët.

• Independence between ESL components

This work considers two assumptions of independence between ESL components. The first assumption is the independence between SLR and present-day ESL parameters, meaning that the current ESL components will be kept constant in the future. The second assumption is the independence between storm tide and wave run-up, i.e., the occurrence of an extreme storm tide event will have no effect on the chance of extreme wave run-up occurring [71]. While these assumptions simplify our understanding of ESL, it is crucial to note that they are not always valid. Recent studies have shown that the predicted SLR could significantly influence each of the parameters composing present-day ESL [81,82,83]. Regarding the independence between storm tide and wave run-up, while the tide component of the storm tide is of astronomical origin and the wave is of meteorological origin, and hence independent, it is essential to note that the storm surge component of the storm tide may be dependent on the wave, as both are driven by the same meteorological processes. Therefore, further studies should investigate these dependence aspects to generate more practical ESL scenarios.

5. Conclusions

The present paper focused on determining the evolution of extreme sea level elevation at the three neighborhoods located along Port-Bouët Bay. The analysis is based on validated modelled historical data from 1981 to 2020, which was built within the framework of this study to overcome the limitations related to data scarcity in this study area. These modelled datasets served to estimate the frequency and magnitude of present-day ESLs. The recent IPCC projections of sea level rise were also integrated to provide future ESLs till the end of the century.

The study results suggest that ESLs have been more intense and frequent these recent years, making the estimated present-day ESL to be relatively high. It was also observed that the current elevations will keep increasing in the future as the mean sea level continue to rise under climate change effects. The analysis also permits us to understand the role played by local coastal morphologies (beach slopes) on the spatial distribution of ESL, as demonstrated by the steepest slopes having the largest ESL elevations and the lowest slopes having the lowest ESL estimates.

These findings were possible by the application of two complementary statistical methods: the Generalized Pareto Distribution fitted to the peaks-over-thresholds approach (GPD-POT) and the joint probability method (JPM). The GPD-POT helped to estimate return values associated with return periods of individual ESL components. The JPM were used in addition to combining the individual likelihood of occurrences to obtain their joint return period/values.

The analysis and results presented in this study provide valuable insights into the current and growing threat of coastal inundation resulting from the increasing magnitude of ESLs along the neighborhoods of Port-Bouët Bay. The obtained ESL information can be useful in decision-making, planning, and emergency response efforts related to coastal flooding and others impacts. As such, these estimates can benefit a wide range of stakeholders, including coastal communities, infrastructure planners, environmental organizations, and researchers. Coastal communities can use this information to prepare for and mitigate the impacts of ESLs, while engineers can use the data to design resilient coastal protection structures. Environmental organizations can advocate for policies that address sea level rise, and researchers can make informed decisions.

Nevertheless, we strongly encourage comprehensive assessments of coastal risks associated with these ESL scenarios, particularly those related to coastal inundation and proactive actions to manage these risks in Port-Bouët Bay. Lastly, we suggest updating these projections when more consistent and accurate in situ data become available.