Next Article in Journal
Analysis of the Effect of Uncertainty in Rainfall-Runoff Models on Simulation Results Using a Simple Uncertainty-Screening Method
Previous Article in Journal
Risk Evaluation Model of Life Loss Caused by Dam-Break Flood and Its Application
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

The Analysis of Short-Term Dataset of Water Stable Isotopes Provides Information on Hydrological Processes Occurring in Large Catchments from the Northern Italian Apennines

1
Scientific High School Aldo Moro, 42124 Reggio Emilia, Italy
2
Geoinvest Srl, 49122 Piacenza, Italy
3
ARPAE Regional Agency for Environmental Protection, 42122 Reggio Emilia, Italy
*
Author to whom correspondence should be addressed.
Water 2019, 11(7), 1360; https://doi.org/10.3390/w11071360
Submission received: 17 May 2019 / Revised: 15 June 2019 / Accepted: 25 June 2019 / Published: 30 June 2019
(This article belongs to the Section Hydrology)

Abstract

:
This study discusses a dataset of water stable isotopes from precipitation (4 rain gauges) and surficial water (9 rivers) from the northern Italian Apennines, an area in which clay-rich bedrocks widely outcrop and the runoff response to precipitation events is very rapid. The dataset has been compiled starting from existing data that had previously been published in the literature and consists of monthly values of stable isotopes oxygen-18 (18O) and deuterium (2H) lasting over the period from January 2003 to December 2006 (precipitation) and from January 2006 to December 2007 (surficial water). For this period, mean residence times estimated by means of a sine-wave fitting technique make evident the significant differences over time spent by water molecules within the 9 catchments. Moreover, isotopic compositions of rivers deviated from those of precipitations revealing the influence of some catchment characteristics in differentiating the isotopic composition in rivers. Further correlations between mean residence times of river water and selected catchment characteristics reveal the role of orography and bedrocks in delaying the water molecules during their flow-paths. In addition, time series and cross–correlation analyses indicate a certain control by the main watershed divide on the isotopic composition of river waters, which is reflected in a progressive isotopic variation with longitude. The study shows that, despite using a short-time dataset (2-years for surficial water) of sparse stable isotopes can provide remarkable indications for depicting hydrological processes in large catchments made up of clay-rich bedrocks.

1. Introduction

Catchments are complex hydrological systems in which the understanding and quantification of the several processes participating in the final runoff is a challenge [1]. This is because the hydrological processes are heterogeneous and their full understanding is not achieved even on small scales (i.e., that of the slope). Thus, extrapolating hydrological processes occurring at the slope scale to a larger scale (i.e., that of catchment) is difficult. Among others, oxygen (18O and 16O) and hydrogen isotopes (2H and 1H) compose the water molecules and have commonly been considered powerful tools for gaining information on the hydrological processes for a long time [2,3]. At the catchment scale, they have been extensively used to understand the runoff mechanisms by collecting isotopic data from the different components participating in the hydrological cycle. The number of components to be monitored depends on the type of catchment. For example, in glacierized catchments, the isotopic data of the snow and glacier melt water, in addition to those of precipitation and runoff, are essential [4,5]. In catchments located at lower altitudes, it is common to collect water isotopes from precipitation, runoff, soil moisture and groundwater [6,7,8].
The above-mentioned studies have mainly involved small-catchments (in which areas are generally less than 100 km2) as the complexity and heterogeneity of hydrological processes are often less there. In larger catchments, the use of water isotopes may be less effective as the number of components with different isotopic compositions generally become higher and their signal unclear. For example, the runoff in large catchments can be fed by aquifers with different isotopic signals while rainfall may not have a unique isotopic value if it takes place at different altitudes (i.e., altitude effect).
Moreover, when considering the annual seasonality of isotope data, isotopic variations in river discharges are generally damped due to the times spent by water molecules travelling along the catchments (i.e., transit times) and to the contribution of groundwater, which has a relatively constant isotopic signal. As a result, the longer the paths are (and/or the larger is the contribution of groundwater to rivers) the more the variation of the isotopes in rivers is damped. In the case of large catchments characterized by long transit times or catchments in which the quota of groundwater in discharge is prevalent, the isotopic variation of isotopes in river water may even be nil [9].
In addition, several processes can take place on a catchment scale modifying the former isotopic signal of precipitation collected in rain gauges. These processes are usually enhanced in large catchments. Among others and without claiming to be exhaustive, evapotranspiration occurring within the soils during specific periods of the year or anthropogenic impacts due to reservoirs, irrigation and water treatment plants are responsible for river water that is typically enriched in heavy isotopes compared to the precipitation. The same effect (i.e., river water enriched in heavy isotopes compared to precipitation) is related to sublimation processes acting on the snowpack accumulating in the uppermost parts of the catchments [10,11].
Due to the problems outlined above, few studies in the literature have used water isotope time-series to depict hydrological processes from large catchments worldwide (see, for instance and without claiming to be exhaustive: [12,13,14,15,16,17,18,19]). Consequently, even fewer studies have compared the temporal and spatial behaviour of water isotopes from rivers at a regional scale. A recent paper by [20] has analysed the temporal and spatial patterns of water isotopes in nine large catchments from Germany and their relationships with isotopic composition in precipitation and some selected catchment characteristics. Firstly, the Authors confirmed that the use of the time series of water isotopes provides information on the hydrological processes in large catchments. Secondly, water isotopes in the river can deviate from the isotopic compositions of precipitations collected in the large catchments because of natural and anthropogenic processes such as evapotranspiration from the soils and river network as well as evaporation from reservoirs. Furthermore, they verified statistical associations between the average isotopic values recorded in catchments and some specific physiographic characteristics.
In this study, we have combined several existing datasets concerning both rains and rivers from the northern Italian Apennines, an area in which the runoff response to precipitation events is usually rapid in all the rivers due to the widely outcrop of clay-rich bedrocks. A final dataset was created and includes monthly isotopic data collected in 4 rain gauges (period January 2003–December 2006) and 9 rivers (period January 2006–December 2007). The data are integrated with the corresponding monthly discharges in river gauges where water isotopes were collected.
The purpose of this work is to verify which information can be gained by this short-time (2 years) isotopic dataset consisting of grabbed water samples from large catchments made up of clay-rich bedrocks. Despite the use of short-time (<3 years) and scattered isotopic series in rivers being uncommon, their temporal fluctuations should be enhanced by the poorly permeable bedrocks so that some hydrological information can be preserved. In this the paper, we will provide an overview of the spatial and temporal variability of isotopic composition in precipitations and river waters together with estimates of residence times in catchments. Finally, time series analysis of the isotopes in rivers is tested and possible statistical associations with selected catchment characteristics will be investigated.

2. Study Area

The study area extends over 6261 km2 in the northern Italian Apennines and includes 9 catchments between the Trebbia River and the Savio River (Figure 1). Elevation decreases toward NE direction and ranges from the 2165 m a.s.l of Mt. Cimone to approximately 40 m a.s.l. of the Savio River gauge. As reported in [21], the mean annual rainfall distribution over the period 1990–2015 exceeds 2200 mm/y near the main watershed divide and progressively decreases to about 900 mm/y in foothills. The rainfall distribution during the year is characterized by a marked minimum in the summer season and two maxima during autumn (the main one) and spring. At the end of the winter season and in the vicinity of the main watershed divide, the cumulative annual snow cover can reach 2–3 m. Potential evapotranspiration ranges from about 500 mm/y up to 650 mm/y in the lowlands and is mainly active during the summer months. All the 9 rivers included in this study originate from the main watershed divide and flow toward NE. Six rivers (namely, Trebbia, Nure, Taro, Enza, Secchia and Panaro) are tributaries of the Po River while the other 3 (Reno, Lamone and Savio) enter the Adriatic sea after the sampling sites. Catchment areas are between 193 km2 (Lamone) and 1300 km2 (Secchia), while flow lengths range from 28 km (Enza) to 85.2 km (Secchia). Mean annual discharges over the period 2006–2016 are included between 8.4 m3 s−1 (Savio) and 30.4 m3 s−1 (Secchia).
Recently, [22] grouped all the bedrocks outcropping in the northern Italian Apennines into six main hydrogeological complexes (namely: clay, marl, flysch, foreland flysch, ophiolite, limestone). Among them, the most represented (in terms of areal coverage) are those composed of clay-rich materials (clay, marl and flysch hydrogeological complexes; see Figure 1). These are considered as impermeable or poorly permeable bedrocks leading to a runoff response of rivers that closely follows the rainfall distribution during the year (pluvial discharge regime). Rivers originating from the most elevated parts of the main watershed divide (Secchia, Panaro) are characterized by a nival-pluvial discharges regime as they are influenced by the melting of snow cover accumulated during the winter months in the upper parts of their catchments.
Due to this, marked low-flows take place in the summer-beginning of autumn periods (August, September and October) while floods usually occur in autumn (October and November) and spring (March and April). With reference to intense precipitation events occurring in the wettest periods, the time-lag of peak discharges is usually less than 24 h.

3. Dataset

3.1. Isotopic Data

Stable oxygen and hydrogen isotope data from precipitation located in the study area were collected from an unpublished master’s thesis [23] and from selected published papers [24,25]. The result of these activities is a final dataset consisting of monthly isotopic data from 4 rain gauges (maximum time window: January 2003–December 2006) located at different altitudes along the northern Apennines (Figure 1). Monthly data of the oxygen and hydrogen isotopes from 9 rivers in the study area are derived from [26]. In the case of both precipitation and river waters, isotopic analyses were carried out by using Isotope Ratio Mass Spectrometry (IRMS). The results are reported as a deviation of the sample from the standard (Vienna Standard Mean Oceanic Water: V-SMOW) and presented in the δ-notation as per mil (‰), where δ = [(RS/RV-SMOW) − 1] × 1000; RS represents either the 18O/16O or the 2H/1H ratio of the sample and RV-SMOW is 18O/16O or the 2H/1H ratio of the V-SMOW. In all the cases, instrument precision (1σ) was in the order of ±0.05‰ for δ18O and 0.7‰ for δ2H.
Starting from δ18O and δ2H values, the deuterium-excess (dE; [27]) was calculated for each sample as:
dE (‰) = δ2H − 8 δ18O
Corresponding uncertainty on dE were assessed as ±0.7‰. All the data are reported in the form of Supplementary material.

3.2. River Discharge and Catchment Data

River discharge data from the 9 gauges in which [26] sampled water for isotopic analyses were provided by [28]. They consist of daily discharge data from 1 January 2006 to 31 December 2007. It is worth nothing that all the 9 catchments are characterized by a low degree of urbanization and no presence of remarkable artificial reservoirs and/or diversions.
Following [20], 7 catchment characteristics (or descriptors) related to catchment area (A), elevation (H), precipitation (P), flow length (F) and specific mean annual runoff (q) were used to identify associations with water isotopes and the estimated residence times. As a further catchment characteristic, the specific river runoff exceeded for 95% of the observation period 2006–2007 was taken into account (q95). The latter is a low flow index that is used worldwide for the regionalization procedure and can be estimated even from a relatively short time series of daily runoffs [29,30].
A summary of the 8 descriptors is given in Table 1 and Table 2. In more detail, average annual precipitation (P) was calculated for the corresponding reference period (2006–2007) for 42 stations. The average values for catchments were estimated by spatial interpolation via ordinary kriging [31]. Area (A), elevation (Hmin, Hmax, Hmean) and Flow length (F) were derived from a 5 × 5 m gridded digital terrain model that was created by the digitalization and linear interpolation of contour lines represented in the regional topography map at a scale of 1:5000.

4. Methodology

4.1. Comparison between Isotopic Compositions in Precipitations and Surficial Water

Averages were computed for isotopic time series of precipitation and river discharge. In detail, the statistical analyses of isotope ratios (δ) and dE were carried out by both (1) arithmetic mean MA and (2) weighted mean Mw as follows:
M A = i = 1 N δ i N
M w = i = 1 N δ i W i i = 1 N W i
where Wi and δi are the amount of precipitation (or river discharge) and the measured monthly isotopic compositions in the i-th month of the series, respectively. N is the number of observations.
Furthermore, linear regression models were used to compare the relationship between δ18O and δ2H in water from precipitation and rivers; the linear relationships are referred to as Meteoric Water Lines (MWLs) and River Water Lines (RWLs). The predictive performances of the regression models were tested by means of Pearson’s correlation coefficient coupled with the t-distribution test.

4.2. Estimation of Mean Residence Times (MRTs) in Rivers

Mean Residence Time (MRT) represents the average time (in months) spent by water molecules within a catchment. Among others, MRTs in catchments can be estimated by comparing the δ18O time series in precipitation and river waters [9]. The method (sine-wave fitting) exploits the seasonal changes of δ18O values in precipitations and rivers, which commonly approximate input (precipitation) and output (rivers) sinusoids with different amplitudes. In fact, amplitudes of sinusoids from river waters are usually reduced as compared with those affecting precipitations because of the time spent by the water molecules travelling along their catchments. It is worth noting that the longer the path spent by the water molecules, the more dampened the δ18O output signal is. Therefore, if we consider the same input signal, smoother output sinusoids indicate catchments with higher MRTs.
From a mathematical point of view, input (precipitation) and output (rivers) sinusoids of the observed δ18O time series can be calculated following the procedures reported in [32,33,34]. In greater detail, the input time series of δ18O in precipitation can be approximated with a sinusoid δ18Oin(i) representing the predicted isotopic composition of precipitation at time i (in months):
δ 18 O i n ( i ) = δ 18 O m e a n + A i n cos ( ω i )
in which Ain is the amplitude, ω is the angular frequency (which is equal to 2π/T, where T is the period between two consecutive δ18O peaks, here sets are equal to 12 months), δ18Omean is the estimated mean δ18O of precipitation over the whole period. The output time series δ18Oout(i) in river has a smaller amplitude Aout and a phase lag Φ with respect to the input sinusoid δ18Oin(i):
δ 18 O o u t ( i ) = δ 18 O m e a n + A o u t cos ( ω i +   Φ )
If we process both Equations (3) and (4), MRTs (in months) for each river can finally be obtained as:
MRT = [ A i n A o u t ] 2 1 ω
It should be noted that the method requires the a priori weighting procedure of the former δ18O time series in precipitation and rivers as their isotopic compositions are linked to the precipitation and discharge patterns [9]. Here, we followed the procedure reported in [35], where the input δ18O time series (i.e., δ18Oin(i)) is obtained by weighting the δ18Oi of monthly precipitation with its volume (P, in mm):
δ 18 O i n ( i ) = N P i i = 1 N P i ( δ 18 O i δ 18 O G ) +   δ 18 O G
where N is the number of months for which isotopes have been collected (in our case from January 2003 to December 2006, i.e., 48), while δ18OG is the long-term mean input isotope ratio, which is calculated as:
δ 18 O G = i = 1 N P i δ 18 O i i = 1 N P i
This approach represents an attempt to approximate the real mass precipitation flux contributing to the river discharge. As in the case of the input δ18O time series, also the output δ18Oout time series must consider the river discharge pattern. The output δ18Oout(i) time series for each river was reconstructed by using the above-mentioned Equations (6) and (7), in which Qi (discharge at the i-month) is used instead of Pi.
Goodness of fits (coefficient of determination R2 and p-value) are given for both the input (δ18Oin) and output (δ18Oout) functions. Uncertainties of the MRTs are also reported and were calculated with the error propagation method (see [36]). As pointed out by [1], the above-mentioned MRT estimates may be affected by remarkable aggregation bias. The latter is mainly due to a certain degree of heterogeneity within the catchment, which leads relationships between tracer cycle amplitudes to be largely nonlinear. Moreover, [1] also highlighted a problem of non-stationarity as the ratio between the volume of water stored within the catchment and the average water flux is unlikely to be constant in time (i.e., steady-state condition is often not reached). The bias increases with the size of the catchment (i.e., high spatial heterogeneity and long residence times) and the simultaneous presence of sub-catchments with different travel times, which in turn leads MRTs to be largely underestimated. However, [1] confirmed that Aout/Ain ratios are themselves “aggregated with very little bias and also very small random errors.” In order to further check the possible presence of aggregation bias within the 9 catchments, Aout/Ain ratios were compared with the corresponding MRTs. In case the size of the catchments (here we recall that catchments areas vary of 1 order of magnitude being included between 208 km2 and 1303 km2, see Table 2) do not lead to aggregation bias, Aout/Ain–MRT pairs should line up according to a straight line.

4.3. Reference Isotopic Compositions in Precipitations from the Northern Apennines of Italy

Precipitation waters worldwide show δ18O and δ2H values that are aligned along a regression line (i.e., the “Global Meteorological Water Line”: GMWL; [37]), which is defined as:
δ2H(‰) = 8.0 δ18O + 10.0
Although this relation is valid everywhere, more accurate regression lines can be obtained by selecting δ18O and δ2H values from reduced areas. This is due to the specific isotopic fractionation processes (i.e., vapour pressure and temperature conditions) controlling the precipitation over each area. For instance, by taking into account the precipitation waters alone from the Mediterranean area, [38] obtained the Mediterranean Meteoric Water Line (MMWL):
δ2H(‰) = 8.0 δ18O + 22.0
in which the intercept is slightly higher than that of the GMWL. On the contrary, air vapour originating from Central Europe area is recognized to provide precipitation with lower intercepts than the above-mentioned ones (see the meteoric water line, here abbreviated CEMWL for convenience, recently obtained by [39] for the Balkan zone: δ2H(‰) = 7.4 δ18O + 2.7).
This fact allows us to roughly discriminate precipitation waters originating from the Atlantic ocean (dE close to +10.0) from those of the Mediterranean Sea (mainly the Tyrrhenian sea; higher values of dE) and Central European areas (lower values of dE).
Recently, several authors have dealt with the isotopic composition of precipitations from several areas in the northern Apennines of Italy [40,41,42,43,44,45]. They obtained a number of Local Meteoric Water Lines (LMWLs) in which the slopes are close to that of the GMWL, MMWL and CEMWL (between 7.7 and 9.0) while intercepts are characterized by higher variability (from 6.8 to 21.5).
The high variation affecting the “intercept” parameter confirms the contributions of precipitation coming from either northern/central Europe and the Atlantic Ocean and the Tyrrhenian Sea [43]. It should be pointed out that these three main areas of origin of the vapour masses have been further confirmed by means of the back trajectories calculation with different atmospheric tracers [46].
Isotopic gradients (i.e., the variation of δ18O or δ2H every 100 m of altitude) in the investigated area are not univocal as they vary more significantly than other mountainous areas in Italy. In particular and with reference to δ18O-altitude most of the δ18O gradients are between −0.15 to −0.25‰/100 m [43,47] but at a local scale may be much more negative (see, for instance, [44], 2017: −0.45‰/100 m) or even nil. For the northern Italian Apennines, [43] explained this behaviour as a consequence of the location of the rain collector (for which the gradient is obtained) with respect to the prevailing wind direction of the air masses. In fact, mountain reliefs play a “shadow effect” in which the relatively dry air masses on the leeside of mountain ranges have water that is enriched in 16O with respect to wetter air masses on the windward side [48,49,50]. The higher the range, the more pronounced the isotopic rain shadow, and thus the area located in the vicinity of the main watershed divide may be more affected by this phenomenon. It should be highlighted that leeward and windward sides change depending on the origins of the air masses. In the case of precipitation related to air masses originating from the Tyrrhenian Sea, the sector of the northern Apennines located north eastward of the main watershed divide, that is, our study area, is defined as the leeward side. On the contrary, water vapour from northern and central Europe condenses in the study area that is now considered the windward side.
In addition, [42] reported an anomalous depletion of δ18O in the foothills of the northern Italian Apennines for the period 2002–2004, which led to mean annual isotopic compositions that were considerably more negative than the normal ones characterizing the lower altitudes. This was due to the marked decrease in the spring and summer precipitation so that the mean annual isotopic values were slightly depleted in 18O. As a result of the δ18O—depletion in the lowland areas such as that which occurred in 2002–2004, the isotopic gradients may be nil or even positive. [25] pointed out that the same effect (i.e., change of mean annual δ-values in rain collectors) can be related to changes in the relative proportion of air masses of different origin (Tyrrhenian, Atlantic and Central European origin) that leads the final isotopic composition to be modified.

4.4. Hierarchic Cluster Analysis

The hierarchical cluster analysis was carried out to identify similarities among the isotopic time series from the 9 different catchments. Clustering was done according to unweighted pair-group average (or centroid) method, in which each group consisted of a δ18O time series from a river. The method was based on a step-by-step procedure in which, by the end, δ18O time series were grouped into branched clusters (dendrogram) based on their similarities to one another. Among the several δ18O time series from rivers, the two most similar ones were selected and linked based on the smallest average distance between all δ18O values. Progressively more dissimilar time series were linked at greater distances; at the end, they all were joined to one single cluster. The cophenetic coefficient was used as a measure of similarity between each pair of clusters; being analysed more than 2 time series, the dendrogram was supported by a cophenetic distance matrix. Further details on this methodcan be found in [51].

4.5. Multivariate Analysis of Variance (Manova)

In order to identify similarities between river sample sites and their isotopic time series, a One-way Multivariate Analysis of Variance (MANOVA) was carried out. This is the multivariate version of the univariate ANOVA testing whether two (or more) groups have the same multivariate mean. The analysis is carried out by pairwise comparisons, in which the within-group covariance matrix pooled over all groups participating in the MANOVA. Tests are significant with p < 0.01. Further details can be found in [51].

5. Results

Table 3 summarizes the mean (MA) and weighted (Mw) values of δ18O, δ2H and dE in water from the rain collector. δ18O-MA and δ2H-MA ranges from −8.96‰ (Parma) to −9.67‰ (Langhirano) and from −61.1‰ (Parma) to −67.5‰ (Langhirano). dE is between 9.8‰ (Langhirano) and 10.6‰ (Parma). Distributions of δ18O and dE values from rain collectors were skewed, with mean values that may differ remarkably from the median (see box plots in Figure 2). Weighted averages (δ18O-MW and δ2H-MW) were slightly negative than the aforementioned unweighted means being between −9.48‰ (Parma) and −9.93‰ (Langhirano) and −64.3‰ (Parma) and −69.1‰ (Langhirano), respectively. In all the rain collectors, dE-MW was higher than the corresponding MA values as they ranged from 10.3‰ (Langhirano, Berceto) to 11.6‰ (Parma). The aforementioned variation between MA and MW pair values is due to the fact that precipitation amount as well as its isotopic composition is not uniformly distributed during the year. In detail, and with reference to the Parma rain collector (here considered as representative of the behaviour of others rain collectors, as clearly visible in Figure 2a), δ18O values showed seasonal variations with a minimum in winter (February 2005: −17.64‰) and a maximum in summer (July 2003: −3.02‰).
The same behaviour is noticed for δ2H (not shown in Figure 2 being characterized by the same pattern of the δ18O), which in case of the Parma rain collector showed a minimum of −127.3‰ (February 2005) and maximum of −13.41‰ (July 2003). When δ18O-MW and δ2H-MW were more negative than the corresponding δ18O-MA and δ2H-MA, this meant the specific year was characterized by a larger amount of winter precipitation (i.e., more depleted isotopic values) and/or a lower quantity of summer precipitations. As dE is obtained from the data of δ18O and δ2H, its composition in precipitation varies during the year and is usually higher in the period between the end of the autumn and spring seasons (Figure 2b; November 2004: 20.2‰; January2006: 19.3‰) and lower in summer (June 2006: −2.3‰). As a result, a larger amount of autumn-spring precipitations leads to a higher value of dE.
By considering the river waters (Table 4 and Figure 2c,d), unweighted means are more positive than the corresponding values from rains. In fact, the δ18O-MA are included between –7.46‰ (Taro) and −9.01‰ (Secchia) while δ2H-MA are between −48.2‰ (Taro) and −61.4‰ (Secchia), respectively. dE-MW values are slightly higher than in the rain as they range between 10.5‰ (Secchia, Panaro) and 12.0‰ (Enza). By considering the weighted averages, all δ18O-MW values became more negative as they are between −7–60‰(Taro) and −9.09‰(Secchia). δ2H-MW values were not all more negative than δ2H-MW and are included between −48.5‰ (Trebbia) and −61.1‰(Secchia). The dE-MW values are slightly more positive than the corresponding dE-MA being between −10.4‰ (Panaro) and 12.7‰ (Secchia). As in the case of precipitation, the mismatch between unweighted (MA) and weighted averages (MW) pair values reflects an isotopic composition that was not uniformly distributed during the year (Figure 2c). In fact, water in rivers was much more enriched in δ18O at the end of the summer season (June, July, August) when the flow rates were lower.
Samples from rivers showed a strong attenuation of the δ18O signatures in comparison with precipitation data (Figure 3). In detail, the Lamone and Savio rivers were characterized by a higher dispersion of the δ18O values than the other rivers. On the contrary, dE distributions from river samples were poorly attenuated with box limits (i.e., the 25th and 75th percentiles) that are often in range with those from rain water (Figure 3). This was evident for the Secchia, Panaro, Savio and Lamone rivers.
δ18O and δ2H were correlated, both in the case of precipitation and river water. The δ18O–δ2H relationships are summarized in Table 3 (Meteoric Water Lines MWLs from rain gauges) and Table 4 (River Water Lines RWLs from rivers), in which slopes, intercepts and coefficients of determinations (R2) are reported separately. Slopes obtained interpolating δ18O–δ2H pairs from rain collectors (MWLs) ranged between 6.9 (Berceto) and 7.8 (Lodesana, Langhirano). With the exception of the Berceto rain collector (−1.9), all intercepts were positive, being between 6.9 (Parma) and 7.8 (Lodesana). Slopes from RWLs were lower and ranged from 4.3 (Secchia) to 6.2 (Taro). Intercepts were always negative showing values from −2.0 (Taro) and −22.8 (Secchia). In the case of the Secchia and Panaro rivers, R2 values (0.33 and 0.36, respectively) indicate weaker correlations that are due to 4 monthly samples from January 2017 to April 2017 lying far below the alignment of the others (not shown).
As described in Section 4.2, residence times of the river were estimated starting from a representative input function (Figure 4). A sinusoid was fitted to the δ18O data of precipitation and the corresponding amplitude (Ain = 4.250) together with the goodness of fits are summarised in Table 5. Fit was statistically robust (p < 0.0001) and the coefficient of determination was in the order of those obtained in other studies (see, for instance: [52,53,54]).
Modelled output sinusoids for rivers are also reported in Table 5. With the exception of the Savio River (p value remarkably higher than 0.0001 and R2 = 0.30), the fits were characterized by coefficients of determination (R2) that were higher than 0.48 and always, or close to being, statistically significant (p < 0.0001). Amplitudes of output sinusoids (Aout) differed slightly as they were between 0.397 (Secchia) and 0.974 (Lamone).
MRTs were estimated with Equation (5) and are also reported in Table 5; they ranged between 8 (Lamone) and 21 (Secchia, Trebbia) months. As shown in Figure 5, the relationship between the Aout/Ain ratios and the MRTs for the 9 catchments deviated little from a linear function, indicating a reduced effect of aggregation bias on the MRT estimates due to a large quota of river waters with fast transit time (lower than 1 year). Thus, the consequent underestimation effects did not significantly affect the final MRT values.
Linear relationships between catchment characteristics and the weighted δ18O and dE averages were not significant as p values were always greater than 0.01 (plots are not reported). On the contrary, some catchment characteristics were correlated with the MRTs of rivers (Figure 6). This was the case for the mean annual specific runoff (q; R2 = 0.59 and p < 0.01), the maximum altitude (Hmax; R2 = 0.63 and p < 0.01) and the specific low flow discharge q95 (R2 = 0.44 and p < 0.01).
In Table 6, correlation coefficients between pairs of δ18O time series from the 9 rivers are reported in the form of a correlation matrix (Manova analysis). With the exception of the south-easternmost river of the area (Savio River), the monthly isotopic composition of rivers strongly correlated with each other (R2 > 0.67). However, some pairs of rivers were characterized by higher correlation coefficients: this is the case of Trebbia-Nure and Secchia-Panaro rivers, with R2 equal to 0.95 and 0.93, respectively. The cluster analysis (Figure 7a) among the δ18O time series from rivers demonstrated that the Savio river was associated with none of the other rivers while two main groups of catchments were clearly separated: the first one comprises Enza, Nure, Secchia and Panaro while the other group was Trebbia, Taro, Reno and Lamone. Different results were obtained from the analyses of the dE time series: correlation coefficients were always lower than in case of δ18O time series and more rivers were not associated with each other (see Table 7 and Figure 7b). In detail, the dE time series of Taro, Trebbia and Nure were associated, as well as those of Secchia and Panaro. This was also the case of the dE time series from Lamone and Reno rivers.
River monthly discharges affected the isotopic composition of river water. In detail, with the exception of the Taro river, values of δ18O depleted with increasing flow rates (Figure 8). On the contrary, the dE values in the river water increased with discharge in almost all rivers (Figure 9). In the case of δ18O-discharge relationships, statistical analyses were significant for 4 rivers (Enza, Panaro, Reno, Lamone). In some cases (namely: Lamone and Savio), the relationship looks exponential rather than linear. dE -discharge linear regressions were statistically significant for 5 rivers (Trebbia, Nure, Enza, Lamone). It must be specified that, as in the case of δ18O-discharge relationships, dE-discharge regressions also appeared to be exponential (Reno, Lamone, Savio).

6. Discussion

6.1. Isotopic Comparison between Precipitation and Rivers and Mean Residence Time of Surface Water

For the considered time period (January 2003–December 2006), δ-values from rain collectors were not characterized by altitude-effect, that is, no depletion of δ18O (and δ2H) with altitude were noticed. In particular, while the rain collectors located at lower altitudes (Parma, Lodensana and Langhirano) showed the depletion effect (which leads to a gradient of about −0.22‰/100 m), those located at the highest altitude and in the vicinity of the main watershed divide (Berceto; 800 m a.s.l.) were characterized by enriched values of δ18O. The latter isotopic value should instead characterize lowland areas in normal periods.
As already anticipated in Section 4.3, [43] have confirmed an overall isotopic gradient (of δ18O) between −0.15 to −0.25‰/100 m that may be affected by changes in space and time. As explained by both [43] and [25], this behaviour depends on (i) local topographical effects (leeward/windward with respect to the origin of the air masses) (ii) modification in the precipitation pattern during the year and (iii) changes in the relative proportion of precipitations by air masses of different origin. In the considered time period (January 2003–December 2006), the most negative δ–values were related to a reduction in the mean annual precipitation that occurred over the lowland areas (from value of 900 mm to about 500 mm). This reduction partially took place during the spring and summer months, that is, when precipitations are usually characterized by enriched values of δ18O and δ2H. [42] have confirmed this reduction in the precipitations in lowlands areas at the foothills of the northern Apennines for the period 2004–2006. Nonetheless, [25] highlighted that the proportion of air masses of different origin had changed in the considered time period and precipitations in the vicinity of the main watershed divide were mainly controlled by air masses of Tyrrhenian origin that exhausted their condensation there. This fact is supported by the sampling of a low-yield spring located near the main watershed divide in which δ18O values were in the order of −8.50‰ [25].
In the study area, mean annual values of δ18O (and δ2H as well) in water samples from rivers were less negative than those from rain collectors (Table 3 and Table 4 and Figure 3). This is true, even taking into account the weighted δ-values with monthly amounts of precipitation (in rain collectors) and flow rates (in rivers).
δ18O enrichment noticed in rivers may indicate that these waters have been subjected to evaporative and/or evapotranspirative processes that have led to 16O enrichment. δ18O–δ2H relationships are reported in Table 3 (Meteoric Water Lines MWLs from rain gauges) and Table 4 (River Water Lines RWLs from rivers) together with slopes, intercepts and coefficients of determinations (R2). Being between 4.3 to 6.2, slopes of the RWLs were always lower than those of the precipitations from the northern Apennines (slopes from 6.9 to 7.8). In addition, in this study all intercepts in RMLs were negatives compared to those of MWLs. These values were in the order of those reported by [25] and [55] for other rivers from the northern Apennines and confirmed that these waters have undergone evaporation.
Some authors highlighted that the lowering in slopes from δ18O–δ2H relationships in rivers is due to instream evaporation. Such processes (i.e., instream evaporation) inducing slope lowering in δ18O–δ2H relationships have already been detected in rivers from non-arid environments such as Germany [20], USA [56], India [57] and Australia [58]. Moreover, they were already reported in some rivers from Italy (Arno River: [14]; Reno River: [59]; Po river: [60]).
If so, and as also suggested by [20] and by [57], a statistical association between river flow-lengths (or catchment area) and δ18O should have been found. In detail, such an association should be represented by an increase in δ18O values with flow-lengths or with catchments areas because of the progressive instream evaporation along the catchment.
In the study area, no relationships were found between unweighted (or weighted δ18O values) with flow-lengths (and catchment areas) of rivers. This leads us to suppose that, even if active, instream evaporative/evapotranspirative processes are not mainly responsible for the δ18O enrichment. Moreover, the comparison of the δ18O weighted and unweighted values with other catchment characteristics allows us to state that the process of δ18O enrichment along catchments is not linked to these factors.
Although [24,34,44,45] have recently reported for the area a number of low-yield springs in which water isotopes were not modified by evapotranspiration processes, slopes between 4 and 6 were found in several freshwater springs [61,62,63]. In this case, there is evidence of a pre-infiltrative evaporation/evapotranspiration that acted by modifying the waters before their infiltration towards the aquifer and, subsequently, to the base flow of rivers.
This is in agreement with the increase in dE detected in river water as the contribution of the base flow from springs in the upper part of the catchments (see [25]). In fact, as in the case of δ-values, there were differences in the weighted dE obtained from the water of the rain collectors and those of the rivers. In detail and with the exception of the Panaro river, values are always slightly higher in rivers than in rain collectors.
MRTs are between 8 and 21 months indicating hydrological differences among catchments that are not explained by aggregation biases (i.e., underestimation effects in MRT estimates are reduced). In order to identify such differences, linear relationships between selected catchment characteristics (namely: catchment area, altitude of stream gauge, maximum altitude, mean altitude, precipitation, flow length, specific mean annual runoff and specific runoff exceeded for 95% of the time) and the weighted δ18O and dE averages were carried out. All linear relationships were not significant as p values were always greater than 0.01. On the contrary, 3 catchment characteristics were positively correlated with the MRTs of rivers (Figure 6). This is the case of the mean annual specific runoff q (R2 = 0.59 and p < 0.01), the maximum altitude (Hmax; R2 = 0.64 and p < 0.01) and the specific low flow discharge q95 (R2 = 0.44 and p < 0.01).
As pointed out by [22], in the northern Italian Apennines, the specific low flow discharge q95 is related to the permeability of the geological formations outcropping in the catchment. The higher the permeability of the geological formations outcropping in the basin, the greater the q95 that represents the quota of groundwater released to the river network through springs. It is evident that an increase in q95 indicates a greater quantity of groundwater that makes up the flow of the watercourse; these waters are characterized by longer paths leading to an increase in the MRT. Therefore, the correlation between q95 and MRTS of water, albeit weak, could be due to the groundwater signal in the flow of rivers. The correlation between MRTs and maximum altitude of the catchments is significantly higher (R2 = 0.64); this may be due to a number of concurring causes related to the Hmax descriptor. Firstly, and as already evidenced by [22], the catchment developing from the main watershed divide (i.e., the areas with higher mountain peaks) is characterized by a large number of springs. Secondly, as anticipated in Section 2, the higher areas are characterized by the presence of a snow cover that blocks the water molecules in solid form for a few months until they melt, resulting in an increase of the MRTs of water in the catchment.

6.2. Longitudinal Control on Monthly δ18O and dE Time Series from Rivers

While the value of δ18O is linked to the condensation temperature of the drops in the clouds, dE is linked to the temperature of evaporation of the moistures composing the air masses. It is evident that δ18O and dE time series from rivers may indicate important evidence of similarity in the mechanisms of condensations (orographic precipitations etc.) and in the main source area of the air masses over the catchments. By following [64], we firstly plotted the δ18O–δ2H values from rivers together with the three main meteoric lines GMWL, MMWL and CEMWL reported in Section 4.3 (Figure 10). All the isotopic pairs were included with MMWL and CEMWL, while most of the samples lay above the GMWL. This means that all water were of Tyrrhenian, Central European and Atlantic origin. With more detail, some samples lay above CEMWL while others in the vicinity of the MMWL indicating that some waters were more probably related to air masses from Central Europe and the Tyrrhenian Sea, respectively. Among the others, samples close to the CEMWL came mainly from the central and easternmost part of the study area (namely: Secchia, Panaro, Lamone, Savio; Figure 10).
With reference to δ18O time-series, all rivers were clustered to one main group with the only exception of the Savio river (Figure 7a). Two sub-groups made up the main group: the first one included Lamone, Trebbia, Taro and Reno rivers while the second one contained Secchia, Panaro, Nure and Enza rivers. It should be pointed out that some neighbour catchments are more similar than distant ones. For instance, this is the case of Taro-Trebbia, Nure-Enza and Secchia-Panaro. This is also confirmed by looking at correlation values among the different δ18O time-series (correlation matrix in Table 6), from which it again emerges that all the rivers (with the exception of Savio) are strongly associated.
The same analyses carried out with dE indicated that all rivers were clustered to one main group with the exception of the Savio river (Figure 7b). Three sub-groups namely Trebbia-Taro-Nure (western sector of the study area), Secchia-Panaro (central sector of the study area) and Reno-Lamone (eastern sector of the study area), were evident and clearly separated. These three sub-groups are evidenced in the correlation matrix (Table 7).
These results lead us to suppose a certain control made by the orography (and thus by the main watershed divide) over the δ18O and dE time-series. In fact, and as already reported in Section 3.2, the main watershed divide is oriented NW-SE (Figure 1) and is close to the Tyrrhenian sea in the westernmost portion of the study area. This may allow air masses originating from the Tyrrhenian Sea to reach the upper parts of the catchments (Trebbia, Taro, Nure, Secchia and Panaro) in the study area (leeward site), having already suffered a significant enrichment of 16O and 1H in the windward site. Moving eastward, this phenomenon is remarkably reduced as a quota of precipitation from the Tyrrhenian Sea is progressively in favour of air masses from Central Europe. For the considered time-period (2006–2007), the control made by air masses of Tyrrhenian origin over precipitation in the vicinity of the main watershed divide is also evidenced by the isotopic composition of water from the Berceto rain gauge.
Furthermore, we cannot neglect the role of sublimation processes acting on snow cover during the winter and spring seasons. In fact, sublimation occurring during sunny days can modify the former isotopic composition of the superficial snow layers allowing the release of a vapour phase from the solid skeleton to the atmosphere. In this case, the final snow cover does not preserve the isotopic composition of the original snowfall from which it was derived [10,11]. As anticipated in Section 2, Secchia and Panaro rivers (i.e., the central sector of the study area with the highest reliefs) are characterized by nival-pluvial discharges due to the melting of the snow cover in the upper part of the catchments during the spring months. In the δ18O and dE time series there are evidences of sublimation in water samples from January 2017 to April 2017 that lead to weaker correlation in the corresponding RMWLs (as already seen in Table 4).

7. Conclusions

This study shows that a short-time series of scattered water isotopes can be useful in order to depict hydrological processes occurring in large catchments composed of clay-rich bedrocks. We believe that this is an important point given both the expensive costs of the more recent in-continuous and in-field sampling/analysis techniques for water isotopes and for the large number of short-time series made of scattered isotopes already available from large and clay-rich catchments worldwide. The analyses of our short-term datasets consisting of monthly values of water stable isotopes from 4 rain collectors (4-years) and 9 large catchments (2-years) have revealed further insights into the hydrological behaviour of large catchments from the northern Italian Apennines. Firstly, the role of the orography has been confirmed (and thus the main watershed divide) in selecting precipitations originating from the Tyrrhenian Sea, which were progressively reduced, moving eastward in favour of those from the Central Europe.
Secondly, MRTs in catchments from northern Italian Apennines are positively correlated with three descriptors, namely maximum altitude as well as mean annual specific runoff and (even if weakly) specific low flow discharge q95. These descriptors recall the effect of the snow cover (maximum altitude) and bedrocks (specific low flow discharge q95) on delaying the water molecules along the catchments. These statistical associations are further confirmed by verifying the almost linear relationship between the amplitudes ratios of isotopic signals from precipitation and catchments with the MRTs estimates, which in turn indicate a low degree of aggregation biases due to the large quota of water with a fast travel time.
Although a short-time isotopic dataset allowed us to obtain a certain amount of information for these catchments, further efforts should focus on investigating the spatial and temporal patterns observed in precipitations and, in particular, the role of the main watershed divide in selecting the air masses that originated elsewhere. This is a crucial point that could allow us to better define the input isotopic signal and, therefore, further improve the calculation of surficial water MRTs. Moreover, a large number of rain collectors (uniformly distributed over the northern Apennines) would allow us to verify the presence of local isotopic anomalies and of any areas where depletion of the isotopic composition with the altitude (altitude-effect) exists. In addition, isotopic data from snow cover accumulated in the central sector of study area (at least in the vicinity of the main watershed divide) should be acquired to reduce uncertainties related to sublimation processes during winter and spring months. Isotopic river monitoring should be implemented with the collection of a longer time series (at least 5 years) to reduce errors in the calculation of the MRTs. In view of the fact that discharge in rivers from the northern Apennines might change suddenly after rainfall (especially during the wet seasons), water sampling for the isotopic analyses should be performed at least hourly during the floods in order to take into account the bulk samples rather than the grab ones.

Supplementary Materials

The following are available online at https://www.mdpi.com/2073-4441/11/7/1360/s1.

Author Contributions

Conceptualization, F.C. and G.M.; methodology, F.C.; formal analysis, F.C.; data curation, F.C. and A.D.; writing—original draft preparation, F.C.; writing—review and editing, F.C., A.D., G.M.

Funding

This research received no external funding.

Acknowledgments

The Authors would like to thank the three anonymous reviewers for their thoughtful and detailed comments on an early version of this manuscript.

Conflicts of Interest

The Authors declare no conflict of interest.

References

  1. Kirchner, J.W. Aggregation in environmental systems—Part 1: Seasonal tracer cycles quantify young water fractions, but not mean transit times, in spatially heterogeneous catchments. Hydrol. Earth Syst. Sci. 2016, 20, 279–297. [Google Scholar] [CrossRef]
  2. Clark, I.; Fritz, P. Environmental Isotopes in Hydrogeology; CRC Press Lewis Publishers: Boca Raton, FL, USA, 1997. [Google Scholar]
  3. Kendall, C.; McDonnell, J.J. Isotope Tracers in Catchment Hydrology; Elsevier: Amsterdam, The Netherlands, 1998; p. 839. [Google Scholar]
  4. Penna, D.; Engel, M.; Mao, L.; Dell’Agnese, A.; Bertoldi, G.; Comiti, F. Tracer-based analysis of spatial and temporal variations of water sources in a glacierized catchment. Hydrol. Earth Syst. Sci. 2014, 18, 5271–5288. [Google Scholar] [CrossRef] [Green Version]
  5. Rücker, A.; Boss, S.; Kirchner, J.W.; von Freyberg, J. Monitoring snowpack outflow volumes and their isotopic composition to better understand streamflow generation during rain-on-snow events. Hydrol. Earth Syst. Sci. Discuss. 2019, 23. [Google Scholar] [CrossRef]
  6. Burns, D.A.; McDonnell, J.J.; Hooper, R.P.; Peters, N.E.; Freer, J.E.; Kendall, C.; Beven, K. Quantifying contributions to storm runoff through end-member mixing analysis and hydrologic measurements at the Panola Mountain Research Watershed (Georgia, USA). Hydrol. Process. 2001, 15, 1903–1924. [Google Scholar] [CrossRef]
  7. Cras, A.; Marc, V.; Travi, Y. Hydrological behaviour of sub-Mediterranean alpine headwater streams in a badlands environment. J. Hydrol. 2007, 339, 130–144. [Google Scholar] [CrossRef]
  8. Marc, V.; Didon-Lescot, J.F.; Michael, C. Investigation of the hydrological processes using chemical and isotopic tracers in a small Mediterranean forested catchment during autumn recharge. J. Hydrol. 2001, 247, 215–229. [Google Scholar] [CrossRef]
  9. Leibundgut, C.; Maloszewski, P.; Külls, C. Tracers in Hydrology; Wiley: New York, NY, USA, 2011. [Google Scholar]
  10. Stichler, W.; Schotterer, U.; Froehlich, K.; Ginot, P.; Kull, C.; Gaeggeler, H.; Pouyaud, B. Influence of sublimation on stable isotope records recovered from high-altitude glaciers in the tropical Andes. J. Geophys. Res. 2001, 106, 22613–22620. [Google Scholar] [CrossRef]
  11. Sokratov, S.A.; Golubev, V.N. Snow isotopic content change by sublimation. J. Glaciol. 2009, 55, 823–828. [Google Scholar] [CrossRef] [Green Version]
  12. Frederickson, G.C.; Criss, R.E. Isotope hydrology and residence times of the unimpounded Meramec River Basin, Missouri. Chem. Geol. 1999, 157, 303–317. [Google Scholar] [CrossRef]
  13. Yi, Y.; Gibson, J.J.; Hélie, J.F.; Dick, T.A. Synoptic and time-series stable isotope surveys of the Mackenzie River from Great Slave Lake to the Arctic Ocean, 2003 to 2006. J. Hydrol. 2010, 383, 223–232. [Google Scholar] [CrossRef]
  14. La Ruffa, G.; Panichi, C. Caratterizzazione Chimico-Isotopica Delle Acque Fluviali: Il Caso del Fiume Arno; Istituti Editoriali e Poligrafici Internazionali: Pisa, Italy, 2000; p. 101. (In Italian) [Google Scholar]
  15. Lu, B.; Sun, T.; Wang, C.; Wang, J.; Dai, S.; Kuang, J. Temporal and Spatial Variations of d18O along the Main Stem of Yangtze River, China. In Monitoring Isotopes in Rivers: Creation of the Global Network of Isotopes in Rivers (GNIR); IAEA-TECDOC-1673; International Atomic Energy Agency: Vienna, Austria, 2012; pp. 211–219. [Google Scholar]
  16. Marchina, C.; Natali, C.; Bianchini, G. The Po River Water Isotopes during the Drought Condition of the Year 2017. Water 2019, 11, 150. [Google Scholar] [CrossRef]
  17. Martinelli, L.A.; Victoria, R.L.; Silveira Lobo Sternberg, L.; Ribeiro, A.; Zacharias Moreira, M. Using stable isotopes to determine sources of evaporated water to the atmosphere in the Amazon basin. J. Hydrol. 1996, 183, 191–204. [Google Scholar] [CrossRef]
  18. Rank, D.; Wyhlidal, S.; Schott, K.; Jung, M.; Heiss, G.; Tudor, M. A 50 Years’ isotope record of the danube river water and its relevance for hydrological, climatological and environmental research. Acta Zool. Bulg. 2014, 66, 109–115. [Google Scholar]
  19. Vitvar, T.; Aggarwal, P.K.; Herczeg, A.L. Global network is launched to monitor isotopes in rivers. Eos Trans. Am. Geophys. Union 2007, 88, 325–326. [Google Scholar] [CrossRef]
  20. Reckerth, A.; Stichler, W.; Schmidt, A.; Stumpp, C. Long-term data set analysis of stable isotopic composition in German rivers. J. Hydrol. 2017, 552, 718–731. [Google Scholar] [CrossRef]
  21. Cervi, F.; Nistor, M.M. High Resolution of Water Availability for Emilia-Romagna Region over 1961–2015. Adv. Meteorol. 2018, 2018, 2489758. [Google Scholar] [CrossRef]
  22. Cervi, F.; Blöschl, G.; Corsini, A.; Borgatti, L.; Montanari, A. Perennial springs provide information to predict low flows in mountain basins. Hydrol. Sci. J. 2017, 62, 2469–2481. [Google Scholar] [CrossRef]
  23. Sole, F. Gli Isotopi Dell’ossigeno e Dell’idrogeno Come Traccianti nelle Acque Meteoriche per L’alimentazione di Acque Superficiali e Freatiche Nell’area Parmense. Master’s Thesis, Università Degli Studi di Parma: Parma, Italy, 2006, unpublished. p. 202. (In Italian). [Google Scholar]
  24. Cervi, F.; Corsini, A.; Doveri, M.; Mussi, M.; Ronchetti, F.; Tazioli, A. Characterizing the recharge of fractured aquifers: A case study in a flysch rock mass of the northern Apennines (Italy). In Engineering Geology for Society and Territory; Springer: Cham, Switzerland, 2015; Volume 3, pp. 563–567. [Google Scholar]
  25. Iacumin, P.; Venturelli, G.; Selmo, E. Isotopic features of rivers and groundwater of the Parma Province (Northern Italy) and their relationships with precipitation. J. Geochem. Explor. 2009, 102, 56–62. [Google Scholar] [CrossRef]
  26. Martinelli, G.; Chahoud, A.; Dadomo, A.; Fava, A. Isotopic features of Emilia-Romagna region (North Italy) groundwaters: Environmental and climatological implications. J. Hydrol. 2014, 519, 1928–1938. [Google Scholar] [CrossRef]
  27. Dansgaard, W. Stable isotopes in precipitation. Tellus 1964, 16, 436–468. [Google Scholar] [CrossRef]
  28. ARPAE-EMR. Regional Agency for Environmental Protection in Emilia-Romagna Region: Annali Idrologici. 2014. Available online: https://www.arpae.it/sim/ (accessed on 31 May 2017).
  29. Laaha, G.; Blöschl, G. Low flow estimates from short stream flow records—A comparison of the methods. J. Hydrol. 2005, 306, 264–286. [Google Scholar] [CrossRef]
  30. Laaha, G.; Blöschl, G. A comparison of low flow regionalization methods-catchment grouping. J. Hydrol. 2006, 323, 193–214. [Google Scholar] [CrossRef]
  31. Isaaks, E.H.; Srivastava, R.M. An Introduction to Applied Geostatistics; Oxford University Press: New York, NY, USA, 1989. [Google Scholar]
  32. Maloszewski, P.; Rauert, W.; Stichler, W.; Herrmann, A. Application of flow models in an alpine catchment area using tritium and deuterium data. J. Hydrol. 1983, 66, 319–333. [Google Scholar] [CrossRef]
  33. Pearce, A.J.; Stewart, M.K.; Sklash, M.G. Storm runoff generation in humid headwater catchments: 1. Where does water come from? Water Resour. Res. 1986, 22, 1263–1272. [Google Scholar] [CrossRef]
  34. Stewart, M.K.; McDonnell, J.J. Modeling base flow soil water residence times from deuterium concentrations. Water Resour. Res. 1991, 27, 2681–2693. [Google Scholar] [CrossRef]
  35. Bergraann, H.; Sackl, B.; Maloszewski, P.; Stichler, W. Hydrological investigations in a small catchment area using isotope data series. In Proceedings of the 5th International Symposium on Underground Water Tracing, Athens, Greece, 22–27 September 1986; Institute of Geology and Mineral Exploration: Athens, Greece, 1986; pp. 255–272. [Google Scholar]
  36. Reddy, M.M.; Schuster, P.; Kendall, C.; Reddy, M.B. Characterization of surface and ground water δ18O seasonal variation and its use for estimating groundwater residence times. Hydrol. Process. 2006, 20, 1753–1772. [Google Scholar] [CrossRef]
  37. Craig, H. Isotopic variations in meteoric waters. Science 1961, 133, 1702–1703. [Google Scholar] [CrossRef]
  38. Gat, J.; Carmi, I. Evolution of the isotopic composition of atmospheric waters in Mediterranean Sea area. J. Geophys. Res. 1970, 75, 3039–3078. [Google Scholar] [CrossRef]
  39. Bădăluță, C.A.; Perșoiu, A.; Ionita, M.; Nagavciuc, V.; Bistricean, P.I. Stable H and O isotope-based investigation of moisture sources and their role in river and groundwater recharge in the NE Carpathian Mountains, East-Central Europe. Isot. Environ. Health Stud. 2019, 55, 161–178. [Google Scholar] [CrossRef]
  40. Cervi, F.; Ronchetti, F.; Doveri, M.; Mussi, M.; Marcaccio, M.; Tazioli, A. The use of stable water isotopes from rain gauges network to define the recharge areas of springs: Problems and possible solutions from case studies from the northern Apennines. Geoing. Ambient. Min. 2016, 149, 19–26. [Google Scholar]
  41. Longinelli, A.; Selmo, E. Isotopic composition of precipitation in Italy: A first overall map. J. Hydrol. 2003, 270, 75–88. [Google Scholar]
  42. Longinelli, A.; Anglesio, E.; Flora, O.; Iacumin, P.; Selmo, E. Isotopic composition of precipitation in Northern Italy: Reverse effect of anomalous climatic events. J. Hydrol. 2006, 329, 471–476. [Google Scholar] [CrossRef]
  43. Giustini, F.; Brilli, M.; Patera, A. Mapping oxygen stable isotopes of precipitation in Italy. J. Hydrol. Reg. Stud. 2016, 8, 162–181. [Google Scholar] [CrossRef] [Green Version]
  44. Deiana, M.; Mussi, M.; Ronchetti, F. Discharge and environmental isotope behaviours of adjacent fractured and porous aquifers. Environ. Earth Sci. 2017, 76, 595. [Google Scholar] [CrossRef]
  45. Deiana, M.; Cervi, F.; Pennisi, M.; Mussi, M.; Bertrand, C.; Tazioli, A.; Corsini, A.; Ronchetti, F. Chemical and isotopic investigations (δ18O, δ2H, 3H, 87Sr/86Sr) to define groundwater processes occurring in a deep-seated landslide in flysch. Hydrogeol. J. 2018, 26, 2669–2691. [Google Scholar] [CrossRef]
  46. Brattich, E.; Hernández-Ceballos, M.A.; Cinelli, G.; Tositti, L. Analysis of 210Pb peak values at Mt. Cimone (1998–2011). Atmos. Environ. 2015, 112, 136–147. [Google Scholar] [CrossRef]
  47. Cervi, F.; Borgatti, L.; Martinelli, G.; Ronchetti, F. Evidence of deep-water inflow in a tectonic window of the northern Apennines (Italy). Environ. Earth Sci. 2014, 72, 2389–2409. [Google Scholar] [CrossRef]
  48. Kerr, T.; Srinivasan, M.S.; Rutherford, J. Stable water isotopes across a transect of the southern Alps, New Zealand. J. Hydrometeorol. 2015, 16, 702–715. [Google Scholar] [CrossRef]
  49. Scholl, M.A.; Giambelluca, T.W.; Gingerich, S.B.; Nullet, M.A.; Loope, L.L. Cloud water in windward and leeward mountain forests: The stable isotope signature of orographic cloud water. Water Resour. Res. 2007, 43. [Google Scholar] [CrossRef]
  50. Moran, T.A.; Marshall, S.J.; Evans, E.C.; Sinclair, K.E. Altitudinal gradients of stable isotopes in lee-slope precipitation in the Canadian Rocky Mountains. Arct. Antarct. Alp. Res. 2007, 39, 455–467. [Google Scholar] [CrossRef]
  51. Davis, J.C. Statistics and Data Analysis in Geology; John Wiley & Sons: New York, NY, USA, 2001. [Google Scholar]
  52. DeWalle, D.R.; Edwards, P.J.; Swistock, B.R.; Aravena, R.; Drimmie, R.J. Seasonal hydrology of three Appalachian forest catchments. Hydrol. Process. 1997, 11, 1895–1906. [Google Scholar] [CrossRef]
  53. McGuire, K.J.; DeWalle, D.R.; Gburek, W.J. Evaluation of mean residence time in subsurface waters using oxygen-18 fluctuations during drought conditions in the mid-Appalachians. J. Hydrol. 2002, 261, 132–149. [Google Scholar] [CrossRef]
  54. Soulsby, C.; Tetzlaff, D.; Rodgers, P.; Dunn, S.; Waldron, S. Runoff processes, stream water residence times and controlling landscape characteristics in a mesoscale catchment: An initial evaluation. J. Hydrol. 2006, 325, 197–221. [Google Scholar] [CrossRef]
  55. Dadomo, A.; Martinelli, G. Aspetti di idrologia isotopica in Emilia-Romagna. Convegni Lincei Accad. Naz. Lincei 2005, 216, 157. [Google Scholar]
  56. Hogan, J.; Phillips, F.; Eastoe, C.; Lacey, H.; Mills, S.; Oelsner, G. Isotopic tracing of hydrological processes and water quality along the Upper Rio Grande, USA. In Monitoring Isotopes in Rivers: Creation of the Global Network of Isotopes in Rivers (GNIR); International Atomic Energy Agency: Vienna, Austria, 2012; p. 258. [Google Scholar]
  57. Kumar, A.; Sanyala, P.; Agrawal, S. Spatial distribution of δ18O values in river water in the Ganga River Basin: Insight into the hydrological processes. J. Hydrol. 2019, 571, 225–234. [Google Scholar] [CrossRef]
  58. Hughes, C.E.; Stone, D.J.M.; Gibson, J.J.; Meredith, K.T.; Sadek, M.A.; Cendon, D.I.; Hankin, S.I.; Hollins, S.E.; Morrison, T.N. Stable water isotope investigation of the Barwon-Darling River system, Australia. In Monitoring Isotopes in Rivers: Creation of the Global Network of Isotopes in Rivers (GNIR); International Atomic Energy Agency: Vienna, Austria, 2012; p. 258. [Google Scholar]
  59. Carlin, F.; Magri, G.; Cervellati, A.; Gonfiantini, R. Use of environmental isotopes to investigate the interconnections between the Reno River and groundwater (Northern Italy). In Isotope Ratios as Pollutant Source and Behaviour Indicators; IAEA-SM-191/6; International Atomic Energy Agency: Wien, Austria, 1975; pp. 179–194. [Google Scholar]
  60. Marchina, C.; Natali, C.; Fazzini, M.; Fusetti, M.; Tassinari, R.; Bianchini, G. Extremely dry and warm conditions in northern Italy during the year 2015: Effects on the Po river water. Rend. Lincei 2017, 28, 281–290. [Google Scholar] [CrossRef]
  61. Bortolami, G.; De Vecchi Pellati, R.; Ricci, B.; Zuppi, G.M. Indagine Idrogeochimica Preliminare delle Sorgenti Alimentate da Circuiti Profondi della zona Collinare Compresa tra il corso del Po e del Tanaro (Piemonte); Atti I Semin. Inform. Prog. Final. Energetica, Sottopr. En. Geot.; C.N.R.: Roma, Italy, 1979. (In Italian) [Google Scholar]
  62. Iacumin, P.; Venturelli, G.; Burroni, B.; Toscani, L.; Selmo, E. The S. Andrea Bagni waters (province of Parma): Origin, mixing with high-salinity waters and inferences on climatic microvariations. Developments in Aquifer Sedimentology and Ground Water Flow Studies in Italy. Mem. Descr. Carta Geol. D’Italia 2007, 76, 219–227. [Google Scholar]
  63. Doveri, M.; Mussi, M. Water isotopes as environmental tracers for conceptual understanding of groundwater flow: An application for fractured aquifer systems in the “Scansano-Magliano in Toscana” area (Southern Tuscany, Italy). Water 2014, 6, 2255–2277. [Google Scholar] [CrossRef]
  64. Vespasiano, G.; Apollaro, C.; De Rosa, R.; Muto, F.; Larosa, S.; Fiebig, J.; Mulch, A.; Marini, L. The Small Spring Method (SSM) for the definition of stable isotope—Elevation relationships in Northern Calabria (Southern Italy). Appl. Geochem. 2015, 63, 333–346. [Google Scholar] [CrossRef]
Figure 1. Sketch map of the Emilia Romagna region together with locations in which precipitation (letters from a to d) and surficial (numbered from 1 to 9) water samplings have been carried by previous studies. Hydrogeological complexes are reported following [22]; GC: clay; GM: marl; GF: flysch; GFF: foreland flysch; GL: limestone; GO: Ophiolite. For further details, see Table 2 (river gauges) and Table 3 (rain collectors).
Figure 1. Sketch map of the Emilia Romagna region together with locations in which precipitation (letters from a to d) and surficial (numbered from 1 to 9) water samplings have been carried by previous studies. Hydrogeological complexes are reported following [22]; GC: clay; GM: marl; GF: flysch; GFF: foreland flysch; GL: limestone; GO: Ophiolite. For further details, see Table 2 (river gauges) and Table 3 (rain collectors).
Water 11 01360 g001
Figure 2. Monthly variation of δ18O and dE in precipitation ((a,b); Parma: black continuous line; Lodesana: red continuous line; Langhirano: black dashed line; Berceto: red dashed line) and rivers ((c,d); values are represented as aggregated means together with 10th and 90th percentiles, reported as dashed black lines).
Figure 2. Monthly variation of δ18O and dE in precipitation ((a,b); Parma: black continuous line; Lodesana: red continuous line; Langhirano: black dashed line; Berceto: red dashed line) and rivers ((c,d); values are represented as aggregated means together with 10th and 90th percentiles, reported as dashed black lines).
Water 11 01360 g002
Figure 3. Box plot for δ18O (a) and dE (b) of all water samples. Letters from a to d and numbers from 1 to 9 recall precipitation and surficial water, respectively. The whiskers represent the 10th and 90th percentiles, the box limits indicate the 25th and 75th percentiles and the line within the box marks is the median. For further details, readers are referred to Table 2 (river gauges) and Table 3 (rain collectors).
Figure 3. Box plot for δ18O (a) and dE (b) of all water samples. Letters from a to d and numbers from 1 to 9 recall precipitation and surficial water, respectively. The whiskers represent the 10th and 90th percentiles, the box limits indicate the 25th and 75th percentiles and the line within the box marks is the median. For further details, readers are referred to Table 2 (river gauges) and Table 3 (rain collectors).
Water 11 01360 g003
Figure 4. Fitted regression models to δ18O data collected in rain collector (a; 41 months from to January 2003 to December 2006) and rivers (19; 24 months from January 2006 to December 2007). Statistics (uncertainties and sinusoidal fittings parameters) concerning each fitted regression are reported in Table 5.
Figure 4. Fitted regression models to δ18O data collected in rain collector (a; 41 months from to January 2003 to December 2006) and rivers (19; 24 months from January 2006 to December 2007). Statistics (uncertainties and sinusoidal fittings parameters) concerning each fitted regression are reported in Table 5.
Water 11 01360 g004
Figure 5. Relationship between Aout/Ain ratio and MRT for water samples from rivers. For further details, see Table 2 (river gauges).
Figure 5. Relationship between Aout/Ain ratio and MRT for water samples from rivers. For further details, see Table 2 (river gauges).
Water 11 01360 g005
Figure 6. Dependence between mean residence time and the 8 catchment characteristics.
Figure 6. Dependence between mean residence time and the 8 catchment characteristics.
Water 11 01360 g006
Figure 7. Dendrogram for δ18O (a) and dE (b) time series in rivers. For further details, see Table 2 (river gauges).
Figure 7. Dendrogram for δ18O (a) and dE (b) time series in rivers. For further details, see Table 2 (river gauges).
Water 11 01360 g007
Figure 8. Dependence between δ18O and discharge separately for each river. For further details, see Table 2 (river gauges).
Figure 8. Dependence between δ18O and discharge separately for each river. For further details, see Table 2 (river gauges).
Water 11 01360 g008
Figure 9. Dependence between dE and discharge separately for each river. For further details, see Table 2 (river gauges).
Figure 9. Dependence between dE and discharge separately for each river. For further details, see Table 2 (river gauges).
Water 11 01360 g009
Figure 10. δ18O–δ2H plot for water samples from rivers (Trebbia: black line; Nure: black plus; Taro: black cross; Enza: red filled dot; Secchia: black filled circle; Panaro: empty circle; Reno: black diamond; Lamone: empty square; Savio: empty triangle) together with meteoric water lines (GMWL: black dashed line; MMWL: blue dashed line; CEMWL: red dashed line).
Figure 10. δ18O–δ2H plot for water samples from rivers (Trebbia: black line; Nure: black plus; Taro: black cross; Enza: red filled dot; Secchia: black filled circle; Panaro: empty circle; Reno: black diamond; Lamone: empty square; Savio: empty triangle) together with meteoric water lines (GMWL: black dashed line; MMWL: blue dashed line; CEMWL: red dashed line).
Water 11 01360 g010
Table 1. Catchment descriptors included in the analysis.
Table 1. Catchment descriptors included in the analysis.
AcronymVariableUnitsMinimumMeanMaximum
ACatchment areakm21936961303
HminAltitude of stream gaugem43171421
HmaxMaximum altitudem115817842165
HmeanMean altitudem526754944
PPrecipitation mm92410901304
FFlow lengthkm20.955.585.2
qSpecific mean annual runoffL s−1 km−22.215.036.3
q95Specific runoff exceeded for 95% of the timeL s−1 km−20.01.01.7
Table 2. Selected catchment characteristics for the 9 rivers considered in this study.
Table 2. Selected catchment characteristics for the 9 rivers considered in this study.
River (Code)Catchment Characteristics
Hmin (m a.s.l.)Catchment Area (km2)Hmax (m a.s.l.)Hmean (m a.s.l.)q (L s1 km2)q95 (L s1 km2)P (mm)Flow Length (km)
Trebbia (1)257655173593819.41.6130465.5
Nure (2)421208175394415.71.3109420.9
Taro (3)135124617997128.41.099071
Enza (4)231430201680222.40.792428
Secchia (5)471303212169436.31.099985.2
Panaro (6)21258421659397.41.7101740.9
Reno (7)601056194563910.81.1106381.8
Lamone (8)13519311585932.20.5125728.5
Savio (9)4358613615269.80.0116777.4
Table 3. Meteoric water lines (RWLs) along with unweighted (MA) and weighted average (MW) compositions of the precipitation at each rain collector.
Table 3. Meteoric water lines (RWLs) along with unweighted (MA) and weighted average (MW) compositions of the precipitation at each rain collector.
Rain Gauge (Altitude in m a.s.l.)CodeMeteoric Water Lines (MWLs)Unweighted Average (MA)Weighted Average (MW)
SlopeInterceptR2δ18Oδ2HdEδ18Oδ2HdE
Parma (65) a7.66.90.98−8.96−61.110.6−9.48−64.311.6
Lodesana (150)b7.87.80.99−9.44−65.510.0−9.68−67.010.4
Langhirano (220)c7.87.60.99−9.67−67.59.8−9.93−69.110.3
Berceto (800)d6.9−1.90.99−9.13−63.010.1−9.77−67.910.3
Table 4. River water lines (RWLs) along with unweighted (MA) and weighted average (MW) compositions of the river water.
Table 4. River water lines (RWLs) along with unweighted (MA) and weighted average (MW) compositions of the river water.
River (Code)River Water Lines (RWLs)Unweighted Average (MA)Weighted Average (MW)
SlopeInterceptR2δ18Oδ2HdEδ18Oδ2HdE
Trebbia (1)5.00−10.80.63−7.55−48.611.8−7.62−48.512.5
Nure (2)6.09−5.30.97−8.57−57.311.2−8.74−58.011.9
Taro (3)6.19−2.00.90−7.46−48.211.5−7.60−49.511.3
Enza (4)5.40−9.90.82−8.41−55.312.0−8.69−56.812.7
Secchia (5)4.28−22.80.33−9.01−61.410.5−9.09−61.111.6
Panaro (6)4.80−17.50.36−8.77−59.710.5−8.91−60.910.4
Reno (7)5.63−6.70.80−7.80−50.711.7−8.15−52.512.7
Lamone (8)5.45−7.80.95−7.72−49.911.9−8.26−52.813.3
Savio (9)5.00−13.730.82−8.10−54.310.6−8.66−57.511.8
Table 5. Summary of the results for residence times. Estimates are reported (in months) along with uncertainties and sinusoidal fittings parameters.
Table 5. Summary of the results for residence times. Estimates are reported (in months) along with uncertainties and sinusoidal fittings parameters.
Sinusoid (Code)R2Ap-ValueResidence Time (Months)Uncertainties (Months)
Precipitation (a)0.844.250<0.0001n.a.n.a.
Trebbia (1)0.520.4020.000421±9
Nure (2)0.570.512<0.000116±10
Taro (3)0.660.705<0.000112±7
Enza (4)0.510.4710.000518±11
Secchia (5)0.480.3970.001621±15
Panaro (6)0.570.4980.000117±9
Reno (7)0.540.5050.000217±10
Lamone (8)0.760.974<0.00018±4
Savio (9)0.300.7940.039210±9
Table 6. Correlation matrix reporting δ18O associations among the time series in rivers. A progressively more intense green colour is associated with a higher correlation coefficient. * not significant at p = 0.01.
Table 6. Correlation matrix reporting δ18O associations among the time series in rivers. A progressively more intense green colour is associated with a higher correlation coefficient. * not significant at p = 0.01.
TrebbiaNureTaroEnzaSecchiaPanaroRenoLamoneSavio
Trebbia-
Nure0.95-
Taro0.770.86-
Enza0.680.780.85-
Secchia0.870.890.860.84-
Panaro0.710.790.870.890.93-
Reno0.800.830.850.820.860.87-
Lamone0.670.710.770.810.750.810.80-
Savio0.45 *0.30 *0.40 *0.30 *0.38 *0.33 *0.40 *0.43 *-
Table 7. Correlation matrix reporting dE associations among the time series in rivers. A progressively more intense green colour is associated with a higher correlation coefficient. * not significant at p = 0.01.
Table 7. Correlation matrix reporting dE associations among the time series in rivers. A progressively more intense green colour is associated with a higher correlation coefficient. * not significant at p = 0.01.
TrebbiaNureTaroEnzaSecchiaPanaroRenoLamoneSavio
Trebbia-
Nure0.74-
Taro0.770.64-
Enza0.590.620.59-
Secchia0.22 *0.31 *0.17 *0.40 *-
Panaro0.28 *0.32 *0.17 *0.44 *0.86-
Reno0.680.640.510.35 *0.01 *0.05 *-
Lamone0.690.730.690.570.04 *0.14 *0.73-
Savio0.580.510.41*0.36 *0.23 *0.31 *0.38 *0.56-

Share and Cite

MDPI and ACS Style

Cervi, F.; Dadomo, A.; Martinelli, G. The Analysis of Short-Term Dataset of Water Stable Isotopes Provides Information on Hydrological Processes Occurring in Large Catchments from the Northern Italian Apennines. Water 2019, 11, 1360. https://doi.org/10.3390/w11071360

AMA Style

Cervi F, Dadomo A, Martinelli G. The Analysis of Short-Term Dataset of Water Stable Isotopes Provides Information on Hydrological Processes Occurring in Large Catchments from the Northern Italian Apennines. Water. 2019; 11(7):1360. https://doi.org/10.3390/w11071360

Chicago/Turabian Style

Cervi, Federico, Andrea Dadomo, and Giovanni Martinelli. 2019. "The Analysis of Short-Term Dataset of Water Stable Isotopes Provides Information on Hydrological Processes Occurring in Large Catchments from the Northern Italian Apennines" Water 11, no. 7: 1360. https://doi.org/10.3390/w11071360

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop