the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
| 29 Sep 2026
Extreme high-water anomalies in southwest coastal Bangladesh: model uncertainty, cyclone-event checks and non-stationarity
Abstract. Extreme high-water estimates are needed for coastal flood-risk assessment and infrastructure design, but short, incomplete and strongly seasonal tidal-river records make them difficult to constrain. We analysed daily maximum water levels through 2025 at five Bangladesh Water Development Board stations in southwest coastal Bangladesh using a station-specific leave-one-year-out Monthly Climatological High-Water Reference (MCHR) derived from monthly means of daily maxima. Annual maxima were modelled with Normal, Gumbel-I and generalized extreme-value distributions fitted by L-moments (GEV-LM), while declustered peaks over threshold were analysed with a generalized Pareto distribution (POT-GPD). After 90 % annual completeness screening, 25–34 annual maxima remained per station. At the 100-year return period, the four stationary approaches differed by 15.6–21.5 %. POT estimates were stable across threshold and declustering choices at four stations, whereas SW244 showed an 18.75 % maximum deviation from its baseline estimate. Combined trend, change-point and likelihood evidence supported time-varying GEV location at four stations after false-discovery-rate adjustment and trend-free prewhitening sensitivity; likelihood evidence remained significant at all four, while SW243 remained stationary. Replacing the mean MCHR with a median reference changed 100-year GEV-LM and POT-GPD estimates by at most 2.79 % and 3.48 %, respectively. Cyclone-date checks showed that Aila, Amphan, Yaas and Remal coincided with annual maxima or unusually high anomalies at multiple stations. Extreme high-water estimates from short tidal-river records are most informative when model choice, sampling, non-stationarity, parameter uncertainty and data quality are reported alongside the point estimate instead of being collapsed into a single deterministic return level.
- Preprint
(2305 KB) - Metadata XML
-
Supplement
(436 KB) - BibTeX
- EndNote
Status: final response (author comments only)
-
RC1: 'Comment on egusphere-2026-5598', Anonymous Referee #1, 27 Oct 2026
-
AC1: 'Reply on RC1', Sujoy Biswas, 31 Oct 2026
Title: Extreme high-water anomalies in southwest coastal Bangladesh: model uncertainty, cyclone-event checks and non-stationarity
Dear Anonymous Referee #1,
We thank the reviewer for the careful evaluation of our manuscript and for the constructive comments. The suggestions helped us identify several aspects of the methodology and interpretation that required further clarification. We have revised the manuscript and Supplementary Material accordingly. The main revisions include clarification of the data quality-control procedure, additional sensitivity analyses, improved statistical diagnostics, and a more careful discussion of uncertainty and physical interpretation.
Our point-by-point responses are provided below.
Major Comments
Comment 1
Reviewer comment:
Lines 48–51: The limitation of conventional extreme value analysis is somewhat solved by the introduction of MEVD (Marani and Ignaccolo, 2015; Caruso and Marani, 2022), which uses ordinary events to model extremes, which is not discussed in the paper.
Response:
We agree that MEVD is relevant to the estimation of extreme water levels, particularly when information from ordinary events is available. Unlike conventional AMS-based approaches, MEVD uses information about ordinary events and their occurrence within a given period to estimate extreme values.
Our analysis was designed to compare AMS-based distributions and POT-GPD using daily-maximum water-level records from five stations. The main focus was to examine how model choice, sampling decisions, non-stationarity and data quality influence long-return-period estimates.
We have therefore added a discussion of MEVD in the Introduction, including its relevance to coastal extreme-water-level estimation and its comparison with conventional AMS and POT-GPD approaches. We have also clarified that MEVD was not fitted in the present study because it was outside the scope of the selected analysis.
A direct comparison of MEVD, AMS and POT-GPD using the same records and appropriate validation is identified as a direction for future research.
Changes made: Introduction, Section 1.
Comment 2
Reviewer comment:
Lines 98–99: Please explain the cleaning of raw data, especially what is meant by "suspicious observations"? How did you decide whether to keep or remove the observation?
Response:
We agree that the original description of data cleaning was not sufficiently clear. We have rewritten the quality-control (QC) procedure to explain how observations were flagged and evaluated.
The records were first sorted chronologically and checked for duplicate dates, non-numeric entries, repeated values requiring inspection, and isolated observations that differed substantially from neighbouring observations or the wider station record. These checks were used to flag observations for further review; they did not lead to automatic deletion.
Flagged high values were checked against adjacent daily observations and a separate comparison series obtained from the Water Resources Planning Organization (WARPO), based on BWDB records. The only isolated high value requiring a manual cross-check was the SW244 observation on 26 June 1996. The original series reports 5.88 m, whereas the WARPO comparison series reports 3.39 m.
As no original BWDB correction sheet or station-history record was available to establish 3.39 m as an official correction, we now describe it as a comparison-based QC value. The primary analysis uses 3.39 m, while the original 5.88 m is retained in the audit record. We have also included sensitivity analyses using the raw value and excluding 1995–1996.
No other isolated high observation was manually changed or removed.
We have additionally clarified that the 90% annual completeness criterion was applied to AMS selection, whereas the separate 90% monthly completeness criterion was used for the MCHR reference. Missing observations were not imputed.
The revised Figure 2 also presents the analytical workflow, including the QC and completeness-screening steps.
Changes made: Section 2.1, Figure 2 and Supplementary Table S1.
Comment 3
Reviewer comment:
Section 2.4: Why are different parameter estimation methods used for different models? This makes comparative analysis less reliable.
Response:
We agree that using different parameter-estimation methods can influence the comparison among the models.
In the primary analysis, Normal parameters were estimated using the sample mean and standard deviation, Gumbel-I parameters using the method of moments, GEV parameters using probability-weighted moments/L-moments, and POT-GPD parameters using the method of moments applied to declustered exceedances.
These methods were selected considering the model structure and the relatively short AMS records, which contain only 25–34 annual maxima per station. In particular, L-moments were used for GEV because of their established application to extreme-value distributions with relatively small samples.
However, we acknowledge that this approach does not provide a controlled comparison of distribution forms alone. The primary four-method comparison also reflects estimator choice and, for POT-GPD, a different sampling framework.
To examine estimator sensitivity, we added a separate common-MLE comparison for the three AMS models using the same accepted annual-maximum observations at each station. The results are reported in Supplementary Table S11.
We have retained the original estimation methods for the primary analysis and clarified that the resulting spread should be interpreted as a combined model, estimator and sampling sensitivity measure rather than as evidence of distribution choice alone. POT-GPD is evaluated separately because it uses threshold exceedances rather than annual maxima.
Changes made: Section 2.4, Section 4.1 and Supplementary Table S11.
Comment 4
Reviewer comment:
Lines 221–222: Justify the use of different numbers of bootstrap samples for different models for uncertainty analysis?
Response:
We agree that the number of bootstrap replicates needs to be justified by examining the stability of the resulting uncertainty intervals.
We therefore compared the 95% interval endpoints at increasing replicate counts for the GEV models.
For the stationary GEV-MLE model (M0), the upper interval endpoint at SW244 changed from 3.675 m with 1,000 resamples to 4.075 m with 5,000 resamples. At 3,000 resamples, the endpoint was 4.068 m. We therefore increased the M0 bootstrap to 5,000 non-parametric resamples for all stations.
For the non-stationary GEV-MLE model (M1), the largest change in either interval endpoint between 1,000 and 2,000 replicates was 0.035 m. Accordingly, 2,000 parametric bootstrap samples were retained for M1.
POT-GPD retained 2,000 non-parametric resamples of the retained exceedances as a practical setting. However, a formal replicate-count convergence test was not performed for POT-GPD. We have explicitly acknowledged this limitation in the revised Supplementary Material.
A separate Poisson-count sensitivity analysis was also retained to account for uncertainty in the number of POT peaks. This is an event-rate sensitivity analysis rather than a replicate-count convergence test.
The revised uncertainty intervals and convergence assessment are presented in Figure 9 and Supplementary Tables S9 and S10.
Changes made: Section 2.6, Section 3.6, Figure 9 and Supplementary Tables S9–S10.
Comment 5
Reviewer comment:
Lines 261–262: Try to explain the physical reason for the evident episodic departures.
Response:
We have expanded the explanation of the episodic positive departures in the revised manuscript.
These departures may reflect a combination of storm-related winds and atmospheric pressure, unusual tidal conditions, tide–surge interaction, rainfall, river discharge and local estuarine hydraulics. Their relative contributions may differ among stations and individual events.
However, the available records consist of daily maximum water levels. Continuous predicted astronomical-tide data, event-specific meteorological forcing and river-flow observations are not available for the full analysis period. Therefore, these physical mechanisms cannot be separated for individual departures.
We have consequently avoided assigning specific physical causes to the observed anomalies. Instead, we interpret them as seasonally adjusted high-water departures and have clarified this limitation in the Results section.
Changes made: Section 3.2 and Section 4.2.
Comment 6
Reviewer comment:
Lines 295–303: It is not clear how the statistical metrics are calculated. What data is used to compute the metrics? Does it use bootstrap results?
Response:
We have clarified the calculation of the statistical metrics and the observations used in the analysis.
For Normal, Gumbel-I and GEV-LM, the accepted annual maxima were sorted in ascending order and assigned Gringorten plotting positions:
Fᵢ = (i − 0.44)/(n + 0.12),
where i is the rank and n is the number of accepted annual maxima.
The fitted model quantiles were evaluated at the same empirical probabilities. RMSE and MAE were then calculated by comparing the observed ordered values with the corresponding fitted quantiles.
For POT-GPD, Q–Q RMSE was calculated separately using the empirical and fitted GPD quantiles of the retained declustered exceedances at matched plotting positions. Because AMS and POT-GPD use different samples, their goodness-of-fit metrics are not directly ranked against each other.
Importantly, these metrics were calculated from the original accepted observations, not from bootstrap samples. Bootstrap resampling was used separately to estimate the 95% uncertainty intervals of the return levels.
We have also clarified that the multi-model/method spread was calculated from the estimated return levels using the maximum and minimum estimates relative to their mean. It was not calculated from bootstrap replicates.
Changes made: Section 2.5 and Section 3.3.
Comment 7
Reviewer comment:
Also, GEV-LM achieved the lowest average RMSE and MAE. This result could be influenced by the parameter estimation methods. Have you tried all models with the same parameter estimation technique?
Response:
We agree that the parameter-estimation method may influence the observed model performance. To examine this, we fitted Normal, Gumbel-I and GEV using maximum likelihood estimation (MLE) with the same accepted AMS observations and Gringorten plotting positions at each station.
The common-MLE comparison also showed that GEV-MLE had the lowest average in-sample RMSE and MAE, with values of 74.9 mm and 60.9 mm, respectively. The corresponding values were 85.5 mm and 65.5 mm for Normal-MLE and 97.7 mm and 73.8 mm for Gumbel-I-MLE.
The best-performing model nevertheless varied among individual stations. These results indicate that the competitive average performance of GEV is not explained solely by the use of L-moments.
POT-GPD was not included in this common-MLE comparison because it uses threshold exceedances rather than annual maxima. We have retained this distinction in the revised manuscript.
We also emphasize that these are in-sample goodness-of-fit results and do not establish superior prediction of long return periods.
Changes made: Section 3.3 and Supplementary Table S11.
Comment 8
Reviewer comment:
It is apparent from Fig. 4 that Gumbel-I (orange line) gives the highest estimates. Not sure about the x-limit; maybe 100 years. In contrast, in Fig. 6, at SW243, POT-GPD provides the highest return-level estimates. This discrepancy should be addressed, and the description should be modified accordingly (lines 312–320).
Response:
We rechecked Figures 4 and 6 together with the numerical values reported in Table 3.
The two figures represent different aspects of the analysis. Figure 4 compares the empirical annual maxima with the fitted Normal, Gumbel-I and GEV-LM distributions at common Gringorten plotting positions. Its horizontal axis represents the empirical return-period range of the accepted AMS observations. It is therefore an observed-range goodness-of-fit assessment, not a comparison of 100-year extrapolated estimates.
Figure 6 presents the stationary return-level estimates from all four approaches, including POT-GPD, at 10-, 25-, 50- and 100-year return periods.
The ranking of the estimates depends on the station and return period. In particular, the 100-year estimates in Table 3 show that Gumbel-I gives the highest estimate at SW243 (2.285 m), whereas POT-GPD gives 2.002 m. Across the five stations, POT-GPD gives the highest 100-year estimate only at SW107.2 (2.129 m); Gumbel-I gives the highest estimate at the other four stations.
We have clarified the distinction between observed-range model fit and long-return-period extrapolation in the Results section. The corresponding figure caption and description have also been revised to avoid suggesting that the model ranking is the same across all return periods.
Changes made: Section 3.3, Section 3.4, Figures 4 and 6, and Table 3.
Comment 9
Reviewer comment:
Lines 318–320 state the influence of an individual extreme on tail behavior under data-scarce settings. However, Gumbel-I uses AMS, and POT-GPD uses peaks over a threshold, which picks more peaks. Whether peaks selected for POT-GPD involve the specific peak on 18 April 1996 or are filtered out? The predictive performance of the two models should also be evaluated based on these factors.
Response:
We checked the SW244 AMS and POT-GPD samples to confirm whether the 18 April 1996 observation was retained.
The 18 April 1996 anomaly (2.080819 m) is the accepted annual maximum for 1996 and is also retained as the fourth peak in the baseline POT-GPD sample. The baseline POT analysis uses the 95th-percentile threshold (u = 0.677901 m) and 3-day declustering. Therefore, the observation was not filtered out of either analysis.
We then examined its influence through targeted sensitivity analyses.
Removing the 1996 annual maximum reduced the Gumbel-I 100-year estimate from 2.1364 m to 1.8194 m.
For POT-GPD, removing only the 18 April peak reduced the 100-year estimate from 1.9058 m to 1.4634 m. Removing all four retained POT peaks from 1996 reduced the estimate to 1.3922 m.
We also examined the effect of holding out these observations. When the 1996 AMS maximum was excluded, the fitted Gumbel-I model implied a return period of approximately 274 years for the held-out observation. When all four 1996 POT peaks were removed, the fitted GPD upper endpoint was approximately 1.869 m, below the observed 18 April anomaly.
These results show that the 1996 observation and associated event cluster have a substantial influence on the fitted upper tail, particularly for POT-GPD. However, these are targeted sensitivity diagnostics rather than a full independent predictive-validation study. We have clarified this distinction in the revised manuscript and have not used these results to claim that either model has superior predictive performance.
The original observation remains in the primary analysis because its influence on the estimated tail does not, by itself, establish that the observation is invalid.
Changes made: Section 3.4 and Supplementary Tables S12a–S12b.
Comment 10
Reviewer comment:
Of course, the threshold for POT-GPD varies across stations, as evidenced for SW244 station. What could the physical drivers possibly be contributing to this large deviation?
Response:
The POT threshold was selected separately for each station using the 95th percentile of its anomaly series. Therefore, differences in threshold values reflect differences in the station-specific anomaly distributions.
The sensitivity analysis showed that SW244 behaved differently from the other stations. Across the 12 threshold and declustering combinations, its 100-year estimates ranged from 1.765 m to 2.263 m, giving a maximum absolute deviation of 18.75% from the baseline estimate of 1.906 m. The other four stations showed considerably smaller deviations.
Local tidal amplification or damping, tide–surge interaction, channel geometry, freshwater discharge, rainfall, cyclone exposure and estuarine morphology may contribute to differences in high-water behaviour among nearby stations.
However, the available daily-maximum records do not allow us to isolate the contribution of any individual physical mechanism. We have therefore clarified that the threshold is a statistical cut-off rather than a direct measure of a physical driver. The large deviation at SW244 is interpreted as evidence of sampling and tail sensitivity, not proof of a specific physical cause.
Changes made: Section 3.5 and Section 4.2; Supplementary Table S7 and Figure S1.
Comment 11
Reviewer comment:
From Fig. 9, it is clear that POT-GPD estimates are less uncertain than the other two approaches. The unusually wide uncertainty band at SW244 should be better explained.
Response:
We have expanded the discussion of uncertainty at SW244 to distinguish between parameter uncertainty and sensitivity to sampling choices.
At SW244, the stationary GEV-MLE model (M0) gives a 100-year estimate of 2.005 m with a 95% interval of 1.153–4.075 m. The non-stationary GEV-MLE model (M1) gives a conditional estimate of 2.499 m with an interval of 1.524–6.312 m. The event-rate-aware POT-GPD estimate is 1.906 m, with an interval of 1.286–2.557 m.
The M1 interval is particularly wide because SW244 has only 25 accepted annual maxima, and the fitted GEV tail parameters and temporal trend are uncertain. The 100-year conditional quantile also extends beyond the available observation period. M0 shows a similar sensitivity to bootstrap replicate count, as reflected in the convergence assessment.
Although the POT-GPD interval is narrower than the two GEV intervals, it is wider at SW244 than at the other stations. This is consistent with the greater sensitivity of the SW244 POT estimates to threshold and declustering choices, as well as uncertainty in the fitted tail and the number of retained peaks.
We have also clarified that the event-rate-aware POT interval accounts for uncertainty in the number of retained peaks. The narrower POT interval should not be interpreted as evidence of superior predictive performance because AMS and POT-GPD use different samples and resampling assumptions.
Changes made: Section 3.6, Figure 9 and Supplementary Tables S9–S10.
Minor Comments
Comment 1
Reviewer comment:
Line 21: Do you mean sea-level rise can create an increase in the frequency of storm surges, or does it increase the magnitude of storm surges? The sentence lacks clarity.
Response:
We have revised the sentence to clarify the role of sea-level rise. Sea-level rise raises the background coastal water level. Consequently, a storm surge of the same meteorological magnitude can produce a higher total water level and greater inundation hazard.
The revised statement does not imply that sea-level rise necessarily increases the frequency or meteorological magnitude of storm surges.
Changes made: Introduction, Section 1.
Comment 2
Reviewer comment:
Line 209: What are M0 and M1? Models?
Response:
We have clarified the definitions of M0 and M1 in the Methods section.
M0 represents the stationary GEV-MLE model, in which location, scale and shape remain constant over time. M1 represents the non-stationary GEV-MLE model, in which only the location parameter varies linearly with centred calendar year, while scale and shape remain constant.
These definitions are now provided before the presentation of the stationarity results.
Changes made: Section 2.5 and Section 3.6.
Comment 3
Reviewer comment:
Line 228: What is the reviewer summary workbook? Do you mean supplementary materials?
Response:
Yes, the reference was intended to describe the Supplementary Material and accompanying reproducibility files. We have removed the ambiguous expression and now refer directly to the Supplementary Material and reproducibility files.
We have also clarified that the verification script checks the reported derived results against the distributed materials; it does not reconstruct the complete analysis from restricted raw BWDB observations.
Changes made: Section 2.6 and Code availability.
Comment 4
Reviewer comment:
Change the word ‘Missingness’; use technical English.
Response:
We have replaced "missingness" with "proportion of missing observations" or "missing observations", as appropriate, to maintain consistent technical terminology.
Changes made: Sections 2.1 and 3.1.
Comment 5
Reviewer comment:
The Q-Q agreement is evaluated using R², which is not the best metric. It would be better to use RMSE.
Response:
We agree that RMSE provides a more direct measure of the difference between empirical and fitted Q–Q quantiles.
We have therefore replaced the previous GPD Q–Q R² diagnostic with Q–Q RMSE in Table 2 and revised Figure 5(c) accordingly.
The resulting Q–Q RMSE ranges from 15.5 mm at SW20 to 38.8 mm at SW244. This indicates that SW20 has the lowest discrepancy between empirical and fitted GPD quantiles, whereas SW244 has the largest discrepancy.
The corresponding Results text and figure caption have also been revised to maintain consistency.
Changes made: Section 3.1, Section 3.3, Table 2 and Figure 5(c).
Comment 6
Reviewer comment:
The conclusion needs to be modified. It should properly convey the main takeaway from this study. Also, add future directions.
Response:
We have revised the Conclusion to provide a clearer summary of the main findings and their implications.
The revised conclusion emphasizes that uncertainty in extreme high-water anomaly estimates arises from several linked sources, including model choice, sample definition, temporal change, parameter uncertainty and data quality. The four stationary approaches differ by 15.6–21.5% at the 100-year return period, and no single model provides the highest estimate at every station.
We have also clarified the station-dependent evidence of non-stationarity. The combined trend, change-point and likelihood results support a time-varying GEV location at SW107.2, SW20, SW241 and SW244, whereas the stationary model remains the primary choice at SW243.
Particular attention is given to SW244 because of its short AMS record, relatively large POT Q–Q RMSE, high threshold and declustering sensitivity, wide uncertainty intervals and unresolved data-history uncertainty.
The conclusion also clarifies that the estimated return levels represent seasonally adjusted high-water anomalies, not conventional storm-surge residuals, flood depths or absolute engineering design elevations.
Finally, we have included future research directions involving a comparison of MEVD with AMS and POT-GPD using suitable out-of-sample validation, longer and higher-frequency water-level records with predicted astronomical tides, authoritative gauge-history records, and integration of extreme-value analysis with hydrodynamic inundation modelling and population and infrastructure exposure data.
We believe these revisions provide a more balanced interpretation of the findings and clarify both the potential engineering relevance and the limitations of the present analysis.
Changes made: Section 5, Conclusions.
We thank the reviewer for the constructive comments. We believe that addressing these points has improved the clarity of the methodology, the interpretation of the results and the overall presentation of the manuscript.
Citation: https://doi.org/10.5194/egusphere-2026-5598-AC1
-
AC1: 'Reply on RC1', Sujoy Biswas, 31 Oct 2026
-
RC2: 'Comment on egusphere-2026-5598', Anonymous Referee #2, 27 Oct 2026
Review of “Extreme high-water anomalies in southwest coastal Bangladesh: model uncertainty, cyclone-event checks and non-stationarity”egusphere-2026-5598Comments:
- Explain the selection of the stations. Why was this choice made?
- The quality control (lines 98 to 105) is a bit confusing. Please revise this section and include a figure to clearly explain the process.
- L107: introduce QC
- L106 to 110 are very difficult to understand
- Too many references are out of date. Please update the bibliography.
- L231: Again, explain the selection of the cyclones’ dates.
- The paragraph with lines 250 to 254 is really complex to understand.
- According to Figure 3, there are several years discarded because of the absence of valid data. Therefore, the series is shorter than expected from the previous sections. You need to clarify this issue.
- Figure 7: The y-axis should have limits of 1 to 3 to allow seeing the complete behavior.
- I assume that “M0 retained constant location, scale and shape, whereas M1 allowed only the location parameter to vary linearly with centred time.” is the definition of M0 and M1; please confirm and clarify in the text (also in the results section)
- I can’t understand why you selected a Cyclone with unavailable data.
- The section regarding Cyclone analysis is poor and should be improved or removed altogether, because it does not provide relevant, useful information.
- I expect a deeper evaluation of the heights in terms of the consequences for the population and the infrastructure.
Citation: https://doi.org/10.5194/egusphere-2026-5598-RC2 -
AC2: 'Reply on RC2', Sujoy Biswas, 31 Oct 2026
Title: Extreme high-water anomalies in southwest coastal Bangladesh: model uncertainty, cyclone-event checks and non-stationarity
Dear Anonymous Referee #2,
We thank the reviewer for the careful assessment of our manuscript and for the constructive suggestions. The comments helped us identify several areas where the methodology, data-processing procedures, figures and interpretation needed further clarification. We have revised the manuscript accordingly, with particular attention to the station-selection rationale, quality-control procedure, cyclone-event analysis, completeness screening and engineering implications.
Our point-by-point responses are provided below.
Comment 1
Reviewer comment:
Explain the selection of the stations. Why was this choice made?
Response:
We selected the five stations primarily because suitable daily-maximum water-level records were available for the analysis period and because they cover three major tidal-river systems in southwest Bangladesh.
SW107.2 (Rayenda) represents the Baleswar-Haringhata system, SW20 (Amtali) represents the Barisal-Buriswar system, and SW241 (Khulna), SW243 (Chalna) and SW244 (Mongla) represent three locations within the Rupsa-Pasur system.
This selection allows us to examine differences between tidal-river systems as well as variations among neighbouring stations within the Rupsa-Pasur system. We have clarified the selection criteria in the revised manuscript and also acknowledge that these five stations do not constitute a statistically representative sample of all coastal gauges in Bangladesh. The analysis is therefore limited to the stations for which suitable records were available.
Changes made: Section 2.1, Study area and water-level data.
Comment 2
Reviewer comment:
The quality control (lines 98 to 105) is a bit confusing. Please revise this section and include a figure to clearly explain the process.
Response:
We agree that the original description of the quality-control procedure was not sufficiently clear. We have rewritten this section and revised Figure 2 to make the data-screening process easier to follow.
The revised procedure begins with chronological sorting of the records, followed by checks for duplicate dates, non-numeric entries, repeated values requiring inspection and isolated observations that differ substantially from neighbouring values or the wider station record.
Importantly, flagged observations were not automatically removed. They were checked against adjacent daily observations and a second comparison series compiled from BWDB records and supplied by our supervisor.
The revised description also explains how the SW244 observation on 26 June 1996 was evaluated and why a separate sensitivity analysis was necessary. In addition, we have clarified the difference between the 90% annual completeness criterion used for AMS and the 90% monthly completeness criterion used for the MCHR reference.
Figure 2 now provides a clearer overview of the complete analytical workflow, including data screening, anomaly construction, extreme-value modelling and subsequent sensitivity analyses.
Changes made: Section 2.1 and revised Figure 2.
Comment 3
Reviewer comment:
L107: introduce QC
Response:
Thank you for pointing this out. We have introduced the abbreviation at its first occurrence by writing "quality control (QC)" in the revised manuscript. The abbreviation QC is used consistently thereafter.
Changes made: Section 2.1.
Comment 4
Reviewer comment:
L106 to 110 are very difficult to understand.
Response:
We agree that the original explanation of the SW244 observation was difficult to follow. We have rewritten this paragraph to explain the issue more directly.
The source series reports 5.88 m for 26 June 1996, whereas the comparison series supplied by our supervisor reports 3.39 m for the same date. The two series agree on nearly all other dates. However, no original BWDB correction sheet or station-history record was available to establish 3.39 m as an official correction.
For this reason, the revised analysis uses 3.39 m as a comparison-based QC value while retaining the original 5.88 m in the audit record. We have also included two sensitivity scenarios: one using the raw 5.88 m observation and another excluding 1995–1996.
We have made it clear that this decision represents an unresolved data-quality uncertainty rather than proof that the original observation was physically impossible. This clarification also helps explain why SW244 requires particular caution when interpreting its extreme-value estimates.
Changes made: Sections 2.1 and 4.3; Supplementary Table S1.
Comment 5
Reviewer comment:
Too many references are out of date. Please update the bibliography.
Response:
Thank you for this recommendation. We have reviewed the literature used in the manuscript and strengthened the discussion of relevant recent developments in coastal extreme-water-level analysis.
In particular, we have improved the discussion of seasonal variability, changing storm-surge behaviour and coastal water-level processes by placing greater emphasis on relevant studies, including D’Arcy et al. (2023), Palmer et al. (2024) and Barroso et al. (2025).
We have also expanded the Introduction to discuss the metastatistical extreme-value distribution (MEVD), with reference to Marani and Ignaccolo (2015) and Caruso and Marani (2022). This provides a clearer explanation of alternative approaches to extreme-value estimation and helps establish the scope of the present AMS and POT-GPD comparison.
We have retained older references where they remain important for explaining the original statistical methods, including extreme-value distributions, trend tests and L-moment estimation. These references are still relevant to the methodological framework and were not removed simply because of their publication date.
The revised Introduction and Discussion now provide a clearer connection between the existing literature and the specific research questions addressed in this study.
Changes made: Introduction, Discussion and References.
Comment 6
Reviewer comment:
L231: Again, explain the selection of the cyclones’ dates.
Response:
We have clarified that the cyclone dates were selected from published India Meteorological Department reports rather than from peaks observed in the station records.
Five major Bay of Bengal cyclone dates were specified for the event check: Sidr (15 November 2007), Aila (25 May 2009), Amphan (20 May 2020), Yaas (26 May 2021) and Remal (26 May 2024).
For each station, we examined the largest observed high-water anomaly within one day before or after the documented cyclone date (D−1 to D+1). We also repeated the comparison using a wider D−2 to D+2 window to check whether the selected event peak was sensitive to the timing window.
An annual-maximum match was recorded only when the accepted annual maximum occurred within the specified event window. Missing observations were recorded as unavailable rather than reconstructed.
We have also clarified that this analysis examines the correspondence between independently documented cyclone dates and observed high-water anomalies; it does not independently establish the magnitude of meteorological storm surge.
Changes made: Section 2.7 and Supplementary Table S8.
Comment 7
Reviewer comment:
The paragraph with lines 250 to 254 is really complex to understand.
Response:
We agree that the original paragraph required clearer wording. We have simplified the explanation of missing observations during the cyclone-season months.
The proportion of missing daily observations during April–May and October–November differed by less than one percentage point from the full-record proportion at each station. This suggests that missing observations were not disproportionately concentrated during these months.
However, this does not mean that every cyclone event can be evaluated. An event occurring entirely during a period with missing observations cannot be assessed from the available records.
We have rewritten the paragraph to make this distinction clearer without implying that missing observations have been reconstructed or that all cyclone-related events are represented in the dataset.
Changes made: Section 3.1, Data screening and extreme samples.
Comment 8
Reviewer comment:
According to Figure 3, there are several years discarded because of the absence of valid data. Therefore, the series is shorter than expected from the previous sections. You need to clarify this issue.
Response:
Thank you for pointing out this potential source of confusion. We have clarified that the gaps shown in Figure 3 do not necessarily represent entire years excluded from the analysis.
Two separate completeness criteria were applied in this study. A year was accepted for the annual maximum series (AMS) when at least 90% of its expected daily observations were available. In contrast, an individual year-month was included in the MCHR reference calculation only when at least 90% of its daily observations were available.
Therefore, a month that failed the monthly completeness criterion was excluded from the climatological reference calculation, but this did not automatically exclude the entire year from AMS.
No missing observations were interpolated or reconstructed. The distinction between these two screening procedures is now explained in the Methods section and explicitly stated in the caption of Figure 3.
After annual completeness screening, the accepted AMS sample sizes were 30, 32, 34, 32 and 25 years for SW107.2, SW20, SW241, SW243 and SW244, respectively.
We believe this clarification resolves the apparent difference between the analysis periods and the gaps displayed in Figure 3.
Changes made: Sections 2.1, 2.2 and 3.2; revised Figure 3 caption and Table 1.
Comment 9
Reviewer comment:
Figure 7: The y-axis should have limits of 1 to 3 to allow seeing the complete behavior.
Response:
Thank you for this suggestion. We have revised Figure 7 to use a y-axis range of 1.0–3.0 m so that the baseline POT-GPD estimates and their sensitivity ranges can be viewed more clearly.
The underlying return-level estimates and sensitivity calculations remain unchanged. Only the presentation of the figure has been adjusted.
Changes made: Revised Figure 7, with a y-axis range of 1.0–3.0 m.
Comment 10
Reviewer comment:
I assume that “M0 retained constant location, scale and shape, whereas M1 allowed only the location parameter to vary linearly with centred time.” is the definition of M0 and M1; please confirm and clarify in the text (also in the results section).
Response:
Yes, this interpretation is correct. We have clarified the definitions of both models in the revised Methods and Results sections.
M0 represents the stationary GEV-MLE model, in which the location, scale and shape parameters remain constant over time.
M1 represents the non-stationary GEV-MLE model, in which only the location parameter varies linearly with centred calendar year, while the scale and shape parameters remain constant.
We have now introduced these definitions more clearly before presenting the stationarity results so that the distinction between the two models is easier to follow.
We have also clarified that the return levels calculated from M1 at the 2025 covariate value are conditional estimates used to examine the effect of non-stationarity. They should not be interpreted as forecasts of an event expected once in the next T years under a continuously changing distribution.
Changes made: Sections 2.5 and 3.6.
Comment 11
Reviewer comment:
I can’t understand why you selected a Cyclone with unavailable data.
Response:
Thank you for raising this point. Sidr was included in the original event list because its date is documented in the published cyclone records. The event dates were specified independently of the observed station peaks.
However, after checking the available records, we found that none of the five stations had observations within the specified D−1 to D+1 event window for Sidr.
We therefore excluded Sidr from Figure 10 and from the interpretation of station responses. It is retained only in Supplementary Table S8 as unavailable so that the data limitation remains transparent.
We could not establish the reason for the absence of observations within this event window, and we have not attributed it to any specific cause or attempted to reconstruct the missing data.
Changes made: Sections 2.7 and 3.7; Figure 10 and Supplementary Table S8.
Comment 12
Reviewer comment:
The section regarding Cyclone analysis is poor and should be improved or removed altogether, because it does not provide relevant, useful information.
Response:
We appreciate this comment and have reconsidered how the cyclone-event analysis should be presented. Rather than treating it as an independent validation of storm-surge magnitude, we have shortened and clarified the section to define its limited purpose.
The revised analysis examines whether independently documented cyclone dates coincide with unusually high observed anomalies at the five stations. It also allows us to describe differences in the response among neighbouring tidal-river stations.
The results show that Aila coincided with the 2009 annual maximum at SW243. Amphan and Yaas coincided with annual maxima at SW20 and SW243, while Remal coincided with annual maxima at SW107.2, SW20, SW241 and SW244. Yaas produced anomalies at or above the 99th percentile at four stations, whereas Remal produced anomalies at or above the 95th percentile at all five stations.
These results provide descriptive evidence that major cyclone periods can coincide with unusually high water, although the correspondence varies among stations.
We have also retained the D−2 to D+2 timing sensitivity check, which did not change any evaluable peak or annual-maximum match.
Importantly, we now state more explicitly that the event check and the extreme-value analysis use the same BWDB water-level observations. Therefore, this comparison is not independent validation of storm-surge magnitude, nor does it separate meteorological surge from astronomical tide or other local influences.
The revised Figure 10 presents only the four cyclones with evaluable station-event data, while Sidr is documented as unavailable in Supplementary Table S8.
We have retained this shorter analysis because it provides useful event context and illustrates spatial differences in observed high-water anomalies, while avoiding claims that cannot be supported by the available data.
Changes made: Sections 2.7, 3.7 and 4.2; revised Figure 10 and Supplementary Table S8.
Comment 13
Reviewer comment:
I expect a deeper evaluation of the heights in terms of the consequences for the population and the infrastructure.
Response:
We agree that the practical implications of the estimated high-water levels need to be explained more clearly, particularly in relation to coastal infrastructure and flood-risk assessment.
We have expanded Section 4.3 to clarify how the estimated return levels may be used and what additional information is required before they can be translated into actual flood impacts.
An important distinction is that the estimated values in this study represent high-water anomaly increments above a station-specific seasonal reference. They are not flood depths, ground-relative water levels or absolute design-water elevations.
The results can therefore help identify stations where more detailed flood-hazard assessment may be necessary. However, they cannot directly quantify the number of people or infrastructure assets exposed to flooding or estimate potential damage.
We have also added the relevant Bangladesh National Building Code (BNBC 2020) return-period context. For flood and surge loading, the code identifies 100-year values for essential and hazardous facilities and 50-year values for other structures. We clarify that the anomaly estimates from this study are useful for comparative engineering screening but should not be treated as absolute BNBC design-water or surge elevations.
A more detailed assessment of population and infrastructure consequences would require datum-consistent absolute water levels, local ground and embankment elevations, relative sea-level changes, drainage and wave information, and hydrodynamic inundation modelling linked to population and infrastructure exposure data.
We have therefore expanded the engineering discussion while avoiding unsupported estimates of inundation depth, affected population or infrastructure losses.
Changes made: Section 4.3, Engineering implications and limitations, and Section 5, Conclusions.
We thank the reviewer again for the constructive suggestions. We believe that the revisions have improved the clarity of the data-processing procedures, the presentation of the statistical analysis and the interpretation of the results. We have also clarified the limitations of the cyclone-event assessment and the extent to which the estimated high-water anomalies can be used for engineering and flood-risk screening.
Citation: https://doi.org/10.5194/egusphere-2026-5598-AC2
Viewed
| HTML | XML | Total | Supplement | BibTeX | EndNote | |
|---|---|---|---|---|---|---|
| 264 | 115 | 46 | 425 | 50 | 49 | 55 |
- HTML: 264
- PDF: 115
- XML: 46
- Total: 425
- Supplement: 50
- BibTeX: 49
- EndNote: 55
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
In this paper, titled “Extreme high-water anomalies in southwest coastal Bangladesh: model uncertainty, cyclone-event checks and non-stationarity”, the authors performed a comparative assessment of asymptotic extreme value analysis approaches to evaluate the likelihood of extreme high-water levels using tidal-river records. The authors sought to explore how different sources of uncertainty interact in observational analyses of inland tidal-river stations. The paper is well structured, and the results are clearly presented, with all figures and tables relevant. Despite this, the manuscript needs to address the following concerns before being considered for publication.
Major comments
1. Lines 48-51: The limitation of conventional extreme value analysis is somewhat solved by the introduction of MEVD (Marani and Ignaccolo, 2015; Caruso and Marani, 2022), which uses ordinary events to model extremes, which is not discussed in the paper
2. Lines 98-99: please explain the cleaning of raw data, especially what is meant by "suspicious observations"? How did you decide whether to keep or remove the observation?
3. Section 2.4: Why are different parameter estimation methods used for different models? This makes comparative analysis less reliable.
4. Lines 221-222: Justify the use of different numbers of bootstrap samples for different models for uncertainty analysis?
5. Lines 261-262: Try to explain the physical reason for the evident episodic departures.
6. Lines 295-303: It is not clear how the statistical metrics are calculated. What data is used to compute the metrics? Does it use bootstrap results?
7. Also, GEV-LM achieved the lowest average RMSE and MAE. This result could be influenced by the parameter estimation methods. Have you tried all models with the same parameter estimation technique?
8. It is apparent from Fig. 4 that Gumbel-I (orange line) gives the highest estimates. Not sure about the x-limit; maybe 100 years. In contrast, in Fig. 6, at SW243, POT-GPD provides the highest return-level estimates. This discrepancy should be addressed, and the description should be modified accordingly (lines 312-320).
9. Lines 318-320 state the influence of an individual extreme on tail behavior under data-scarce settings. However, Gumbel-I uses AMS, and POT-GPD uses peaks over a threshold, which picks more peaks. Whether peaks selected for POT-GPD involve the specific peak on 18 April 1996 or are filtered out? The predictive performance of the two models should also be evaluated based on these factors.
10. Of course, the threshold for POT-GPD varies across stations, as evidenced for SW244 station. What could the physical drivers possibly be contributing to this large deviation?
11. From Fig. 9, it is clear that POT-GPD estimates are less uncertain than the other two approaches. The unusually wide uncertainty band at SW244 should be better explained.
Minor comments
1. Line 21: Do you mean sea-level rise can create an increase in the frequency of storm surges, or does it increase the magnitude of storm surges? The sentence lacks clarity.
2. Line 209: What are M0 and M1? Models?
3. Line 228: What is the reviewer summary workbook? Do you mean supplementary materials?
4. Change the word ‘Missingness’; use technical English.
5. The Q-Q agreement is evaluated using R2, which is not the best metric. It would be better to use RMSE.
6. The conclusion needs to be modified. It should properly convey the main takeaway from this study. Also, add future directions.