Modelling Prospective Flood Hazard in a Changing Climate, Benevento Province, Southern Italy

: The change of the Earth’s climate and the increasing human action (e.g., increasing impervious areas) are inﬂuencing the recurrence and magnitude of ﬂooding events and consequently the exposure of urban and rural communities. Under these conditions, ﬂood hazard analysis needs to account for this change through the adoption of nonstationary approaches. Such methods, showing how ﬂood hazard evolves over time, are able to support a long-term plan of adaptation in hazard changing perspective, reducing expected annual damage in ﬂood prone areas. On this basis, in this paper a reevaluation of ﬂood hazard in the Benevento province of southern Italy, is presented, providing a reduced complexity methodological framework for near future ﬂood hazard prediction under nonstationary conditions. The proposed procedure uses multiple nonstationary probability models and a LiDAR-derived high-resolution inundation model to provide present and future ﬂood scenarios in the form of hazard maps. Such maps are derived using a spatialization routine of stage probability across the inundation model that is able to work at di ﬀ erent scales. The analysis indicates that, overall, (i) ﬂood hazard is going to decrease in the next 30 years over the Benevento province and (ii) many areas of the Calore river ﬂoodplain are going to be subject to higher return level events. Consequently, many areas would require new guidelines of use as the hazard level decreases. Limitations of the analysis are related to the choice of the probability model and the parameter estimation approach. A further limit is that, currently, this method is not able to account for the presence of mitigation measurements. However, result validation indicates a very high accuracy of the proposed procedure with a matching degree, with a recently observed 225-years ﬂood, estimated in 98%. On this basis, the proposed framework can be considered a very important approach in ﬂood hazard estimation able to predict near future evolution of ﬂood hazard as modulated by the ongoing climate change.


Introduction
The change of the Earth's climate and the increasing human impact are influencing the water cycle with effects on river discharge [1]. Throughout the world, the ongoing warming associated with local change in land use and increasing urbanization is modifying both the rainfall pattern and the response of catchments (in terms of discharge modulation) in a complex manner [2][3][4][5]. As a consequence, extreme river discharge change is highly variable from a region to another. In China, for example, rivers of the northern semi-arid area are experiencing an increase in river discharge greater than those in the southern subtropical humid area due to a more accelerated warming and wetting [6]. In the United States, increasing trends in many flood flow series have been identified [7]. In northern Europe, the ongoing climate change is generally increasing annual maximum and decreasing minimum river flow [8,9]. In southern Europe a robust decreasing trend in water availability has been detected [9]. Accordingly, in Italy, maximum river discharge has a consistent negative trend over the whole national territory. This downward trend has been explained as an effect of the increasing temperature [10].
Changes in average and extreme river discharges, as a consequence of variation in rainfall distribution, influence the recurrence and magnitude of flooding events and consequently the exposure of urban and rural communities [11]. While many areas are experiencing an increase in exposure to flooding, due to an increase in flood frequency (Southeast Asia, Peninsular India, eastern Africa and the northern half of the Andes; [12,13]), other areas are dealing with a consistent decrease of these events and related increase of potential drought hazard [14][15][16]. In both cases, flood hazard analyses need to account for this variability to reliably support land planning and policy decision oriented to protect communities from floods.
Hazard analysis in flood perspective can be completed through a variety of approaches based on different data types [17][18][19][20][21][22]. While geomorphological observations can support flood hazard mapping in ungauged basins, fluvial stage records can support stationary or nonstationary probability analyses [23][24][25][26]. For instance, Guerriero et al. [27] used stationary probability analysis of multiple fluvial-stage time series to map flood hazard across the Benevento Province. To produce maps form probability models, they used a flood inundation model derived by LiDAR data (e.g., numerical floodplain models; [28]) as a basis. This analysis was completed under the assumption of stationarity avoiding consider the potential presence of trends in available time series.
In presence of significant trend and/or shift in long-term fluvial stage records, consequence of the ongoing climate change and associated anthropogenic changes in river basins, flood hazard analysis needs to account for non-stationarity, showing how flood hazard evolve over the reference time period [26,[29][30][31][32][33]. A simple way to account for non-stationarity of hydrological records and predict near future hazard is to enhance standard extreme value analysis by time dependent modelling of function parameters [32,34]. Other methods use single or two-way (nonstationary-stationary-nonstationary) time series transformation and extreme value analysis to derive nonstationary models [25,35].
These nonstationary approaches have been recognized as one of the current needs for an improved analysis and optimization of water resources risk assessment [36]. Such methods should support also flood hazard analysis in order to prevent underestimation or overestimation of the current and future probability of a disastrous extreme flood events and related investment planning. However, they require at least long hydrological records to reliably predict future conditions in a rapidly changing world [30,37]. In Europe, major rivers have been monitored since at least 18th century. In Italy, the first monitoring stations have been installed around 1918 and their number has consistently risen, reaching hundreds of units. The resulting dataset can be of great importance in flood hazard assessment [38], having the potential to support flood hazard evaluation also in consideration of the ongoing climate change.
On this basis and with the aim of providing a reduced complexity methodological framework for near future flood hazard prediction and improving the understanding of flood hazard in the Benevento province of southern Italy, a nonstationary frequency analysis of data from multiple hydrometric stations has been completed producing multiple flood hazard scenarios. Especially, nonstationary generalized extreme value analysis, enhanced by time varying fitting parameters, and a flood inundation model, were used to derive multiple flood hazard scenarios that describe how perspective flood hazard is going to change in relation to the ongoing global warming. It is important to notice that the study area has been already investigated by [27] assuming stationarity of hydrologic time series. The analysis presented here, and associated framework, represents an evolution of that proposed by [27] modified as to account for non-stationarity of time series. In this way, the analysis provides a new original perspective about current and near future flood hazard across the central sector of the Benevento province.

Study Area
The study area extends crosswise the Benevento province in southern Italy ( Figure 1). It comprises a segment of the Calore river floodplain and three segments of the Ufita, Tammaro and Sabato rivers floodplains. The considered river segments are approximately 50 km, 1.5 km, 5 km and 11 km long, respectively. The less developed floodplains segments belong to the tributaries of the Calore river. A 6 km long segment of the Volturno river is also part of the study area. These rivers are consistently characterized by a single active channel and eight hydrometric stations actively monitor river stages. Of these stations, only five, installed before 1999, have a time history characterized by a time-span longer than 40 years that can be considered significant in nonstationary frequency analysis perspective ( Table 1). The stations selected for the analysis consisted of a staff gauge and water stage data were acquired using the observational method. Since AD 2000, visual staff gauges were discontinued and replaced with continuous recording stations equipped with ultrasonic water level sensors. It is expected that measurements completed with the new system are consistent with measurements completed using visual inspection of staff gauge. Exceptions are the Benevento and the Chianche stations that were re-located in AD 2000; in these cases, a geometric correction based on field measurements and discontinued and working hydrometers measurements was adopted.
Water 2020, 12, x FOR PEER REVIEW 3 of 17

Study Area
The study area extends crosswise the Benevento province in southern Italy ( Figure 1). It comprises a segment of the Calore river floodplain and three segments of the Ufita, Tammaro and Sabato rivers floodplains. The considered river segments are approximately 50 km, 1.5 km, 5 km and 11 km long, respectively. The less developed floodplains segments belong to the tributaries of the Calore river. A 6 km long segment of the Volturno river is also part of the study area. These rivers are consistently characterized by a single active channel and eight hydrometric stations actively monitor river stages. Of these stations, only five, installed before 1999, have a time history characterized by a time-span longer than 40 years that can be considered significant in nonstationary frequency analysis perspective ( Table 1). The stations selected for the analysis consisted of a staff gauge and water stage data were acquired using the observational method. Since AD 2000, visual staff gauges were discontinued and replaced with continuous recording stations equipped with ultrasonic water level sensors. It is expected that measurements completed with the new system are consistent with measurements completed using visual inspection of staff gauge. Exceptions are the Benevento and the Chianche stations that were re-located in AD 2000; in these cases, a geometric correction based on field measurements and discontinued and working hydrometers measurements was adopted.  This area has been repeatedly affected by river floods and local relict fluvial scarps, represent the signature of such historic floods. The most recent event, a destructive flood induced by the overflow of the Tammaro and Calore rivers, occurred in 2015. It consistently damaged the town of Benevento and rural areas located along the Calore river valley [22,[39][40][41]. During this event the river  This area has been repeatedly affected by river floods and local relict fluvial scarps, represent the signature of such historic floods. The most recent event, a destructive flood induced by the overflow of the Tammaro and Calore rivers, occurred in 2015. It consistently damaged the town of Benevento and rural areas located along the Calore river valley [22,[39][40][41]. During this event the river stage of the Calore river increased of several meters reaching approximately 9.5 m and 10 m at the hydrometric stations of Benevento and Solopaca, respectively.

Materials
The nonstationary flood hazard analysis across the considered sector of the Benevento province was completed on the basis of fluvial-stage time series characterized by a significant length of their records. In this case, considering the selected nonstationary approach to the problem, only records longer than 40 years were selected. Such records consist of daily measurements completed by visual inspection of staff gauge, before AD 2000, and using a continuous recording system equipped with ultrasonic water level sensors, after AD 2000. Available daily fluvial stage records were used to derive annual maxima time series characterized by a sample size between 42 values and 68 values ( Table 1). Of the eight stations installed in the study area, only five were considered for the analysis, three installed along the lower course of the Calore river, one along the Sabato rivers, and one along the Volturno river ( Figure 1

Nonstationary Frequency Analysis and Reference Flood Scenario
Nonstationary frequency analysis of selected fluvial stage records began testing time series for trends. In this perspective, a linear regression function was used e.g., Hess et al. [42] and Kundzewicz et al. [43]. This method was selected because of its simplicity and the presence of discontinuities in analyzed time series that prevent the use of more complex nonparametric test, like the commonly used Mann-Kendall test. In the assumption that a time series ((Yi) = Y1, Y2, . . . , Yn) is representative of a stochastic process in which a deterministic linear trend component is present (i.e., non-stationarity), the linear modelling function has the form: where t is the time, a is the intercept (i.e., location parameter), b is the trend coefficient, and ε t represents the estimation error. Function parameterization was automatically completed using the least squares method in MATLAB™ environment and significant trends have been detected for b 0 and a p-Value lower than 0.1 that confirm rejection of the null hypothesis of no trend in each time series e.g., Hess et al. [42].
Once detected a consistent trend in all of the available records (i.e., a downward trend), nonstationary frequency analysis of annual maxima time series was completed by adopting the approach of using standard extreme value analysis enhanced by time dependent modelling of function parameters e.g., Coles [34]. Especially, standard extreme value analysis was completed using a Generalized Extreme Value parametric distribution function GEV [44] to describe available annual maxima time series (Table 1; station positions are in Figure 1). This function has the following form: where ξ is the shape parameter, σ is the scale parameter and µ is the location parameter. To account for non-stationarity, the location parameter µ is considered to change with time and its variation is represented by a linear trend model e.g., Coles [34], Delgado et al. [45] and Gül et al. [46] in which the trend µ trend is considered equal to the trend of the time series: Parameters estimation was completed using the gev.fit function of the ismev package implemented in R (http://cran.r-project.org/package=ismev; Table 2) and a matrix of covariates corresponding to the modeled time (including discontinuities) was used for nonstationary condition simulation. Models parameterization was based on the maximum likelihood estimator (MLE) that provide function parameters and associated standard error (Table 2). To describe and quantify uncertainty in the estimation of flood magnitude, related to the characteristics of time series (i.e., sample size and continuity [47,48]) the 95% confidence intervals were derived on the basis of the estimated standard error. Since the location parameters varies with time, such estimation was completed considering the initial location parameter µ 0 as reference. To estimate the degree of fit of probability models, residual analysis was completed considering a probability models in which location parameter is estimated such that R2~max (i.e., µ = µ 0 + µ t × sample size/2). Such analysis indicates a very high goodness of fit of the probability models (R2 ≥ 0.94; Figure 2).
Following the method described by Guerriero et al. [27], developed probability distributions were used as basis to prepare multiple probability of exceedance databases representative of each selected scenario. Each database, containing probability models for each considered station, was developed using a specific location parameter calculated in reference to the specific flood hazard scenario to be modeled (corresponding to a defined time period in year; i.e., the t parameter of Equation (3)). Particularly, a past flood hazard scenario corresponding to the flood hazard in AD 2000, a current flood hazard scenario corresponding to the flood hazard in AD 2020 and two near future scenarios corresponding to the flood hazard in AD 2030 and AD 2050 were considered.
Since flood hazard analysis requires the definition of a boundary condition, in other words the identification of a design or extreme event [49], a 500-years flood has been considered. This choice is consistent with the guideline of the Italian Ministry of Environment suggesting that a 300-years to 500-years event should be considered as reference flood scenario in hazard evaluation perspective. Nonstationary models indicate that the 500-years event has an estimated magnitude that varies in time and space as reported in Table 3.  2. Time series (a,d,g,j,m) recorded at hydrometric stations considered for the analysis, nonstationary Generalized Extreme Values (GEV) probability models (b,e,h,k,n) indicating initial conditions (i.e., µ = µ0) and best fit condition (such that R2~max; i.e., µ = µ0 + µt × sample size/2) and variation in time (years) of the location parameter of the nonstationary GEV models (c,f,i,l,o). Dark dashed lines of GEV models graph depict 95% confidence interval and provides an overview of the uncertainty associated to estimated flood magnitude at specific return period in the initial condition (i.e., µ = µ0). As expected, low frequency flood has much significant uncertainty associated to magnitude estimation.

Figure 2.
Time series (a,d,g,j,m) recorded at hydrometric stations considered for the analysis, nonstationary Generalized Extreme Values (GEV) probability models (b,e,h,k,n) indicating initial conditions (i.e., µ = µ 0 ) and best fit condition (such that R2~max; i.e., µ = µ 0 + µ t × sample size/2) and variation in time (years) of the location parameter of the nonstationary GEV models (c,f,i,l,o). Dark dashed lines of GEV models graph depict 95% confidence interval and provides an overview of the uncertainty associated to estimated flood magnitude at specific return period in the initial condition (i.e., µ = µ 0 ). As expected, low frequency flood has much significant uncertainty associated to magnitude estimation.

Flood Hazard Mapping Accounting for Multiple Probability Model
For converting databases of probability of exceedance in high resolution flood hazard maps, the procedure and code proposed by Guerriero et al. [27] was used. This procedure consists in deriving a flood inundation model of the study area from topographic data and subsequently in using an interpolation/substitution routine to spatialize each probability of exceedance database using the inundation model as reference. The adopted routine iteratively extracts each line of the probability database containing the reference fluvial stage and associated probability for each considered station and interpolates probability values across the study area, on the basis of the position of each station creating a probability raster. Once spatialized, the interpolated values of probability of exceedance are iteratively assigned to corresponding cells in terms of X and Y coordinated of the inundation model that contain the modeled fluvial stage.
In this case, the flood inundation model was computed as a difference between the available DSM and the water surface digital elevation model (wDEM). The wDEM was computed through the interpolation of a triangulated irregular network of the water elevation data at cross-sections defined along the river reaches following indication of [50] Figure 1. The obtained inundation model was normalized to the maximum magnitude of the modeled 500-years flood and subsequently used in flood hazard spatialization through the procedure proposed by Guerriero et al. [27].
Due to the presence of multiple probability of exceedance databases corresponding to past, current and near future flood hazard scenarios, multiple flood hazard maps were produced and classified. Flood hazard maps classification was completed identifying specific hazard zones: a (i) very high hazard zone, floodable by a 1-year to 5-year floods, a (ii) high hazard zone, floodable by a 5-year to 30-year floods, a (iii) medium hazard zone, floodable by a 30-year to 100-year floods and a (iv) low hazard zone, floodable by a 100-year to 500-year floods. Flood magnitude for each specific return period was derived from probability analysis considering that the theoretical return period is the inverse of the probability of exceedance.

Validation
To validate results from the analysis, a comparison between the 225-years floods extent predicted by a specifically developed flood hazard scenario corresponding to the AD 2015 and that observed during the event of October 2015 across the lower sector of the Calore river by Guerriero et al. [22] see Figure 1 was completed. Such choice is related to the indication of probability models that the return period of the October 2015 event is around 225 years. To validate the model, the degree of matching of the observed vs. predicted flooded areas across this part of the study area was considered e.g., Molinari et al. [51]. Since Guerriero et al. [27] indicated that this flood hazard mapping approach can be considered an alternative also to hydraulic simulation-based flood hazard analyses, 2D hydraulic simulation of estimated hydrographs was not completed for validation purposes. Table 1 reports result from time series trend estimation using a linear regression fit. As indicated in the table, time series present a downward trend characterized by a slope coefficient that range between −0.040 and −0.010. Results from p-Value analysis indicate statistical significance of trends in the analyzed time series. It is important to notice that time series recorded at the Benevento station is characterized by a p-Value equal to the considered significance threshold. This result might be related to the presence of wider discontinuities in this time series. In addition, the coherence of the trend with that characterizing the further two stations installed along the Calore river, and the spatial variation of trend coefficients (slope) that mimics a consistent decrease down valley, should be considered as indicators of consistence of the significant trend detected by linear regression in the hydrometric record registered at the Benevento station.  (Table 1). Time series reported in Figure 2a,d,j provide an overview of the characteristics of the Calore river, while time series reported in graphs g and m provide an overview in the characteristics of the Sabato and Volturno rivers, respectively. In the case of the Calore river, a tendency to an increase of maximum fluvial stage down-valley is observable with a peak at the Solopaca station where a maximum level of 10 m has been observed. Figure 2a,d,g,j,m report also nonstationary GEV probability models as linear function, indicating time variation of the magnitude of 5-years, 10-years, 30-years, 50-years, 100-years, 200-years and 500-years events (multi-color lines). As expected, such variation is characterized by a linear trend having a slope equals to µ trend . Figure 2b,e,h,k,n show relative frequency distributions and obtained initial GEV models (i.e., µ = µ 0 ). Such models offer an overview of the potential magnitude of the flood at different return periods for the area of study. Especially, for the Calore river, the magnitude of the 10-years flood varies between approximately 4.2 m and 8.3 m and the magnitude of the 100-years flood between approximately 6 m and 11.3 m. As indicated by 95% confidence intervals (i.e., dashed lines), such estimation is characterized by larger uncertainty in comparison with that of the station located along the Sabato and the Volturno rivers. Along the Sabato river, the magnitude of the 10-years and 100-years floods is estimated in 3.4 m and 4.9 m, respectively. Along the Volturno river, the magnitude of these two potential events are estimated in approximately 4.4 m and 5.3 m, respectively. Such graphs also report the probability models calculated considering a location parameter µ such that R2~max. Overall, the R2 is higher than 0.94. Figure 2c,f,i,l show variation of the GEV distributions location parameter µ with time as modeled by a linear trend function in the reference period (i.e., AD 2000 to AD 2050). Overall, the location parameters follow the downward trend characterizing the time series with equal slope coefficient. Figures 3 and 4 depict flood hazard evolution in the study area between AD 2000 and AD 2050. Figure 3a,b depict zonation in flood perspective referred to a past and a present scenario, corresponding to AD 2000 and AD 2020, respectively. Figure 4a,b depict zonation in flood perspective referred to two future scenarios, a first one (a) corresponding to AD 2030 and a second one (b) corresponding to AD 2050. A general increase with time of extent of areas floodable by high return period events is observable from the maps. This increase is consistently related to the decrease of the extent area floodable by low return period events.

Maps reported in
The past flood scenario (i.e., AD 2000) indicates how a large fraction of the floodplains belonging to the study area were classified as a high to very high hazard zones. Zones of medium to low hazard were poorly developed. The industrial area of Benevento was an exception since a significant part of this area was floodable by a 100-years to 500-years event. Moving from the past to the present (i.e., AD 2020) and then to the two future scenarios (i.e., AD 2030 and AD 2050), an increase of the zones floodable by events characterized by 30-years to 100-years and 100-years to 500-years flood events is observable. This is consistent with the downward trend observed in time series. These maps offer a perspective of the evolution of the hydrologic response of the floodplain to floods.
Maps of Figure 5 shows flood hazard at selected locations of the study area. Three zones of particular significance in terms of human influence on river dynamics e.g., Di Baldassarre et al. [52] are considered. The zones are (1) the manufacturing area of Benevento (Figure 5a,d,g,i), (2) a sector its urban and sub-urban area (Figure 5b,e,h,m) and, (3) the lower part of the Calore river floodplain (Figure 5c,f,i,n). From these maps, an increase of the zones floodable by events characterized by 30-years to 100-years and 100-years to 500-years flood events in time is observable.   The past flood scenario (i.e., AD 2000) indicates how a large fraction of the floodplains belonging to the study area were classified as a high to very high hazard zones. Zones of medium to low hazard were poorly developed. The industrial area of Benevento was an exception since a significant part of this area was floodable by a 100-years to 500-years event. Moving from the past to the present (i.e., AD 2020) and then to the two future scenarios (i.e., AD 2030 and AD 2050), an increase of the zones floodable by events characterized by 30-years to 100-years and 100-years to 500-years flood events is observable. This is consistent with the downward trend observed in time series. These maps offer a perspective of the evolution of the hydrologic response of the floodplain to floods.

Validation
Results of validation analysis are depicted in Figure 6a. Especially, a comparison between data by Guerriero et al. [22] representative of the area inundated by the 2015 flood event and the estimated 225-years flood scenarios is presented. The degree of matching between the area floodable by a 225year event estimated by the model with that affected by the 2015 flood, is of 98%. At specific locations, the hazard model tends to slightly overestimate the extent of the flooded area introducing false positive in the estimation. The zones that suffer from major differences in estimation are the confluence between the Calore and the Volturno rivers and where the morphology of the river valley is particularly complex like the zone at the beginning of the area selected for the validation.  Figure 1). Figure (a-c) show the past flood hazard scenario (i.e., AD 2000) across the area covered by polygons a, b and c, respectively, in Figure 1. Figure (d-f) show the present flood hazard scenario (i.e., AD 2020). Figure (g-i) and (l-n) show the future flood hazard scenarios (i.e., AD 2030 and AD 2050, respectively).

Validation
Results of validation analysis are depicted in Figure 6a. Especially, a comparison between data by Guerriero et al. [22] representative of the area inundated by the 2015 flood event and the estimated 225-years flood scenarios is presented. The degree of matching between the area floodable by a 225-year event estimated by the model with that affected by the 2015 flood, is of 98%. At specific locations, the hazard model tends to slightly overestimate the extent of the flooded area introducing false positive in the estimation. The zones that suffer from major differences in estimation are the confluence between the Calore and the Volturno rivers and where the morphology of the river valley is particularly complex like the zone at the beginning of the area selected for the validation.

Discussion
The analysis suggests that the considered time series have a statistically significant downward trend. Such negative trend is consistent with the observed change in river discharge in southern Europe as consequence of the ongoing warming [10], thus it is potentially the signature of the ongoing climate change. In addition, authors of Magliulo et al. [53] indicated that, before AD 2000, the Calore river experienced incision, narrowing and reduction in the fluvial bars area. These changes, consistent with the negative trend in river discharge, were potentially related also to the increased presence of human infrastructures and activities. However, the analysis provides only a statistical overview of

Discussion
The analysis suggests that the considered time series have a statistically significant downward trend. Such negative trend is consistent with the observed change in river discharge in southern Europe as consequence of the ongoing warming [10], thus it is potentially the signature of the ongoing climate change. In addition, authors of Magliulo et al. [53] indicated that, before AD 2000, the Calore river experienced incision, narrowing and reduction in the fluvial bars area. These changes, consistent with the negative trend in river discharge, were potentially related also to the increased presence of human infrastructures and activities. However, the analysis provides only a statistical overview of trend that might be difficult to associate to the physics of the system [54]. For the selected study area, this is connected to the presence of discontinuities in the time series and the absence of reliable information about past human modification especially over the long-term. In this condition, being worthwhile the use of nonstationary approach for hazard analysis [36], the use of frequency analysis of annual maxima, enhanced by time dependent modelling of fit parameters (or other covariates like teleconnections; Coles [34], Gilroy and McCuen [55] and Lima et al. [56], can be recommended because the inference about the nature of the trend is included in the model and does not need to be necessarily interpreted [57]. In the specific case, only the location parameter of the probability distribution was considered to change with time. Although this is the simplest way to account for non-stationarity, Šraj et al. [32] indicated that, in estimation of flood magnitude, this approach has smaller uncertainty in comparison with those that use multiple time-varying parameters. In addition, a trend in a time series, may not be a direct consequence of a deterministic process, since other statistical behaviors such as long-range dependence (LRD) may induce trend-like behaviors [58]. This is the case of long-term climatic oscillation that might be not fully related to the ongoing climate change and considerable a fully stochastic process. On this basis, in the presented analysis, the prediction of future flood hazard was limited to the next 30 years in which the representativeness of results can be speculated.
Results from nonstationary extreme-value analysis indicate that the magnitude of the extreme flood scenario (i.e., 500-years event) is going to decrease by 10% to 30% between AD 2000 and AD 2050. In this perspective, the GEV probability function is suitable to describe the statistical behavior of annual maxima, also under nonstationary conditions [59], but tends to underestimate the intensity of very high return period events [34]. Thus, a different probability model might be selected for improving the analysis and uncertainty related to flood magnitude estimation needs to be evaluated and considered in future predictions [47]. Computed confidence intervals indicated substantial uncertainty in magnitude estimation of low return period events, also in relation to the limited sample number of available time series, the presence of discontinuity, and the use of the MLE estimator in model parameterization. Several authors indicated that the use of the Bayesian approach (MCMC estimator) often provides a better estimate of model parameters (i.e., narrower confidence interval; Šraj et al. [32] and Parkers and Demeritt [60]. In addition, Katz et al. [61] indicated that, with short time series, MLE approach might provide unrealistic estimate of the shape parameter. Uncertainty estimation has, then, to be carefully taken into account since unrealistic parameterization of probability models can lead to an under/overestimation of flood magnitude, especially for low probability events [47]. Having say that, validation analysis indicated a very high accuracy of the hazard model that, associated with its intrinsic high resolution, make hazard estimation more reliable also in presence of urban areas [62]. A key point of the used framework and its application is the possibility of deriving flood hazard maps from nonstationary probability models through a sort of "classification" of a LiDAR derived flood inundation model [27]. This is related to the use of a specific procedure of flood hazard spatialization [27] that, being able to consider multiple probability-models in a high dimensionality framework (i.e., space and time [63]), is able to work at local to regional scale.
The maps indicate that, overall, flood hazard is going to decrease across the study area in the next 30 years and many areas would require new guidelines as the hazard level decreases. This is an important implication derived by the adoption of a nonstationary approach for time series analysis. A comparison between the current flood hazard scenario and the floodplain zonation in flood perspective adopted by the Basin Authority of southern Apennine of Italy (Figure 6b) indicate substantial differences in terms of hazard estimation, especially in the identification of zones floodable by events characterized by a return period between 100 years and 500 years. In addition, the Basin Authority zonation is substantially based on the identification of two zones, a first floodable by a 100-years event and a second floodable by a 300-years event.
Conversely, the produced hazard maps are able to better describe flood hazard assigning to each cell of the inundation model the return period of the event that can involve that specific area, providing distributed information across the study area. The production of multiple scenarios describing also what is likely to happened in the near future extends the descriptive capability of the proposed method, that should be able to support a long-term plan of adaptation for scenarios of hazard change reducing expected annual damage in flood prone areas [64].

Conclusions
In this paper flood hazard across the central sector of the Benevento province in southern Italy was analyzed accounting for non-stationarity of available time series, proposing a reduced complexity framework for flood hazard analysis that account for the ongoing climate change. The procedure, modified from that proposed by Guerriero et al. [27], uses multiple probability models and a LiDAR derived high resolution inundation model to provide present and future flood scenarios in the form of hazard maps. The possibility of producing prospective flood hazard scenarios is related to the non-stationarity of records, consequence of the ongoing climate change and associated anthropogenic changes. In the proposed framework, non-stationarity induced by the ongoing climate change associated with human effects is considered by adopting standard extreme value analysis, enhanced by time dependent modelling of the location parameter. In this way, among the available hydrological records, only time series having a relevant length have been selected (i.e., samples > 40). Flood hazard maps are produced through a routine of spatialization of stage probability across the inundation model that is able to work at different scale.
Hazard zonation maps provide an overview of the flood hazard in this area, already hit by a number of destructive flood events, showing how riverine flood hazard is going to change in the next 30 years. Especially, according to global river discharge projections, in this part of the Europe, flood hazard is going to decrease in the near future with a decrease in magnitude of extreme flood events (i.e., 500-years flood) of between 10% and 30%. Consequently, many areas would require new guidelines of use as the hazard level decreases. For instance, restrictions along floodplains might be re-modulated on a decadal basis and new guidance in the context of the agricultural production adaptation might be provided. This is an important implication derived by the adoption of a nonstationary approach for time series analysis underlining the relevance of the proposed framework. In this way, predictions are able to support a long-term plan of adaptation for both scenarios of hazard change reducing expected annual damage in flood prone areas.
Limitations of the proposed analysis are related to the use of GEV probability model that tends to underestimate the magnitude of low frequency events and the use of MLE approach for model parameterization that might introduce ambiguity in probability model shape. However, result validation indicates a very high accuracy of the proposed procedure with a degree of matching with a recently observed 225-years flood of 98%. On this basis, the proposed framework, representing an evolution of that proposed by Guerriero et al. [27], can be considered relevant in flood hazard estimation due to its ability of predicting near future evolution of flood hazard as modulated by the ongoing climate change.