Use of historical information in extreme-surge frequency estimation: the case of marine flooding on the La Rochelle site in France

Nuclear power plants located in the French Atlantic coast are designed to be protected against extreme environmental conditions. The French authorities remain cautious by adopting a strict policy of nuclear-plants flood prevention. Although coastal nuclear facilities in France are designed to very low probabilities of failure (e.g., 1000-year surge), exceptional surges (outliers induced by exceptional climatic events) have shown that the extreme sea levels estimated with the current statistical approaches could be underestimated. The estimation of extreme surges then requires the use of a statistical analysis approach having a more solid theoretical motivation. This paper deals with extreme-surge frequency estimation using historical information (HI) about events occurred before the systematic record period. It also contributes to addressing the problem of the presence of outliers in data sets. The frequency models presented in the present paper have been quite successful in the field of hydrometeorology and river flooding but they have not been applied to sea level data sets to prevent marine flooding. In this work, we suggest two methods of incorporating the HI: the peaks-over-threshold method with HI (POTH) and the block maxima method with HI (BMH). Two kinds of historical data can be used in the POTH method: classical historical maxima (HMax) data, and over-a-threshold supplementary (OTS) data. In both cases, the data are structured in historical periods and can be used only as complement to the main systematic data. On the other hand, in the BMH method, the basic hypothesis in statistical modeling of HI is that at least one threshold of perception exists for the whole period (historical and systematic) and that during a giving historical period preceding the period of tide gauging, only information about surges above this threshold have been recorded or archived. The two frequency models were applied to a case study from France, at the La Rochelle site where the storm Xynthia induced an outlier, to illustrate their potentials, to compare their performances and especially to analyze the impact of the use of HI on the extreme-surge frequency estimation.


Introduction
France derives over 75 % of its electricity from nuclear energy.Most nuclear power plants in France are not located on the coasts.Only five are located on the Atlantic coast (including the Channel).Like any other installations, nuclear power plants can be subject to external influences and aggressions such as extreme environmental events (river and/or marine flooding, earthquakes, etc.).The Blayais nuclear power plant was partially flooded when storm Martin struck the French coast in 1999.A combination of an exceptional surge (outlier), of a high tide and high waves (induced by strong winds) led to the overflow of the dikes.According to Mattéi et al. (2001), these dikes were not designed for such a concomitance of events.A guide to protection, including some fundamental changes in the assessment of flood risks, has therefore been produced by the Nuclear Safety Authority (ASN, 2013).However, some issues like the frequency estimation of extreme surges remain among the priorities of the Institute for Radiological Protection and Nuclear Safety (IRSN).During the last 3 decades, France has experienced several other violent storms (the Great Storm of 1987, Martin in 1999, Klaus in 2009and Xynthia in 2010) that caused very high surges which can be considered to be outliers.An outlier can be defined as an observation that is well outside of the expected range of values in a sample.

Y. Hamdi et al.: Use of historical information in extreme-surge frequency estimation
Surge frequency estimation is an important step in the analysis of the risk associated with marine flooding.The estimation of the frequency of occurrence of extreme environmental events using probability functions has been a common issue for many decades (e.g., Dalrymple, 1960;Gringorten, 1963;Cunnane, 1987;Chow et al., 1988;Rao and Hamed, 2000;Bernardara et al., 2014).The engineer generally needs to determine the surge of a given return period T , i.e., the surge quantile X T or design surge.Traditional methods for analyzing and estimating the frequency of extreme events have been generally based on available local observations from the systematic record alone.However, it seems that these methods are not really suitable for an extreme events data set containing a surge (outlier) which is much higher than any other observed surges (e.g., Hu, 1987;Ebrary, 1999).Indeed, the presence of an outlier in a sample may lead to poor estimation of the distribution parameters.In order to base our statistical inference on the right tail of the selected distribution, detection and treatment of outliers are key elements to an effective frequency estimation and risk analysis (Barnett and Lewis, 1998;Chebana et al., 2012).During the last 4 decades, several authors (Leese, 1973; US Water Resources Council Hydrology Committee (USWRC), 1982; Stedinger and Cohn, 1986;Condie, 1986;Jarrett, 1990;Salas et al., 1994;Ouarda et al., 1998;Hamdi, 2011;Payrastre et al., 2011Payrastre et al., , 2013) ) have recognized the value of using other sources of information in the frequency estimation of extreme events.Indeed, to obtain both a consistent and accurate estimate of extreme events requires the use of a consistent technique which improves the incorporation of the outlier in the frequency estimation.Regional estimation, in which on-site observed exceptional events may become normal regional extreme observations and do not appear to be outliers any more, was considered by the scientific community to be a serious track to analyze the frequency of occurrence of the surges (Bardet et al., 2011;Bernardara et al., 2011;Weiss and Bernardara, 2013).However, the inter-sites dependency issue in the regional framework must be revised (Bardet and Duluc, 2012).Additional information refers also to historical events which have been experienced before the systematic period.HI may arise from verbal communication from the general public, written records in archives (books, newspapers, damage reports, unpublished written records, etc.) and from high-water marks left by extreme floods for instance.Other sources of HI such as paleoflood data (which can be obtained from the manipulation of certain types of proxy data) can be useful as well.A review of the literature on HI and the role it can play in a frequency analysis has been made by several authors (e.g., Stedinger and Baker, 1987;Salas et al., 1994;Ouarda et al., 1998).The basic reason for working on such a topic arises from the fact that, despite their significant impacts on nuclear-related facilities and on economic and social activities, statistical characterization of extreme storm surges, using HI, has been explicitly handled in only two references in all the literature (Van Gelder, 1996;Bulteau et al., 2014).Van Gelder (1996) used a Bayesian framework to account for known historical sea floods in the estimation of sea dikes design level in the Netherlands.More recently, Bulteau et al. (2014) have also used a Bayesian framework to estimate extreme sea levels for the design of coastal defences.The authors applied the Bayesian method to the site of La Rochelle in France.The probabilistic and statistical treatment of surge data containing outliers is also limited in the literature especially in a local frequency analysis context.
A basic hypothesis in statistical modeling of HI is that at least one threshold of perception (S 0 ) (below which the magnitude of the surge is unknown) exists and that, during a giving historical period preceding the period of tide gauging, only surges above the threshold (large enough to be remembered and/or to leave a mark somewhere on the coastal region) have been recorded or archived.Typically S 0 is the surge at which significant economic damage occurs.We named this method the block maxima method with HI (BMH).Plotting position rules, to calculate observed probabilities based on both historical and systematic information, have been proposed in the literature (Hirsch, 1987;Hirsch and Stedinger, 1987;Guo, 1990).The development proposed by Hirsch (1987) is considered herein.The reader is referred to Hirsch (1987) and to Sect.3.1 of the present paper for more details on the BMH model.Another method of using HI in the frequency estimation of extreme events, which has not been frequently explored in the literature, is based on the use of two kinds of historical data: classical historical maxima (HMax) data, and over-a-threshold-supplementary (OTS) data.In both cases, the data are structured in historical periods and can be used only as complement to the main systematic data.This method is called the peaks-over-threshold method with HI (POTH) in this paper.More details about the OTS and HMax data are provided in Sect.3.2 of the present paper.
Whatever the approach used, parametric methods based on the method of the maximum likelihood (ML) for estimating the distribution parameters have been developed and used (Leese, 1973;Hosking and Wallis, 1986;Guo and Cunnane, 1991).The selection of the ML method was based on its statistical features for large samples and for its ability to easily incorporate any additional data in the estimation process.HI are often imprecise, and we should consider their inaccuracy in the analysis.However, even with important uncertainty, the use of HI is a viable mean to decrease the influence of outliers by increasing their representativeness in the sample (Hosking and Wallis, 1986a;Wang, 1990;Salas et al., 1994;Payrastre et al., 2011).
The objective of the present work is to examine the potential gain in estimation accuracy with the use of Historical events for the two methods described above.The La Rochelle tide gauge located in the French Atlantic coast is used in this paper as a case study.The paper is organized as follows.A brief review of the basic definition and concepts is pre-sented in Sect. 2. The third section presents the case study with systematic and historical data.In the Sect.4, the developed frequency analysis models will be applied to the La Rochelle site.The results of the analysis are further discussed in the same section, before the conclusion and perspectives in Sect. 5.

Basic definitions and concepts
A frequency estimation procedure with HI must include the following steps: (i) hypothesis testing (verification of stationarity, randomness and homogeneity) and detection of outliers; (ii) selection of S 0 (for the BMH method) and the historic threshold (for the POTH-OTS); (iii) compute empirical probabilities corresponding to systematic observations (for the BMH method, only observations below S 0 are taken into account); (iv) the computation of empirical probabilities of historical data (using the Weibull formula for the POTH approach and the exceeding formula for the BMH one).Because some surge events exceeding the perception-threshold S 0 occur among the systematic data, these events are virtually removed from the systematic period and are treated as historical data (Bayliss and Reed, 2001) in the BMH frequency model; (v) the selection of a parameters-estimation method (ML in our case) and estimating a theoretical distribution function (the GEV (generalized extreme value) for the BMH model and the GP for the POTH one) to fit the observations; (vi) calculating from the theoretical distribution the T -year values of interest (say 100 or 1000 years).We analyze systematic observed skew surges that are available in the period from 1941 to 2010.Data sets were extracted from the predicted tide levels and corrected observations (Bardet and Duluc, 2012).The selected extremes for a particular year were obtained from statistically independent surges.It should be noted that the same raw data set analyzed by Bardet and Duluc (2012) is used in this paper.
A common problem in statistical modeling is the existence of gaps within the observation period (due to damage or failure of the measurement system, human errors, strikes, wars, etc.).Whatever the HI approach (BMH, POTH) to be used, we must take into account these missing periods in the analysis.Worth noting is that it happens quite often that failure in the measuring stations occurs because of the storm, thus creating a non-independent gap.Then we should ensure that the gaps occur independently from a measured variable.By means of classical extreme value distributions we are often able to analyze the statistical behavior of the maximum of a sequence of independent identically distributed (iid) random variables.Such a sequence usually represents values measured on a regular time scale.In our application the extremesurge data used for the BMH frequency model are the annual maximum (AM) surges.As it is the maximum of a block of values, it is often denoted as block maximum (Fig. 1).
Figure 1.A sketch of a systematic record.To the left: the block maxima data set (for the BMH method).To the right: the POT data set (for the POTH method).In POT only surges X i > u s are modeled through exceedances Y i = X i − u s .The distribution F (X) is known only for the upper part x > u s .
On the other hand and in the context of a point process model, surge events occur at successive random times T i when a random variable is observed.We consider that only large surge values are of interest (values exceeding a sufficiently high threshold u s ).We can assume that the times T i corresponding to large enough surges X i should be described by a homogenous Poisson process for which the number N of surges on a time interval of length w s follows a Poisson distribution with mean λw s .The λ parameter is the rate, and it is generally given in number of events per year.Surges X i are assumed to be iid random variables with continuous distribution F X (x) and density f X (x).The return period in such a case is given by T = 1/λ (1 − F X (x)).The annual maximum surges (for the BMH model) and the peaks-over-threshold (POT) surges (for the POTH model) are used in this paper as systematic data sets.The construction of the probability plot is a key step of a frequency analysis.There are several formulas that can be used to calculate the observed probability of an event in the systematic period.Based on several studies (e.g., Alam and Matin, 2005;Makkonen, 2006) the Weibull plotting position rule was used (p emp = i/(n + 1)).

Historical information
HI refers to events which have been experienced before the systematic gauging period.Although the frequency analysis using HI often lacks precision in the data used, it has been shown to improve the frequency estimation of extreme events.As mentioned earlier in this paper, there are two methods (BMH and POTH) to integrate HI with systematic data in a frequency model.Whatever the method used and comparing with the traditional frequency models, there are two fundamental differences when using HI.
i.A different class of plotting position formulas can be used for the empirical probabilities.For instance, the exceedance formula is used when applying the BMH frequency model.The exceedance formulas are based on the fact that the frequency of the whole surges ex-  (Weibull, 1939) is used, but the probabilities of HI and of the systematic data are separately calculated.
ii.A second fundamental difference when using HI in a frequency model is that an adapted MLE methodology should be used.To take into consideration available information on historical events, this methodology leads to use additional terms (depending on the method and the information type) in the likelihood expression.
The development of the plotting position formula for the BMH model and derivations of the likelihood equations for the frequency models with historical information are addressed hereafter.
Typically, historical surges are stated in terms of a large event (a huge submersion for instance) which took place in a particular year without knowing the characteristics (rank or magnitude) of the events that took place in other years.Also worth noting is that the completeness is not guaranteed.In fact, we only know extreme surges that have caused damages and thus participating in high levels.Extreme historical surges associated with low or medium tides are often not known.In this paper, surges are described in four types: (i) exact historical surges (magnitudes can be accurately estimated); (ii) events above a threshold (lower bound); (iii) an upper bound over a given period: this bound was never exceeded; (iv) surges whose magnitudes are known with some uncertainty and bounded by upper and lower limits (range).On the basis of these HI types, a joint likelihood function of the historical and systematic data (G|θ ) will be expressed and discussed in more detail in the rest of the paper.With the tools exploited in this study, both BMH and POTH models can use exact values as HI.However, the second and the last HI types can only be used with the BMH frequency model whereas the third type can only be used with the POTH model (particular case of OTS historical information without events above the threshold).

The block maxima method with HI (BMH)
In applying the BMH method one may assume that HI are available in relation to one or more fixed thresholds of perception.Typically these thresholds of perception (S t ) are surges at which important economic damage occurs.A first difficulty with this method is to quantify the magnitude of historical surge events (as it is often subject to important uncertainty) and to estimate the empirical frequencies from a sample of both systematic data and HI.A second difficulty lies in the hypothesis of the HI completeness.Indeed, historical surge data are supposed to be exhaustive and all the events exceeding S t are reported.In other words, for the entire historical period, we assume that the thresholds have not been exceeded except for the available historical data.The selected statistical distribution for the BMH model is the generalized extreme value (GEV).Indeed, the annual maxima data converge asymptotically to this distribution function.
Mathematical derivations to calculate plotting positions of systematic and historical data for the BMH model using more than a threshold of perception is detailed in Appendix A and illustrated by Fig. 2. The ML estimators are obtained by maximizing the likelihood function over the parameter space θ (or usually simpler, by maximizing the logarithm of this function).Under the assumption that the surges are iid, the global likelihood function of the whole data sample is any function L(G| θ ↔ ) which is the product of the likelihood functions of the particular types of events and information.It can then be expressed as Details about each term of this equation are presented in Appendix B.

The peaks-over-threshold method with HI (POTH)
The historical data (HMax and OTS data) are structured in historical periods and can be used only as complement to the main POT systematic data (Fig. 3).The OTS data are overa-threshold historical data recorded on periods with known durations and known exceedances.The periods are assumed to be potentially disjoint from the systematic period and other historical periods.For each period with known duration w OTS , we must have a threshold u h and all observations exceeding this threshold.The historical threshold u h cannot be smaller than the main threshold u s (used for the POT systematic data).A HMax data period corresponds to a time interval of known duration w HMax during which historical n k -largest values are available.Periods are also assumed to be potentially disjoint from the systematic period and other historical periods.The selected statistical distribution for the POTH model is the generalized pareto (GP) distribution.Indeed, the exceedances law converges to this distribution function.
Calculating empirical probabilities of historical data and constructing an empirical distribution function (in a POTH context) is the same as that of systematic data.A classic empirical formula can be used.Otherwise, as mentioned earlier in this paper, it was shown that the Weibull plotting position formula p = i/( n + 1) is more adequate than the other commonly used rules.In this formula, n is a prediction of the number of events on the historical period for HMax (since this number is unknown).A natural choice is n = λw HMax (λ is the events rate on the systematic period).For an OTS period with duration w h and with no surge observations, the  never exceeded threshold u h is shown as a horizontal segment with return period up to w h .To estimate the distribution parameters by using the maximum likelihood technique in the POTH model, let us assume a set of POT systematic observations X sys, i with two sets of historical OTS and HMax surges X OTS,i and X HMax,i .The global log-likelihood can be expressed as Details about each term of this expression are presented in Appendix C.

Data and models settings
The frequency analysis is performed at the La Rochelle site, which is located on the French Atlantic coast.Independent surges time series were obtained from the predicted tide levels and corrected observations (Bardet et al., 2012), provided by the French Oceanographic Service SHOM (Service Hydrographique et Océanographique de la Marine) for the tide gauge named "La Pallice".One of the most important features of the station La Rochelle is the fact that the region, in which this station is located, has experienced important storms during the last 2 decades (the storm Martin in 1999 and Xynthia in 2010).Figure 4 displays the geographic location of the La Rochelle site.

Systematic data and statistical tests
In the case of the BMH model, the annual maximum surges available in the period from 1941 to 2010 constitute the systematic record and all events occurred before 1941 are considered to be HI.The choice of systematic data for the POTH approach is different.Indeed, to use all available data separated by periods of missing data, we concatenated all the surges gauged since 1941 to form one systematic record.It was shown in previous researches (e.g., Bardet et al., 2011;Hamdi et al., 2014) that a threshold u s equal to 41 cm is an adequate choice that gives a rate λ equal to 4 events year −1 .BMH and POTH systematic data from 1941 to 2010 are represented by bar plots in Figs. 5 and 6 respectively.
Stationarity, homogeneity and randomness of time series are required conditions in a frequency analysis (Rao and Hamed, 2001).The Kendall test for stationarity (Mann, 1945), the Wald-Wolfowitz test for randomness (Wald and Wolfowitz, 1943) and the Wilcoxon test for homogeneity (Wilcoxon, 1945) are the three non-parametric tests used as a prerequisite for frequency analysis.The La Rochelle station passed these tests at the significance level of 5 %.The Grubbs-Beck test (GBT) for the detection of outliers is also  used in the present work.The results of this test show clearly that the 2010 event is undoubtedly an outlier.The reader can be referred to Grubbs and Beck (1972) for more details on this test.

Settings of the BMH and POTH frequency models
Some historical known events experienced by the La Rochelle site are available.We chose in this part to consider only the information collected with a tide gauge before 1941.This information was collected as part of the PhD thesis of Gouriou (2012) at the University of La Rochelle and corresponds to a 12-year surges record (between 1863 and 1874) at the Old Harbour of La Rochelle.Only the four highest values are considered in the present work (Table 1).Through Garnier and Surville (2010) or Breilh et al. (2014), we could identify that some other extreme events that happened in the La Rochelle area, during gaps in the systematic period exploited in the present work.In particular, an event occurred in 1957 is comparable in terms of consequences to the Xynthia event.But these events have not been measured by tide gauges or the estimated levels are not clearly associated with La Rochelle site.However, this is worrying for us since ex-  haustiveness in HI has been assumed in the developed BMH model and in the POTH model with OTS data.A possible solution is to choose a threshold of perception S t (in the case of the BMH method) or a historical threshold u h (in the case of the POTH method) as high as possible.This would allow us to get closer to the completeness of the HI.Care must be taken when using the available HI (Table 1), it is not necessary to use all the available historical data and depending on the selected threshold, one or more historical data can be excluded from the analysis.This is indeed the case here, only the historical surge of 111 cm induced by the event that took place in 1866 was used.Some elements of justification will be presented in the next section.The other data are considered quite small and taking a threshold of perception below these values further exacerbate the problem of the assumed completeness.
We introduced in a previous section some parameters characterizing data sets (containing both HI and a systematic record) that we apply here on the La Rochelle records.

Settings of the BMH frequency model
As mentioned in Sect.4.1, events occurred before 1941 are considered to be HI.The threshold of perception is not easy to estimate.In addition to the fact that it was not exceeded during all the historical period except when information is available, the value of the threshold of perception must meets the following criteria: (i) should be as high as possible; (ii) should be exceeded by the selected historical values introduced in the estimation; (iii) may result in serious economic damage and disruption.In the case of marine surges, the last criterion is not always taken into account for the simple reason that high surges can occur at the same time with low tide and may not cause damage.
The 2010 event (Xynthia) is exceptional (comparing to the other observed systematic surge peaks) and affected not only the French coasts but also western Europe.The surge value (153 cm) caused by the Xynthia event has largely exceeded twice the average of the other surge records (61 cm).This event is considered to be historic in the BMH frequency model.This changes its empirical probability but not the value of the overall likelihood.Consequently, the theoretical probabilities remain unchanged.All other systematic surges seem to be ordinary extreme ones with relatively low return periods.The choice of the threshold of perception becomes then a goal easy to achieve; we simply need to choose a threshold that allows us to consider only the Xynthia event to be historic among the systematic ones.It should range between the Xynthia event surge (153 cm) and the maximum value of all the other systematic observations (90 cm).A value of S t = 110 cm with a total of two exact values (n k = 2) is then utilized in the present paper.Only one of these values is systematic (n e = 1) and there are no lower bounds (lb = 0) and ranges (n r = 0) that can be used for this case study.Figure 5 illustrates well these model settings.

Settings of the POTH frequency model
The parameters characterizing a POTH model data set are u h and w h .In order to be consistent with the BMH model assumptions and to compare the two models results, the HI is used in this work as OTS data.As shown in Fig. 6, HI complements then the systematic record with these OTS data recorded on one historical period (1866-1940) with known duration (w OTS = 75 years) and one exceedance (n k = 1).As was the case for the BMH frequency model and to ensure exhaustiveness of HI, a threshold as high as possible should be chosen.A value of u h = 110 cm is then utilized.Unlike the BMH frequency model, systematic events that exceeded the historical threshold are not considered to be historical in the POTH model.Consequently their empirical probabilities will be different (compared to their equivalents in the BMH model) but theoretically the likelihood remains insensitive to the fact that these observations are considered to be systematic or historic.Worth noting is that theoretically, with the proposed BMH and POTH models settings, the portion of the likelihood related to the addition of the HI must converge to the same value and the impact of the use of the HI must be the same for the two concepts.
The POTH frequency model is applied with the Renext tool developed by the French Institute for Radiological Protection and Nuclear Safety (IRSN).The Renext software proposes fits and diagnostics for the POT method.The software is programmed in the R environment (open-source software for statistical computing: http://www.r-project.org/).

Results and discussion
In this section we report the results of the BMH and POTH frequency models for the extreme storm surges frequency es-timation and analysis applied to the La Rochelle station.All the developments and simulations were carried out within the R environment.The packages "evd" (Stephenson, 2012) and "Renext" (Deville, 2013), which have specifically been developed for extreme events frequency analysis, were used.As shown in Figs.5-6, the surge frequency analysis was performed considering only one scenario for each method: 30 systematic surges and 1 exact historical observation for the BMH frequency model and 106 systematic POT surges and 1 OTS data for the POTH model.Whatever the frequency model settings, the historical period considered in these two cases starts when the first historical data was observed (75 years long).For each model, the evaluation of the impact of additional HI (by comparing quantiles with no historical data and those obtained using the whole 1866-2010 data set, including both systematic and HI), constitutes the main result of this paper, and it is highlighted later in the present section.These results are presented in form of probability and diagnostic plots (Figs.7-8) and in form of tables summarizing the numerical estimates of the return levels of interest and uncertainty parameters.
Fitting the GEV/GP distributions to the data yields the following parameters, the absolute and relative widths ( CI and CI/S T , respectively) of the confidence intervals of the representative surge quantiles which are presented in Tables 2  and 3. Since high return periods are needed for the safety of nuclear installations, the main focus was set to the 100, 500 and 1000-year surge quantiles.The uncertainty corresponding to the return level is closely related to that of the model parameters.The Delta method was used to calculate the variance of the return level estimates.The distribution functions (GEV for the BMH model and GP for the POTH model) was fitted to data. Figure 7-8 consists of 4 subplots.In each figure are displayed the fitting and the Q-Q plots (top left panel and top right panel, respectively) excluding the HI.The graphics with HI are presented in the bottom.The solid line in the middle is the fitted distribution function, while the dashed lines are the lower and the upper bounds of the 70 % confidence intervals of surge quantiles.In the plots at the bottom (with HI), the systematic data are represented by empty circles, while the HI is depicted by a solid red circle.In the next sub-sections, we present the results of each frequency model separately.However, we can afford to analyze first in the present sub-section the results when no HI is included.It can be seen that the fit with no historical data included, presented in Figs. 6 and 7, is poor at the right tail of the distribution and more particularly at the outlier.It is observed that the POT model gives lower surge estimates for the La Rochelle station than the AM model.Quantiles corresponding to 100, 500 and 1000-year return periods given by the POT model are approximately 6 cm lower than those obtained by the AM method.The corresponding confidence intervals are very large when the annual maxima data are used, and it is much narrower in the case the POT frequency model is used.This was expected because the robustness of the re-  sults greatly depends on the amount of available information and since the POT method generally uses more data than the AM one.Also note that this difference between the AM and the POT models should decrease, and the estimated values become fairly close when additional data (given by historical information in our case) is used.This also has been proven in a theoretical way by the scientific community (e.g., Langbein, 1947;Chow et al., 1988;Hamdi et al., 2013).

Results of the BMH frequency model
In a classical frequency analysis (with systematic data only), some goodness-of-fit criteria and tests (e.g., χ 2 , Kolmogorov-Smirnov and Anderson-Darling tests) can be applied.As soon as any HI is used, testing the goodness of the fit becomes a more complicated business.The way additional HI affects the quantile and uncertainty estimates can even though be analyzed using a visual inspection of diagnostic plots (Fig. 7) and the numerical results presented in Table 2.As shown in Fig. 7, the use of HI appreciably changes not only the empirical probabilities and the estimated quantiles but also the uncertainty associated to these estimates.The observed probabilities of surges that exceeded the threshold of perception were calculated using the exceedance formula.It can be seen that this formula allows better positioning of the outliers in the systematic data and gives to these events more reasonable return periods.It can then be concluded that the fit with empirical frequencies of large surges (at the right tail of the distribution and especially at the outlier) is more satisfactory when using HI.
Including HI reduces the relative width of the confidence bounds.This finding is underpinned by the values of the absolute and relative widths of the confidence intervals of surge quantiles presented in Table 2.This is the case when more data are included in the parameters estimation (if the extra data are consistent with previous ones).The relative widths of the confidence bounds for the BMH frequency model with no HI involved are 2.4 times larger than those obtained with both systematic and historical surges.This holds for all the return periods of interest (100, 500 and 1000 years).This fact is also graphically underpinned by the confidence intervals plots presented in Fig. 7.The BMH model with both systematic and historical data gives return levels systematically lower than those given with systematic data only.

Results of the POTH frequency model
The visual inspection of diagnostic plots is complicated when the POTH frequency model is used.Remember here that the empirical probabilities are calculated with the Weibull plotting position rule (not with the exceedance formula), but in a separate way for the historical information and the systematic data.Their positions are therefore not the same as those given by the BMH model.We can quickly see from Fig. 8 that the visual fit at the outlier is almost the same as that obtained with no HI included.The fit of the large empirical frequencies, at the right tail of the distribution excluding the outlier, is slightly better than that obtained without HI.The improvement is not significant, and as shown in Fig. 7, several observations have remained outside the confidence interval.
The results given in Table 3 indicate that the values of relative confidence intervals widths with no HI involved are 1.65 times narrower than those obtained with both systematic and historical surges.In other words, the user of this method is more confident in its parameter and quantile estimates when using HI.As with the BMH model, return levels given by the POTH frequency model with HI are also lower than those obtained with no historical data included.And what are the shortcomings and the advantages of the two frequency models?
In order to compare the two approaches, the two model settings were chosen to be as close as possible.This is why for the POTH model we have considered one period ending where the systematic period began (the same period as that used in the BMH model).The POTH historical data can be handled differently; one could use for instance the four available historical data (Table 1) over the corresponding period of 12 years without making assumptions about the nongauging period .From Tables 2 and 3, the convergence of the AM quantiles towards the POT ones is clearly observed when the HI is used.100, 500 and 1000-year quantiles obtained by the POTH method are close to those given by the BMH one (with a difference of only 1.8, 1.6, and 1.4 % respectively).
This convergence can be seen as tangible evidence of the robustness of the two models results.Unlike the POTH model, the BMH one gave a relatively better visual fit not only at the outlier but also at the other large empirical frequencies.The confidence intervals obtained by the POTH method are much more advantageous (much tighter), and the relative width is only 20 % for the 100-year return surge (against 60 % for the BMH model).
The biggest disadvantage of the BMH frequency model comes from the completeness assumption.It is assumed that during the historical period the perception threshold was never exceeded except in years when the information is available.This hypothesis can quickly become unsuitable if we have evidence that one or more significant events occurred during the historical period and we cannot get more details on these events The application with one piece of HI (1863-1874) follows on nicely from the methodology proposed in this paper.However, that the single additional surge (111 cm) has made the Xynthia datum seem more consistent with the estimated distributions could just be a coincidence arising as a consequence of the extreme value distributions' shape parameter being difficult to estimate with precision.In addition, the fact that the empirical probabilities are better positioned with the exceedance formula can be a misleading element in our assessment of the fits and the uncertainty.This section deals with the impact of additional extreme value data, which occurred during both the systematic and historical periods, on the robustness and the uncertainties in both the frequency estimation of the Xynthia event and the extrapolated levels corresponding to other desired return periods.
As mentioned in Sect.4.2, additional extreme surges are available.Breilh et al. (2014) have constructed a historical database of storm-induced coastal floods that affected the La Rochelle region over the last 5 centuries.The analysis of the documented events (some flooding events described as "comparable to Xynthia") led the authors to propose observed and modeled maximum water levels possibly reached during these events (Table 4), in particular at La Pallice, where the systematic data were measured.Three of the deducted skew surges (1941, 1957 and 1999) are generated by events that have occurred during the gaps of the systematic period.Fortunately, the developed POTH frequency model has been developed to manage this type of data by considering them to be HI.Moreover some HI was collected from archives by Garnier and Surville (2010) or Gouriou (2012).
It is noteworthy that, unlike the 1863-1875 data, this additional information is not clearly related to the La Rochelle station and not always precise (some events as in 1940 are known with intervals); therefore, care must be taken when using these data in the developed frequency models.Therefore, to complete the data of Breilh et al. (2014), only the high skew surges collected for the 1924 and 1940 events are considered in this work, as presented in Table 4.
The question that rises is the following: would additional extreme HI improve the fits and the uncertainties?Four study cases were considered with the POTH model: (1) frequency estimation with systematic data only (w s = 26 years); (2) the 1999, 1957, 1941 and 1940 events are considered to be HMax data with a length of historical period w HMax equal to 45 years (the 1940-2010 period without the systematic data); (3) the 1924 historical surge is added to the HI considered in case 2 (the length of the historical period becomes 61 years by adding the 1924-1939 missing period) and (4) the 1863-1874 piece of HI (the four highest surges) is considered to be a separate historical period with an effective duration w HMax1 equal to 12 years, and is combined with the five previous historical surges over a historical duration of 110 years.For the 1940 surge which is known with an interval, we considered a mean value with the POTH model.The estimates of the return levels (in cm) of interest and corresponding relative widths (in %) of the 70 % confidence intervals for each case are summarized in Table 5.
This complementary discussion illustrates the important reduction of confidence intervals associated with the use of additional historical data.This reduction clearly appears in the four case studies.Indeed, the results given in columns 3, 4 and 5 of Table 5 indicate that the value of the relative confidence interval width decreases as we add more historical data.However the value of the relative confidence interval width with no HI involved are 1.4-1.5 times narrower than that obtained with both systematic and all of the historical surges.In other words, the user of this method is more confident in its parameters and quantile estimates when using a more exhaustive HI.For instance, the use of only one event (1866) as in Sect.5.2 could be misleading because the small value of the relative confidence interval can be explained by the fact that we consider only two surges exceeded 110 cm over one century, whereas we observed five other extreme surges at least over the same period.With the same consideration, it can be also shown from Table 5 that the return levels given by the POTH frequency model with HI are significantly higher than those obtained with no historical data included or with only the 1866 event for which too much importance was given to the exhaustiveness of a medium historical value.
Another issue that has arisen in the present paper is about the historical data that were not used because they were below the threshold of perception (Table 1).Assuming that these data are of valuable information, we strongly believe that these data must be exploited.Theoretically, nothing prevents us from considering them as systematic information even if it occurred 50 years before the beginning of the systematic period.In other words, systematic data that have exceeded the threshold of perception are considered to be historical information and historical data that have not exceeded this threshold would be considered as systematic information, subject to an adapted length of the whole systematic period.

Conclusion
In this paper, two methods for how to use historical surges into the local frequency analysis have been presented and applied to the La Rochelle site located on the French Atlantic coast.The first method is based on the presence of a perception threshold using annual maxima data (BMH model) and the second method is based on the use of historical data periods and POT surges (POTH model).Two adapted likelihood functions and different class of plotting positions formulas are then built to properly handle the information on historical surges.Several types of HI can be considered in the frequency analysis when using these methods.Unlike the BMH frequency model, systematic events that exceeded the historical threshold are not considered as historical data in the POTH model.Consequently their empirical probabilities are different and fitting to the observations was improved especially with the BMH frequency model, but the likelihood remains the same with similar hypotheses for the historical information.Although the two approaches are physically different, the use of historical information has close impacts when using the two methods with the available data and the applied settings.The comparison of the theoretical distributions fitted to AM and POT extreme events have indicated that, at high return periods (greater than 100 years), the two approaches give comparable predictions of surge magnitudes when historical information is used.This can be interpreted differently; the amount of information is more critical and important than the choice of the approach.
Care is always needed when modeling data with the aim of extrapolating whatever the model used (BMH or POTH) and regardless of whether or not HI is used.As indicated by the results of the case study and reasoning in terms of high return levels and associated uncertainty, we have two models that handle differently the impact of the use of HI on high quantiles estimation.In such an exercise in which extreme values are used with the aim of extrapolating, it is very important to account for the underlying error with which we make our inferences.Despite the improvement of the uncertainty degree resulting from the use of HI in the BMH model, this uncertainty remains high.The distribution parameters and quantiles estimates are qualified with a much more appropriate degree of uncertainty when using the POTH frequency model, because more systematic data are available with this model.On the other hand, how the theoretical distribution is representative of the observations is also an important is-sue.Overall observations are relatively well fitted except at the outlier when the POTH frequency model is used.Unlike this model, the BMH one provides a much better visual fit at large frequencies and tends to similar return levels when HI is used.Indeed, the BMH frequency model has the advantage of allowing the repositioning in history of the systematic observations that exceeded the threshold of perception (the outlier for instance).However, the BMH frequency model needs to assume the exhaustiveness of HI, whereas the POTH method is more faithful to the observations and needs fewer hypotheses with regards to the empirical distributions.
Overall, adding information on historical surges to the local frequency models has reduced the variance of the distribution function parameters estimates, except if the assumption of exhaustiveness is too strong as with an only medium surge over a large historical period.The relative confidence intervals widths decrease as more historical data are used.However additional historical data do not necessarily lead to more robust estimates, the high return levels being very sensitive to the historical values.Then it is absolutely necessary to have historical data the more complete possible, in order to have both robust estimates and reduced uncertainties.
During our research that led to the writing of this paper we noted three perspectives.The first one is related to the kinds of historical information, which are different in the two models (lower and/or upper bounds, threshold of perception, etc.).The BMH and POTH models could be adapted to consider similar HI, and the plotting positions for the POTH model could be calculated with the exceedance formula.In addition, given the importance of the exhaustiveness of the implemented historical data, another kind of HI considering "at least k events above the threshold of perception" would be very interesting for the two models, then partially answering to the exhaustiveness issue.This treatment would be particularly relevant when the number or the values of historical data cannot be accurately estimated, and could be based on a "binomial censored data" approach (Stedinger and Cohn, 1986).The second perspective is relative to the spatial extension of extreme surges data including HI.It can be envisaged that using the two methods presented herein may improve the frequency estimation and uncertainty in a regional framework, by increasing the amount of information.This work is ongoing.Finally, the last perspective is no less important and may be complementary to the first and the second ones.The question that has arisen is about the historical data that were not used because they were below the threshold of perception.Assuming that these data are of valuable information, we think that these data have to be exploited.

Y. Hamdi et al.:
Use of historical information in extreme-surge frequency estimation evaluated at the observed sample, and it is the product of the likelihood functions of the particular types of events and information E i : L(E i |θ ). (B1) The events can be the (s − e) systematic events (x i ), the n k historical surges (y i ) above the threshold of perception S t with exactly known magnitudes, the n low lower bounds historical events y low,i i=1,...n low and the n r range historical events y low , y up i .
-ML estimator for the BMH systematic data: the likelihood function is given by the joint distribution of X 1 , X 2 , . .., X s-e .It can be expressed as f X (x i , θ ) . (B2) -ML estimator for the BMH historical exact data: the likelihood function is given by the joint distribution of n K , Y 1 , . .., Y n k , It can be expressed as f X|X>S t y j , θ , (B3) where y i denotes a surge exceeding the threshold of perception S t and observed in the total period n (as described in Sect.2.3).Since f X|X>S t (y) = f X (y)/p, the above equation simplifies to f X y j , θ .-ML estimator for the OTS historical data: for HI, consider an OTS period with threshold u h , duration w OTS and with given k observed surges X OTS,i for i = 1, 2, . .., k.The events with X OTS > u h form an homogenous Poisson process with a rate λ (1 − F (x, θ )) and their number K follow a Poisson distribution.Conditional on {K = k}, the k observed surges are independent with density f (x)/ (1 − F (u h , θ )) for x > u h .Hence

Figure 2 .
Figure 2. A sketch of systematic records and historical information: a multiple-threshold case.

Figure 3 .
Figure 3. POTH: a sketch of systematic records and historical information (HMax and OTS data).

Figure 4 .
Figure 4. Geographic location of the La Rochelle site.

Figure 5 .
Figure 5. La Rochelle site -systematic and historical data and settings of the BMH frequency model.

Figure 6 .
Figure 6.La Rochelle site -systematic and historical data and settings of the POTH frequency model.

Figure 7 .
Figure 7.The GEV distribution fitted to the AM of peak surges at the La Rochelle station: with no historical information included (top left panel) and the corresponding Q-Q plot (top right panel); with both systematic and historical data (bottom left panel) and the corresponding Q-Q plot (bottom right panel).

Figure 8 .
Figure 8.The GP distribution fitted to the POT surges at the La Rochelle station: with no historical information included (top left panel) and the corresponding Q-Q plot (top right panel); with both systematic and historical data (bottom left panel) and the corresponding Q-Q plot (bottom right panel).
One of the most important questions still remains is about the robustness: what is the impact of the use of the HI on the BMH and POTH models results (for high return levels)?

--
ML estimator for the BMH historical lower bound data: ,i , θ − F y low,i , θ .ML estimator for the POT systematic data: let us assume a set of n POT systematic observations X i and a selected threshold u s and consider w s the total duration.The likelihood L X sys,i |θ isL X sys,i |θ = Pr (N = n) • n i=1 f X sys,i , θ .(C1)For a homogeneous Poisson process with rate λ, the above equation can be re-expressed asL X sys,i |θ = (λw s ) n n! exp (−λw s ) • n i=1 f X sys,i , θ .(C2)The log-likelihood X sys,i |θ is then: X sys,i |θ = n log (λw s ) − log (n!) − λw s + n i=1 log f X sys,i , θ .

Table 2 .
The T -year quantiles (in cm) estimated on the basis of the systematic and historical periods and absolute and relative widths of their 70 % confidence interval (the BMH frequency model).

Table 3 .
The T -year quantiles (in cm) estimated on the basis of the systematic and historical periods and absolute and relative widths of their 70 % confidence intervals (the POTH frequency model).

Table 4 .
Breilh and al. (2014)(2010)Historical skew surges estimated with the flooding levels collected byGarnier and Surville (2010)in the La Rochelle area.bHistoricalskew surges estimated with the maximum water levels modeled byBreilh and al. (2014)at La Pallice. a

Table 5 .
The T -year quantiles (in cm) estimated on the basis of the systematic and historical periods and relative widths (in %) of their 70 % confidence interval (the POTH frequency model).