Introduction

The coastal cities worldwide, and more so in the Mediterranean, are threatened by sea-level rise (SLR), inundation, storm surges, and precipitation variability (Essink et al. 2010; Nicholls et al. 2007). On top of that, land subsidence poses a severe threat to coastal areas by intensifying the level of exposure to inundation or flood hazards (Galloway et al. 1999; Mimura 2021; Ogie et al. 2018). Climate change and land subsidence co-occurring cause instability of low-lying areas, coastal inundation and erosion, and saltwater seepage from the deepest portion of the aquifer with consequent contamination of freshwater resources (Abdoulhalik and Ahmed 2018; Nicholls et al. 2021; Sestini 1996). On this note, coastal cities and the adjacent rural areas often rely on either non-structural or structural measures (engineering structures) or both to mitigate inundation and protect a fragile water infrastructure consisting of drainage channels, dams, levees, dikes, detention reservoirs, seawalls, polders, pumping stations, floodgates, and more (Nazarnia et al. 2020; Ogie et al. 2018).

Climate change and land subsidence have a profound role in the severity of coastal flooding, with far-reaching water management, flooding infrastructure, ecological and socio-economic implications for low-lying coastal cities and polder systems. These two main drivers of coastal flooding have attracted the interest of many scientists and researchers in recent years (Essink et al. 2010; Mimura 2021). Recognizing the role these two drivers have on coastal surface water drainage through time is imperative to strategic water planning and development of coastal urban areas and their drainage waterways. Past and recent studies have greatly improved our understanding of these drivers, but correlation with seawater seepage and drainage is widely unrecognized in many coastal studies.

Coastal regions with low elevations and high populations are at the highest risk of inundation (Mimura 2021; Syvitski et al. 2009). The appealing nature of the coastal areas has accelerated human migration to the area, thus promoting rapid urbanization, tourist resort growth, and intensive farming. In the process, the natural environment and hydrological cycle are modified in many ways, viz drainage, land reclamation, discharge of sewage contaminants into coastal waters, river modification, sand mining, hydrocarbon production, and more. These kinds of human interventions promote saline intrusion of the shallow coastal aquifer and have a long-term impact on freshwater management and water infrastructure (flooding), which could be costly to mitigate in the future (Nicholls et al. 2007, 2021 ).

In the Ravenna coastal basins (Northern Italy), a surface water drainage system and pumping stations are necessary to manage floodwater and lower water table heads (Giambastiani et al. 2007). The initial objective of pumping was land reclamation for allowing agricultural activities, protect infrastructure, accommodate human settlement (Antonellini et al. 2015; Giambastiani et al. 2020; Teatini et al. 2005, 2006). Subsequently, land drainage and pumping, managed by the Land Reclamation Consortium (LRC), allows farmland activities by controlling vadose zone depth and maintaining constant water table depth in the range of 1.5–2.0 m below the ground during the year (Antonellini et al. 2015, 2019).

Relative Sea Level Rise (RSLR), which includes land subsidence, both natural and anthropogenic, and eustacy, has occurred across the entire coastal areas of the Northern Adriatic Sea and has produced relative ground settlements of meters in Ravenna (Carbognin et al. 2011). The area has suffered from an intensive subsidence rate, up to a maximum value of 110 mm/year in the 1970s, during the high economic boom due to deep groundwater natural gas withdrawals (Teatini et al. 2005). More recently, the rates of subsidence have significantly decreased in many areas due to strong mitigation law imposed on groundwater exploitation (Carminati and Martinelli 2002), but they still remain important in some places like the coastal area (Artese et al. 2016; Bitelli et al. 2015; Cerenzia et al. 2016; Teatini et al. 2005). The subsiding land raised the need to pump more water out of the coastal basins to keep pace with the rising groundwater level (Antonellini et al. 2019). At the same time, freshwater availability in the area is becoming scarce due to limited rainfall, severe drainage, saline intrusion, and long periods of drought (Antonellini et al. 2008). The climate and weather variability in the area has intensified seasonal imbalance in groundwater budget, characterized by water deficit in the summer and surplus in winter when drainage is at its maximum peak (Benini et al. 2016; Greggio et al. 2018).

The general aim of our paper is to explore correlations between these aspects and to show how drainage in different coastal basins may be more or less sensitive to land subsidence, water management practices, and climate variability (i.e., droughts, climate extremes, etc.). In this context, our work adds novelty and is important for coastal water management. In order to achieve our objective, we compute local climate extreme indices as part of the necessary analysis required to understand drainage evolution. The study prioritizes the need to better understand drainage evolution in coastal basins under long-term and seasonal conditions. A previous study (Giambastiani et al. 2020) on factors affecting water drainage long-time series in the Quinto coastal basin of Ravenna (Fig. 1) proves that the amount of extra water entering the aquifer is due to vertical seepage from the bottom of the upper sandy portion of the surface aquifer. This extra amount of water is caused by an increase in hydraulic head difference driven by land subsidence rate, while rainfall time series analysis (1971–2017) and land-use change do not show a well-defined trend and play no important role in justifying the increasing drainage trend. More study is required to understand these above-mentioned relationships and the need to understand different coastal basins behavior to drought indices and other important extremes in a similar timeframe prompted the current research work.

Fig. 1
figure 1

Location maps of the study area (a) and the Emilia-Romagna Region (b) the three low-lying coastal basins, along with the main hydrology, pumping stations and surface drainage network are shown in (c)

Full size image

The novelty of our study, however, is showing that there is no single process controlling a change in drainage but there are many and complexly interrelated. Land subsidence, climate variability and extremes (droughts) as well as past management practices all play an essential role. We show that there is a wide range of factors that affect drainage time series and knowledge gained from a single basin cannot be extrapolated to the whole coastal area.

The study-specific objectives are: (1) to examine the behavior of long-term drainage data in the lowest 3 coastal basins of Ravenna and further compare it with the hydroclimatic component, and (2) to identify the main factors contributing to surface water drainage through time in the area by considering a wide range of climate extremes, which include both precipitation and temperature indices. The connection highlighted in our study area provides a clear and complete description of the processes occurring in low-lying coastal areas subjected to drainage, subsidence, and salinization and could serve as a motivation and reference for other similar coastal basins worldwide. Moreover, the findings can help to take management actions for the sustainable use of water resources in coastal areas facing similar problems.

The paper is structured in a first general description of the study area, including the water management and drainage system, hydrogeological, and climate setting (Sect. “General description of the area”). Following, Sect. “Data acquisitions and elaboration” describes the methodology applied for drainage, climate and subsidence time series analysis and modelling, as well as water budget calculation to examine the effect of change in hydrologic cycle component on the water resources of the studied basins. In Sect. “Results and discussion”, the drainage trends are compared and discussed with regards to climate extremes and subsidence trends; finally, the main findings are reported in the conclusions (Sect. “Summary and Conclusions”).

Materials and methods

General description of the area

The city of Ravenna is in a low-lying coastal plain situated at the eastern end of the Emilia-Romagna region (North-Eastern Italy), about 60 km south of the Po River delta (Fig. 1). Ravenna is one of the oldest Italian towns founded in the VII century BC; it harbors great historical monuments that attract the interest of tourists all around the world (Teatini et al. 2005). The presence of artistic monuments and proximity to the sea has promoted a wide range of socio-economic and land-use activities in the area. The city land elevation has a range of 0.5 to − 0.5 m a.s.l. (Greggio et al. 2018). Currently, the Ravenna coastal area is composed of wetlands, lagoons, urbanized areas with industrial facilities, reclaimed agricultural land, rivers, canals, and complex drainage network systems (Giambastiani et al. 2008; Greggio et al. 2012).

The studied coastal area includes three mechanically drained coastal basins: San Vitale Basin (1,080 ha), Rasponi Basin (2,639 ha), and Quinto Basin (9,200 ha). The Quinto basin is the largest low-lying basins in the area, with a surface covering 10% of the Ravenna province, whereas the San Vitale and the Rasponi basins cover approximately 3% and 1%, respectively. The Quinto basin land use includes farmland, pine forests, beaches, and gravel quarries (Giambastiani et al. 2020; Mollema et al. 2012). The Rasponi basin is mostly supporting agricultural activities, while the San Vitale basin is a small-reclaimed area with natural and human-made features such as the pine forests, recreational areas, coastal lagoons, and farmlands (Giambastiani et al. 2007). Together, they cover the entire southern coastal zone of Ravenna.

Water management and drainage system

Water management along the Po river plain has always been a necessity for effective and sustainable farming since 2000 years ago (Antonellini et al. 2015; Greggio et al. 2018). The Romans used hydraulic techniques to distribute water from the Apennines down to the Adriatic Sea (Port of Ravenna) through mountainous, hilly, and lowland territories characterized by different problems and characteristics. The hydraulic interventions on the natural river networks, such as the digging of canals and river dikes construction, had the objective to ensure efficient use of land for farming (Stefani and Zuppiroli 2010). Nowadays, the entire stretch of the Emilia-Romagna plain has more than 30 pumping stations with a network of ditches, canals, and rivers interconnected to one another (Antonellini et al. 2019).

The Ravenna municipality accommodates nine active and three emergency pumping stations. The study area includes three mechanical pumping stations (called \(P1\), \(P2\), and \(P3\) in Fig. 1). The extensive land reclamation and impoldering system in the area emerged in the twentieth century when the existing hydraulic infrastructure was designed and constructed (1920 and 1960) (Vandenbohede et al. 2014). Due to the low topography, the water lifted from the basins is discharged directly into the Adriatic sea by the connected pumps and drainage channels (Giambastiani et al. 2020; Greggio et al. 2018; Stefani and Vincenzi 2005). During the summer period (Jun–Sept), the drainage canals receive surface water from adjoining rivers through a system of connected pipes and sluice gates. This method of water distribution system promotes efficient water management for drainage and irrigation distribution to support rural, industrial, recreational, and urban activities. At the same time, the Land Reclamation Authority manages water scarcity emergencies in summer and drought periods through several water-saving techniques such as the reduction in the withdrawal of surface water, rationalization of water distribution (Greggio et al. 2018), and inter-basin water transfers (IBWT) when and where necessary.

Hydrogeological settings

The geological setting of the studied area is composed of a sequence of layered sediments deposited during the Quaternary and Upper Pliocene, which represent different sedimentary environments (Gambolati et al. 1991; Teatini et al. 2005). The shallow coastal aquifer is composed of Holocene sediments, and the Pleistocene continental clay below forms the impermeable bottom confining layer. The coastal sediment of the upper shallow aquifer was transported and deposited by sea-level fluctuations (Amorosi et al. 1999; van Straaten 1970) due to the barrier-lagoon-estuary natural dynamic behavior (Fig. 2). The shallow coastal aquifer sedimentology includes wedge-shaped dunes and beach deposits intercalated with fine sediment and frequent thin (centimeters) peat layers formed in lagoons, marshes, and alluvial plains. This dynamic nature of the coastline resulted in alternating layers of sediments that correspond to multiple parallel dune belts formation (Giambastiani et al. 2007).

Fig. 2
figure 2

Image modified from Mollema et al. (2012)

Geological cross-section showing the stratigraphy and lithology of the Ravenna shallow coastal aquifer.

Full size image

In particular, the shallow coastal aquifer is composed of two sandy units separated by alternating clayey and sandy loam layers (prodelta deposits). The sandy units have relatively thick medium-grain sand in the upper part of the aquifer (0–10 m thickness) and fine-grained sand (− 21 to − 26 m thickness) in the lower part (Giambastiani et al. 2020; Mollema et al. 2012). Further inland, the two central sand units are laterally connected to gravel deposits, which separate them from the lagoon deposits. The silty-clay basement is at a depth varying from − 15 to − 20 m in the west and from − 25 to − 30 m in the east (present shoreline) (Amorosi et al. 2008; Campo et al. 2017).

Climate of the area

The climate of the Emilia-Romagna plain is peculiar because of the bordering mountains and the sea—the Apennines to the south, the Alps to the north, and the Adriatic Sea to the east. The orography and marine area are significant drivers influencing the spatial distribution of temperature and precipitation in the region (Tomozeiu et al. 2006), particularly for the coastal cities (ISPRA 2013; Tomozeiu et al. 2006).

During the past 48 years, the annual precipitation of the area varies between 300.2 and 1063.3 mm, with an average yearly value of 639 mm/year and the mean annual temperature varies from 10.8 °C to 15.5 °C with an average value of 13.5 °C/year. The annual rainfall was the lowest in 1987 (~ 489 mm) and the highest in 1997 (~ 774 mm). In the past decade, the rainfall over the area reaches its maximum during the autumn.

Data acquisitions and elaboration

This study focuses on the use of different statistical methods to analyse the coastal drainage data (Time Series Analysis, Local Regression Model), DEM elaborations for land subsidence mapping, and indexes computed from daily meteorological data from 1971 to 2018. We also elaborate sub-periods water budgets of the coastal basins. A flowchart diagram of the methodologies is in Fig. 3.

Fig. 3
figure 3

Flowchart of the methodology used to assess the potential threats caused by land subdidence, and climate variability and change in the study area

Full size image

Drainage data

The pumping data of coastal watersheds (San Vitale, Rasponi, and Quinto basins) were obtained from the original bookkeeping records of LRC from 1971 to 2018. It is worth mentioning that data are not available between 1988 and 1991 for both Rasponi and San Vitale basins. The daily hours of operation were recorded in a logbook sheet of each basin and digitized for future elaboration. Computationally, the daily drained water \(\left( {{\text{m}}^{3} /{\text{day}}} \right)\) of each station was estimated by multiplying the flow rate of each pump \(\left( {{\text{m}}^{3} /{\text{h}}} \right)\) by the duration of its operation \(\left( {\text{h}} \right)\). Afterward, the equivalent depth of the water columns \(\left( {{\text{m}}/{\text{day}}} \right)\) in each of the basins was obtained by dividing the total volume of water drained per day by the drained area of the individual basin \(\left( {{\text{m}}^{2} } \right)\).

Subsidence data

The subsidence surveys of each coastal basins were obtained from literature data published by Teatini et al. (2005) from 1972 to 1992 in the following sub-periods: 1972–1977, 1977–1982, 1982–1986, 1986–1992; and from the Regional Agency for Prevention, Environment, and Energy of Emilia-Romagna (Arpae website; www.arpae.it/index.asp?idlivello=1414) from 1992 to 2016 for sub-periods: 1992–2000, 2002–2006, 2006–2011, and 2011–2016. The subsidence data were acquired with topographic and Interferometric Synthetic Aperture Radar (InSAR) techniques calibrated on a 2005 topographic survey and a GPS permanent network (Cerenzia et al. 2016; Teatini et al. 2005). The data were elaborated with the QGIS software to obtain subsidence maps and isokinetics (mm/year) for the different sub-periods.

Climate data and climate extremes

The daily meteorological data—minimum, maximum, and mean of precipitation (\(P)\) and temperature (\(T)\)—were taken from different meteorological stations within the Ravenna area to get the longest time dataset possible. The assumption is that the dataset obtained from each weather station has insignificant differences among individual stations, because the stations are a few km apart. The dataset covers a period from 1971 to 2018 and it is relevant for the three coastal basins studied. The daily dataset for the periods 1971–1975 and 1976–2018 was retrieved from the Arpae periodical Annals (https://www.arpae.it/) and the Arpae DEXT3R website (https://simc.arpae.it/dext3r/), respectively.

The core set of descriptive climate extreme indices used in our study was developed by the joint \({\text{CCl/CLIVAR/JCOMM}}\) Expert Team on Climate Change Detection and Indices (ETCCDI) in 1999. The ETTCDI indices include temperature and precipitation indices. Indices calculation allows straightforward monitoring of trends that can potentially cause stress to coastal structures, humans, and the natural environment in a short timeframe (ISPRA 2013; Zhang et al. 2011). In achieving this, we used 13 temperature and 10 precipitation indices from the list of ETCCDI indices and two additional drought indices (Table 1).

Table 1 Descriptive table of the climate extreme indices selected
Full size table

The climate extremes indices were calculated using the RClimDex software developed by ETCCDI. The RClimDex is a \(R\)-based software (Zhang et al. 2018) that provides a user-friendly interface to compute indices of climate extremes and allows for user-defined input when calculating the indices. The drought indices were calculated using the \({\text{SPEI}}\) and \({\text{SPI}}\) package present in \({\text{R}}\) (https://www.r-project.org/). The \({\text{SPI}}\) and \({\text{SPEI}}\) calculation depends on the probability of events (such as \(P\) and \(T\)) transformed into a standardized series for the desired timescales ranging from 3 to 48 months (Hayes et al. 1999; McKee et al. 1993; WMO 2012). However, the most suitable timescale to monitor long-term hydrological changes is at least 24 months (such as reservoir level, groundwater, streamflow) (Hayes et al. 1999; WMO 2012).

Water budget calculation

Water budget evolution through time in the different coastal basins needs to be quantified to examine the effects of change in hydrologic cycle components on freshwater resources. The water budget analysis, in fact, is vital for assessing changes in freshwater availability within coastal basins. The irrigation and drainage activities in the coastal basins make water budget analysis complex to compute. The topography of the coastal basins is flat and the drainage system mainly collects the runoff water (Mollema et al. 2012). The general formulation of water budget is as follows:

$${\text{Input}} = {\text{Output}} \pm \Delta S$$
(1)

Potential evapotranspiration \(\left( {{\text{PET}}} \right)\) and actual evapotranspiration \(\left( {{\text{AET}}} \right)\) were obtained using the Thornthwaite equation (1948) and the Thornthwaite and Mather Procedure (1955,1957), respectively. The change in soil water storage (\(\Delta S\)) and infiltration \(\left( {{\text{Infl}}} \right)\) were also included in the calculation. The Thornthwaite-Maher procedure estimates for different parameters such as the soil moisture balance, recharge, water deficit (Steenhuis and Van Der Molen 1986).

While solving this equation, also vertical water seepage in the calculation must be considered. The vertical seepage in a polder is the movement of water under the pressure gradient existing between the aquifer bottom and the water table (which is below sea level in our case). For simplicity, here, we assume an upward seepage (\(Qs\)), although a small lateral flow component may be possible in certain conditions (Fig. 4). In general, the change in the amount of water leaving the system (\(\Delta {\text{Pu}}\) and \(\Delta {\text{AET}}\)) is equal to the change in the amount of water entering the system (\(\Delta P\) and \(\Delta {\text{Qs}}\)). In our case, the coastal aquifer is considered a closed system. The \(\Delta P\) and \(\Delta {\text{AET}}\) are assumed homogeneous throughout the coastal basins. Therefore, the change \(\left( {\Delta } \right)\) in the basins groundwater budget variations for different sub-periods (1971–1986, 1987–1996, 1997–2006, and 2007–2018) through time follows Eq. (3) below:

$$P + {\text{ Qs }} = {\text{ ETR }} + {\text{ Pu ,}}$$
(2)
$$\Delta {\text{Qs}} = \left( {\Delta {\text{ETR}} + \Delta {\text{Pu}}} \right){-} \Delta P ,$$
(3)

where, \(\Delta P\) is the change in precipitation (mm/year), \(\Delta {\text{Qs}}\) is the change in water seepage (mm/year), \(\Delta {\text{ETR}}\) is the change in actual evapotranspiration (mm/year), \(\Delta {\text{Pu}}\) is the change in pumping (mm/year).

Fig. 4
figure 4

Conceptual diagram showing the cross-section of the shallow coastal aquifer and the anthropogenic and natural processes at work in the Ravenna coastal basins

Full size image

The annual water budget components are random variables through time.

Timeseries analysis and modelling in R

We performed a time-series analysis of the \(P\), \({\text{PET}}\), \({\text{Infl,}}\) and \(\Delta S\) by applying the locally weighted scatterplot smoothing (LOESS) regression algorithm implemented in \({\text{R}}\) (https://www.r-project.org/). The fitting procedure uses a numerical algorithm that prescribes how the value of the regression surface \(g\left( {x_{i} } \right)\) is computed, which is the estimate of the regression surface \(g\) at a specific value of any point in the space of the predictors \(x\). In the fitting of local regression models, the regression surface properties and the errors require making choices on the specification of properties of the errors and the regression surface. The choices include gaussian or symmetric distribution, constant variance (a priori weights), locally linear or locally quadratic in numeric predictors, neighborhood size, normalization of the scales, dropping squares, and a conditionally parametric subset (Cleveland et al. 2007). However, this study uses identically distributed, Gaussian errors: one numeric predictor \(\left( {x_{i} } \right)\) local regression model. The Gaussian errors have a constant variance \(\sigma^{2}\), and working with one numeric predictor assumes \(x\) is any value along the scale of measurement of the variable. By definition, the local regression model relates both the response and predictors as follows:

$$y_{i} = g\left( {x_{i} } \right) + \varepsilon_{i} ,$$
(4)

where \(\varepsilon_{i}\) is a random error, and \(g\left( {x_{i} } \right)\) is the expected value of \(y_{i}\) (mm).

By letting \(\Delta_{i } \left( x \right) = \left| {x - x_{i} } \right|,\) let \(\Delta_{i} \left( x \right)\) be the values of these distances ordered from the smallest to the largest, and using the tricube weight function.

$$T(u;t) = f(x) = \left\{ {\begin{array}{*{20}l} {1 - \left( {\left( {u/t} \right)^{3} } \right)^{3} } \hfill & {{\text{for}} \,0 \le u < t} \hfill \\ 0 \hfill & {{\text{for}} \, 0 \ge t} \hfill \\ \end{array} } \right.$$
(5)

The smoothness of the LOESS fitting depends on the specification of the neighborhood parameter (\(\alpha > 0\); \(\alpha\) is the angle of inclination). As \(\alpha\) increases, \(g\) becomes smoother. Supposedly, \(\alpha < 0\), and the assumption is that \(g\) will be equal to the number of points in each smoothing neighborhood \(\left( {\alpha n} \right)\) truncated to an integer. Therefore, the defined weight for the \(\left( {x_{i} ,y_{i} } \right)\) is given as follows:

$$W_{i} \left( x \right) = T\left( {\Delta_{i} \left( x \right); \Delta_{\left( q \right)} \left( x \right)} \right)$$
(6)

When \(\alpha \le 1\), we let \(q\) be equal to \(\alpha n\) truncated to an integer. For \(\alpha > 1\), the \(W_{i} \left( x \right)\) are defined in the same way, but \(\Delta_{\left( q \right)} \left( x \right)\) component is replaced with \(\Delta_{\left( n \right)} \left( x \right)\alpha\). The neighborhood weights called \(W_{i} \left( x \right)\) either decrease or remain constant as \(x_{i}\) increases in distance from \(x\). Another parameter is \(\lambda\); a parametric function used to describe the regression surface in either linear or quadratic polynomial forms. For \(\lambda = 1\), the specification is linear, and if \(\lambda = 2\), the specification is quadratic (see Cleveland 1979; Cleveland et al. 2007; Cleveland and Devlin 1988; Loader 2012, for detailed explanation).

Results and discussion

Time-series analysis

The hydroclimatic parameters time-series represent the climate of the coastal basins (San Vitale, Quinto, and Rasponi basins). The LOESS model (Eq. 4) was used to produce average smooth polynomial lines on the annual time series. The changes in the modelled trends (Fig. 5 and Supplementary Material: Fig. S1) indicate the change in individual components over different sub-periods. The sub-period variation of individual components sums up to give a total water budget changing through time.

Fig. 5
figure 5

Annual time series data and modelling of pumping in the following: a San Vitale, b Rasponi, and c Quinto basins from 1971 to 2018 (in mm/year). The blue line is the locally fitted non-parametric curve used to produce a smoothed model surface with the LOESS method; the 95% confidence interval is shown in grey shading

Full size image

The annual precipitation rate in the area is characterised by periods of dryness and wetness (Supplementary Materials: Fig. S1). The driest period recorded in the coastal area is around the mid-1980s, whereas the wettest period is observed around 2000. In the 1970s–1980 s, the annual precipitation in most parts of Italy was negative in trend. The negative trend in rainfall led to long periods of dryness, resulting in severe freshwater shortage (Piervitali et al. 1998). In a similar period (1971–1986), the annual rainfall in the Ravenna coastal basins dropped to ~ 492 mm (a minimum value). Afterward, the amount of rainfall rose to ~ 760 mm in 1998, which is the highest peak of rainfall recorded in the area. Due to natural rainfall variability, its amount dropped again to ~ 593 mm in 2006.

The trend of AET (Supplementary Materials: Fig. S1) is dependent on the level of soil moisture content and the amount of rainfall. Therefore, the AET sinusoidal trend is similar to the P trend. The P and AET trends are marginally increasing with time in the area, while the PET trend (Supplementary Materials: Fig. S1) is significantly increasing due to a rise in average air temperature. Consequently, the high PET and evaporation in the coastal area will lead to increased water loss from open water surfaces, declining soil water content, and low drainage during the warming months.

Apropos to drainage evolution, between 1971 and 1983, the rate of drainage in all coastal basins is characterised by two apparent peaks seen in 1973 and 1980 (Fig. 5a–c). The pumping rate is high in this period, especially in the San Vitale and Rasponi basins. The reason for the high drainage rate is related to the high subsidence rate and incomplete wetlands drainage. In fact, Simeoni and Bondesan (1997) reported that large freshwater marshes and brackish wetlands present at the beginning of the twentieth century were close to being reclaimed in the 1970s. Another major landmark event between 1971 and 1983 was severe land subsidence (> 65 mm/year) that promoted inundation in the coastal basins (Carbognin et al. 2011; Carminati and Martinelli 2002; Gambolati et al. 1991; Teatini et al. 2005). At the same time, Adriatic climate-induced sea storms permanently flooded most parts of the Ravenna area in the late 1970s (Teatini et al. 2005). As a result, water management authorities decided to increase the drainage capacity by increasing the number of pumps to reduce flooding and water stagnation where the topography is low and experienced a high subsidence rate.

After 1983, the pumping rate dropped significantly (a minimum in 1987) in the Rasponi and San Vitale basins and remained low until around 1993. The reason for the sudden decline in drainage in all coastal basins is an extended period of drought between 1985 and 1995 (as shown in Fig. 9). By 1995, the pumping rate in the coastal basins has started to increase again due to a rise in precipitation until 2000. After that, the pumping rate increases steadily with a 22-mm rise in rainfall in a similar period.

Monthly trends

The rainfall amount received by the coastal basins often increases in the autumn (Sep–Nov) and the spring (Mar–May). The amount of water supplied to the root zone by rainfall (\(P \ge {\text{PET}}\)) is the largest during autumn. On the other hand, the amount of rainfall decreases during the summer, particularly in July. The amount of soil moisture stored within the root zone is minimal during the summer due to limited precipitation and high evapotranspiration rates in the area (Supplementary Materials: Fig. S2). The substantial increase in summer temperature further worsens the condition of dryness in the coastal watersheds.

Apropos to drainage, a considerable portion of rain falls during the autumn when drainage is still low. Autumn rainfall has to compensate for water loss during the summer periods (Fig. 6a, 7a, 8a). In the late fall and early winter months, the aquifer eventually gets saturated allowing extra water to reach the drains and the intake canals to the pumps. The rate of daily pumping reaches its peak in winter despite the declining winter rainfall. In summary, the rate of pumping in the three coastal basins depends on the antecedent precipitation and soil–water interaction. Giambastiani et al. (2020) already excluded the influence of the urbanization on drainage series in the Quinto basin and showed no increase or decrease in the correlation between rainfall and the equivalent amount of drained water over periods of minimum (1971–1981) and maximum urbanization (2006–2016).

Fig. 6
figure 6

Monthly boxplot showing: a pumping, b pumping/precipitation ratio, c infiltration, and d change in water storage in the San Vitale basin from 1971 to 2018. The red dot is the mean of the monthly observations

Full size image
Fig. 7
figure 7

Monthly boxplot showing: a pumping, b pumping/precipitation ratio, c infiltration, and d change in water storage in the Rasponi basin from 1971 to 2018. The red dot is the mean of the monthly observations

Full size image
Fig. 8
figure 8

Monthly boxplot showing: a pumping, b pumping/precipitation ratio, c infiltration, and d change in water storage in the Quinto basin from 1971 to 2018. The red dot is the mean of the monthly observations

Full size image

Currently, the coastal basins have a water deficit during the summer and water surplus during the winter when drainage is maximum. Water withdrawal (drainage) and rainfall recharge fluctuations cause stress changes within the unconfined aquifer units, similar to natural subsidence mechanisms (Chaussard et al. 2013). The lowering of water level and the soil moisture conditions determined by the intensity of drainage, have an effect on the degree of subsidence, especially where organic material layers (i.e. peat layers) and drained cultivated peatlands are present due to organic matter oxidation (Gambolati et al. 2006; Hoekstra et al. 2020; Wösten et al. 1997; Zanello et al. 2011). Moreover, Antonellini et al. (2019) showed that the elastic and delayed-elastic components of the natural subsidence in the shallow aquifer of this coastal area are closely related to water table fluctuations. These fluctuations change the effective stress at the daily (in the sandy unconfined portion of the aquifer) and seasonal time scales (in the organic-rich fine sediment of the semiconfined prodelta portion), respectively.

Further observation reveals that the infiltration rate (Infl) and change in water storage (ΔS) become noticeable at the end of the year, resulting in increased baseflow within the drainage canals. Land reclamation authorities often raise the pumping rate in winter to prevent unexpected flooding and waterlogging of the crops where the water table is shallow and close to the topography. However, the Pu/P ratio diagram (Fig. 6b, 7b, 8b) shows that the rate of pumping increases in the San Vitale and Quinto basins during the summer, particularly in July. This could be due to specific management practise applied on the two basins, where extensive agriculture practice and important irrigation occur. Here, many canals have a double use, as drainage and irrigation canals, especially during the irrigation season (spring and summer months). The increase in pumping in July could be explained by the return flow, that is the portion of the irrigation water that travels directly from the irrigated fields as surface runoff, as well as the water which has infiltrated into the soil, and part of its flows downstream.

The high summer Pu/P ratio value in the Rasponi basin is less apparent than other coastal basins (Fig. 7b).

Climate extreme analysis

Climate extreme indices calculation is fundamental to detect early climate anomalies and allow for adaptation strategies and prevent societal damages (Zhang et al. 2011). Climate extreme indices that are statistically significant lead to severe impact. Likewise, an insignificant trend does not mean that they are less likely to occur, because an extreme event may result from the accumulation of moderate rainfall events or interaction with other extreme conditions such as storm surge and exceptional tides (Seneviratne et al. 2012). However, most of the warming extremes indices examined are significantly increasing in trend (Table 2). These include \({\text{SU - 25}}\), \({\text{GSL}}\), \({\text{TXx}}\), \({\text{TX90p}}\), \({\text{TX10p}}\), \({\text{TXn}}\), \({\text{TN10p}}\), and \({\text{TN90p}}\), \({\text{CSDI}}\). Significant temperature trends imply that evaporation and evapotranspiration processes will increase in the area, thus partly affecting the water table elevation and the amount of water reaching the vadose aquifer. On the other hand, all extreme rainfall indices trends are weak and not as significant as the warming extremes are (Table 2). This is in line with previous works carried out on a sub-portion of the study area (Quinto basin, Giambastiani et al. 2020) and on the entire coastal aquifer of the Emilia-Romagna (Giambastiani et al. 2021). These studies show no clear rainfall trend and no change in intensity and extreme rainfall events in the last 50 years (1971–2017) in the Ravenna territory.

Table 2 Climate extreme indices for the low-lying coastal basins of Ravenna (1971 to 2018)
Full size table

Apropos to drainage, all trends of extreme rainfall indices slightly increase with a marginal rise in drainage rate except in the San Vitale basin, where they show a negative trend between 2000 and 2018. Based on rainfall trends and documented flood memories in the area (caused by land subsidence, exceptional tides, and storm surges), the coastal drainage infrastructures in the area have increased in numbers and capacity to efficiently and effectively accommodate the processes conveying water into the basin.

The \({\text{SPI}}\) and \({\text{SPEI}}\) drought indices indicate the relationships between drainage and rainfall deficit. The prolonged drought or rainfall deficit witnessed between 1985 and 1995 is seen to be the main factor influencing the drainage rate in all three coastal basins (Fig. 9). The most prolonged dry spell \(\left( {{\text{CDD}}} \right)\) and the longest wet spell \(\left( {{\text{CWD}}} \right)\) trends are weak and not significant. For this reason, we do not have enough evidence as to whether the \({\text{CWD}}\) and \({\text{CDD}}\) trends will increase or not in the near future. The \({\text{CDD}}\) and \({\text{CWD}}\) are related to the rainfall occurrences, and they can serve as drought indicators, which have great importance to agro-meteorologist and water managers (Casanueva et al. 2014).

Fig. 9
figure 9

Monthly series for drought indices in the low-lying areas: a showing the SPI, and b the SPEI index from 1971 to 2018. Positive SPI values indicate higher than median precipitation, whereas negative values indicate less than median precipitation.

Full size image

Subsidence rate

The rate of subsidence for the sub-period from 1972 to 1980 is the highest (> 65 mm/year) recorded in the past half-century (Supplementary Materials: Fig. S5). This coincides with the development of on-shore and off-shore deep gas exploitation coupled with groundwater withdrawals, which started in the early 1950s and went on until the late 1970s (Carminati et al. 2003).

The subsidence mapping of the low–lying coastal areas shows that the land portion farther from the coastline is more stable with respect to the areas close to the coastline that are sinking faster (~ 20 mm/year), particularly in the Quinto and Rasponi basins. The current rate of land subsidence in the San Vitale basin is smaller than in the other two basins (about 2.5 mm/year). The subsiding part of the coastal basins poses significant threats to coastal towns, hydraulic infrastructures, and an expansion of the drainage network (new pumping stations, more and more powerful pumps) to cope with a more considerable amount of water to be drained. The increase in hydraulic head difference, driven by land subsidence rates, causes vertical seepage from the bottom of the coastal aquifer and an amount of extra water entering the basin (Giambastiani et al. 2020). In this context, the higher sinking rates observed at the coastal edge of the Quinto and Rasponi basins have promoted more seawater influx into their coastal aquifers compared to the San Vitale basin (Supplementary Material: Fig. S5). These two basins correspond to the most saline portion of the Ravenna territory, with the lowest water table depth, as reported in the reported maps developed by Giambastiani et al. (2021) and related to the entire coastal aquifer of the Emilia-Romagna. The pumping rate is increasing in the Rasponi and the Quinto basins due to the constantly increasing subsidence rates. The San Vitale basin has been relatively stable in the past 15 years, as shown by the decline of drainage rate in this period. Cerenzia et al. (2016) reported that Marina di Ravenna, a town in the San Vitale basin, is recently more stable and emphasized that the forthcoming sub-periods subsidence mapping by Arpae will be useful to either deny or accept the trend.

Water budget analysis

The sub-period (1971–1986, 1987–1996, 1997–2006, and 2007–2018) water budget trends show the changes in seawater seepage (ΔQs), rainfall (ΔP), evapotranspiration (ΔET), and pumping (ΔPu) in the coastal basins (Supplementary Materials: Fig. 6a–c). The input components of each sub-period recharge include vertical seepage of seawater (Qs) and rainfall (P), while the output components include evapotranspiration (ET) and pumping (Pu) (Fig. 4). The Qs increases within the coastal aquifer when AET and Pu components are greater than P. On the other hand, the Qs decreases because the AET and Pu components are lower than P.

In the first sub-periods (1971–1986), pumping decreases in all the coastal basins (San Vitale, ΔPu = − 8.3 mm/year; Rasponi, ΔPu = − 7.8 mm/year; Quinto, ΔPu = − 2.1 mm/year) with decrease in rainfall (− 14.6 mm/year). The output components are greater than rainfall recharge, so vertical seepage through the shallow unconfined aquifer increases (San Vitale, ΔQs = 1.5 mm/year; Rasponi, ΔQs = 2.1 mm/year; Quinto, ΔQs = 7.7 mm/year). In this period, vertical seepage increment plays a more important role by contributing more to water table elevation with implications for freshwater availability.

During the next sub-period (1987–1996), pumping in the Quinto (ΔPu = 11.61 mm/year) and San Vitale basins (ΔPu = 15.4 mm/year) increases with an increase in rainfall recharge (31.2 mm/year), while it slightly decreases in the Rasponi basin (ΔPu = − 0.1 mm/year). Despite this, the output components are lower than rainfall in all coastal basins, and so the vertical seepage decreases everywhere (San Vitale, ΔQs = − 1.4 mm/year; Quinto, ΔQs = − 5.2 mm/year; Rasponi, ΔQs = − 17.6 mm/year). In this case, rainfall for this sub-period is enough to counterbalance saline intrusion because the sum of changes in processes that promote salinity (seepage, ET, and pumping) are lower than rainfall in this period.

From 1997 to 2006, pumping in the Rasponi and Quinto basin increases (Quinto, ΔPu = 2.3 mm/year; Rasponi, ΔPu = 0.8 mm/year) as rainfall decreases (− 20.1 mm/year). The current sub-period decline in rainfall is not as drastic as the 1971–1986 drop. In contrast, the pumping in the San Vitale basin decreases (ΔPu = − 6.4 mm/year). The output components are greater than the rainfall input and consequently vertical seepage flow into the unconfined aquifer increases in all basins (San Vitale, ΔQs = 4.5 mm/year; Quinto, ΔQs = 13.2 mm/year; Rasponi, ΔQs = 11.7 mm/year). Since the amount of water leaving the shallow aquifer is higher than the rainfall input, the upward-directed seawater seepage has contributed more to the water table increase in this period.

During the last period (2007 to 2018), pumping increases for both Quinto and Rasponi basin (Quinto, ΔPu = 5.2 mm/year, Rasponi, ΔPu = 9.6 mm/year). On the other hand, pumping continues to decline in the San Vitale basin from the last sub-period (current sub-period ΔPu = − 17.2 mm/year) despite the increase in rainfall. The output components are higher than rainfall in Rasponi, while they are lower for San Vitale and Quinto basins. Consequently, the vertical seepage in the Rasponi basin increases (ΔQs = 1.2 mm/year), while it decreases for the other two basins (Quinto, ΔQs = –3.2 mm/year; San Vitale, ΔQs = –25.7 mm/year). The current declining pumping rate in the San Vitale basin is remarkable, as also observed in the last sub-period 1997–2006. The significant reduction in upward seepage through the shallow unconfined aquifer of the San Vitale basin also explains the current reduction in pumping.

Over the last 48 years, the overall vertical seepage strongly decreases with a significant decrease in pumping in the San Vitale basin (Supplementary Materials: Fig. S6a); the same occurs in the Rasponi basin even if at a lower rate (Supplementary Materials: Fig. S6b). The positive vertical seepage in the Quinto basin, on the other hand, contributes more to water table rise and, as a result, an increase in drainage rate (Supplementary Materials: Fig. S6c). This agrees with the increased hydraulic vertical gradient calculated for the Quinto basin (period 1971–2017) by Giambastiani et al. (2020) and resulting in about 139–92 mm of a total equivalent amount of water due to vertical seepage. The increase or decrease in drainage and seepage can cause a similar change in freshwater availability and water table elevation within the coastal aquifer. However, rainfall and evapotranspiration also affect the water elevation table and change in freshwater availability; for instance, more rainfall recharge counteracts other processes that promote aquifer salinity such as seepage, drainage, and ET.

Summary and conclusions

The coastal areas worldwide are at risk for climate change effects and land subsidence problems. There is a need to adequately inform the coastal communities, coastal managers, and scientists interested in how climate extremes and land flooding could shape the future of freshwater availability and mitigation measures put in place. In this context, the conclusions addressed in this study could serve as a motivation and reference to understand other coastal areas undergoing similar conditions worldwide.

The results of this study show that the drainage behaviour in the coastal basins is influenced by land subsidence, past management practices, and climate variability. Climate change indicators (such as wet periods, drought periods) and land modification indicators (such as reclamation, groundwater and gas extraction) have profoundly controlled surface water drainage evolution. The most important indices to monitor drainage evolution are the drought indices: SPI and SPEI; they allow detecting the precipitation severity in terms of wetness and dryness and permit easy comparison of them with drainage evolution along time scale.

Land subsidence causes an increase in upward-directed vertical seepage of saline water into the coastal aquifer and the drainage network. This leads to an increase in drainage through time that promotes salinization and a decrease in freshwater availability. The vertical seawater seepage could cause significant long-term concern for ground and surface water users in the low-lying coastal basins.

The period characterized by the highest subsidence rate (1950–1970 s) had a more significant influence on drainage than in any other period and promoted more seawater influx into the aquifer. The upward-directed vertical seepage of saline water through the shallow unconfined aquifer in the most subsiding basin (Quinto basin) has enormously contributed to the increase in the drainage rate due to the extra water influx reaching the drainage canals. On the contrary, a decrease in land subsidence rates in the San Vitale basin has decreased seepage from the aquifer bottom, thus contributing significantly to the decline in pumping rate and water table rise.

On the other hand, drainage has proven to be an essential mitigation infrastructure in the low-lying coastal basins to mitigate potential inundation, but with a long-term impact on water resources management and planning. An effort should be made to control the depth of drainage by keeping the water table depth as small as possible, while maintaining agricultural crop requirements. Accurate planning and timely control of the drainage system and the pumping station activity can be used to obtain this control.

The coastal basins are currently characterized by water deficit during the summer and water surplus during the winter, when aquifer gets saturated by autumn rain and drainage is maximum. At the seasonal timescale, the rate of pumping depends on antecedent rainfall and soil–water storage. The abundant autumn rainfall does not correlate with an increase in drainage until the end of November or December, because the soil water storage needs to be restored before a rise in water table triggers drainage. The freshwater availability in the area is becoming scarce due to limited rainfall, persistent water drainage, and long periods of drought. Drainage activities in the coastal basins have prevented rainfall recharge and the vadose soil from receiving enough freshwater to counteract the saltwater intrusion.

The detailed analysis carried out in sub-periods for different mechanically drained coastal basins highlights the temporal (i.e., periods in which seepage does not exist) and spatial variability of the processes, and emerges that an upscaling cannot be done without considering the water budgets of the individual basins. The novelty of our study is showing that there is no single process controlling a change in drainage but there are many and complexly interrelated. Land subsidence, climate variability, climate extremes, as well as past management practices all play an important role. We show that there is a wide range of factors that affect drainage time series and knowledge gained from a single basin cannot be extrapolated to the whole coastal area. Therefore, water management practices in this coastal zone should be tailored to the specific characteristics (hydrologic balance response) of the basin and focus on those processes that mostly affect drainage.