A record of changes in the Gran Sasso groundwater before, during and after the 2016 Amatrice earthquake, central Italy

We performed continuous recordings (May 2015 – January 2017) of hydraulic pressure and electrical conductivity of groundwater in the 190 m-long horizontal S13 borehole drilled next to the deep underground laboratories of Gran Sasso (LNGS-INFN), located in the core of the Gran Sasso carbonate aquifer (central Italy) at a distance of about 39 km south-eastward from the 24 August 2016 Amatrice earthquake (6.0 Mw) epicenter. Using a 3-channel, 24-bit ADC we achieved a sampling rate of groundwater physical properties up to 50 Hz for each channel. We focused on the analysis of data recorded before, during and after the Amatrice earthquake, describing and discussing in detail the evidence for significant hydraulic pressure and electrical conductivity anomalies recorded before the main shock. We identified unambiguous signals in the hydraulic pressure data starting on 19 August, i.e. five days before the 24 August mainshock. A more careful analysis allowed us to detect the inception of a weak change up to 40 days before the Amatrice earthquake and a significant variation in the electrical conductivity data about 60 days before. The data revealed highly dynamic aquifer behaviour associated with the uprising of geogas probably related to the preparation stage of the Amatrice earthquake.

L' Aquila mainshock and inferred that fluids played a key role in the fault failure process [44][45][46][47][48] . Similarly, during the 24 August 2016, M w 6.0 Amatrice earthquake there was evidence of hydrogeochemical and hydrogeological changes before and during the seismic sequence 51,52 , that included about 80,000 events 53 . Several geochemical anomalies were interpreted as reliable seismic precursors for an extensional tectonic setting 51 . We made our observations from a tunnel at 965 m a.s.l. (Fig. 2), in the core of the Gran Sasso aquifer beneath a 1,400 m-thick layer of rocks. The horizontal S13 borehole (Figs 2 and 3) was drilled around the end of the 1980s in a hall excavated in Meso-Cenozoic limestones and dolostones close to the underground laboratories of Gran Sasso (LNGS-INFN). It has a horizontal length of slightly over 190 m and dips gently upward by ~5° (Fig. 3). The S13 borehole intercepts a thrust fault near its end (Figs 2 and 3). The initial 175 m of S13 borehole are tubed with a well casing, while the last 10 m drain inside the upper Triassic dolomite 54,55 .

Hydrogeological Setting
Most of the Abruzzi region is a seismic area characterized by normal faulting earthquakes reaching maximum intensity values up to XI MCS, corresponding to magnitudes close to 7.0; the most destructive known events occurred in 1349, 1461, 1703, 1915 and 2009 56 . The seismic activity of the Abruzzi region is quite like that of most of the central-southern Apennines, with predominant normal faulting earthquakes generated by predominant NW-SE trending fault systems of the Apenninic chain. [57][58][59][60][61][62] . The Gran Sasso chain is in the northern part of the Abruzzi region and corresponds to the main topographic and geological feature of central Apennines: its culmination (Corno Grande, 2,912 m a.s.l.) is the highest peak in peninsular Italy. The LNGS-INFN underground laboratories, with a rock overburden of 1,400 m, are placed in this area together with the S13 borehole site (Fig. 2). No significant historical earthquakes are known in the Gran Sasso chain, but paleoseismological studies reported multiple surface-faulting events; the magnitude of the two most recent strong events (5 th -3 rd cent. BC and 6 th -5 th millennium BC) have been estimated at about 7.0 63 .
During the 6 April 2009, M w 6.3 L' Aquila earthquake and the related seismic sequence, several investigators studied groundwater post-seismic changes (spring discharge, water table and water content of main ions, 222 Rn, CO 2 , water isotopes and Uranium) of the Gran Sasso aquifer [64][65][66][67][68] . Later, the S13 borehole was specifically selected with the goal of continuous monitoring of hydraulic pressure, electrical conductivity and temperature at a very high sampling frequency. This choice was guided by the observation that S13 borehole clearly recorded water table changes induced by the 2009 L' Aquila earthquake 64,65 .
The 1,000 km 2 -wide fractured and fault-partitioned Gran Sasso carbonate aquifer is hosting the nuclear physics underground laboratories of Gran Sasso (LNGS-INFN) and six horizontal perforations of the borehole hall, among which the S13 borehole ( Fig. 3 for details) 54 . This aquifer is a representative carbonate aquifer of the Mediterranean domain and over the past few decades has been the subject of investigations in the framework of several projects concerning hydrogeology, hydrogeochemistry and isotope hydrology [69][70][71][72][73][74][75][76][77][78] . The Gran Sasso aquifer is composed of Meso-Cenozoic basin-to-slope and reef-platform carbonate rocks, which are arranged in a thrust belt geometry formed during the Upper Miocene Apenninc orogeny 69,70 . The thrust belt is subsequenty displaced by Plio-Quaternay normal faults generating by the NE-migration of the Tyrrhenian post-orogenic extensional front 54,55,79 . The Gran Sasso aquifer is hydraulically confined to the north and to the east by the Upper Miocene   54,55 . Groundwater from the core aquifer is partially drained by the two highway tunnels and by the underground laboratories (UL LNGS in Fig. 2); it is important to stress that there are no pumping systems in the whole area, and groundwater is only passively drained from the two sides of the Gran Sasso chain 78 . Since the 1980s, groundwater flowing in the tunnels has been exploited as drinking water. The decrease of spring discharge measured in the period 1979-1990 was interpreted as a transient phenomenon caused by the tunnel excavation. Once the aquifer reached a new steady condition, however, progressive stabilization of discharge was expected 70 Table; HT: Highway Tunnels; UL LNGS: Underground Laboratories; BH: Borehole Hall; CA: Calderone glacier (high elevation water reservoirpreferential recharge area); 1: overflow spring (CP: Capo Pescara spring); 2: preferential groundwater flowpath area; 3: preferential groundwater flowing toward the UL; PR: Preferential Recharge; DR: Diffuse Recharge; S13: monitored horizontal borehole. The hydrogeological relationships in the square are showed into details in (b). (b) and (c) Detailed hydrogeological relationships between Calderone glacier acting as a water reservoir for the carbonate aquifer down below; i: hydraulic gradient; k: hydraulic conductivity (kh: karst horizon; lim: limestone; dol: dolomite). The hydrogeological relationships in the square are showed in detail in (c).  80 . Figure 4 shows the linear behaviour and fit of hydraulic pressure vs its flow of S13 borehole, Fig. 5 shows a simple hydraulic scheme of the apparatus located at the S13 borehole wellhead. The hydrological characteristics of the S13 borehole showed a specific peculiarity: high hydraulic pressure with respect to very low flow rate (Fig. 4). To measure temperature and electrical conductivity we needed to spill some water, but to avoid disturbance in hydraulic pressure measurement we decided to limit as much as possible the spilled flow. To this end, we installed the hydraulic system shown in Fig. 5; a manual hydraulic valve (d) in Fig. 5 permitted to regulate the water flow. The hydraulic pressure sensor was connected directly to the S13 horizontal borehole output, (c) in Fig. 5, while the other two sensors (temperature (e) and electrical conductivity (f)  [54]). The figure shows the characteristic of the 6 horizontal boreholes (S13, S14, S15, S16, S17 and S18). The S13 borehole was monitored in this study. Bottom: geological logs of S13 and S14 boreholes (top scheme and  (Fig. 6). Under these conditions of high pressures and outlet of water from an extremely small opening, the flow rate was not in a laminar regime but in a fully developed turbulence.
A seismic station of the national network 53 of Istituto Nazionale di Geofisica e Vulcanologia (INGV) is located about 250 m from the S13 borehole inside the underground laboratory of Gran Sasso (LNGS-INFN). The seismic station (international code: GIGS) was deployed in the framework of the GINGER experiment [81][82][83] and is equipped with two broadband seismometers, including a Nanometrics Trillium 240 s (see http://iside.rm.ingv. it/iside/standard/info_stazione.jsp?page=sta&sta=2571 for further details) and a Guralp CMG 3 T 360 s. This instrumentation is used both for continuous microseismic monitoring of the Gran Sasso chain and for recording global seismicity.

Data Discussion and Interpretation
During the months of continuous monitoring (May 2015-Jan 2017) the hydraulic pressure signal of the S13 horizontal borehole showed the recharge and depletion cycle of the Gran Sasso aquifer (Fig. 6) and displayed also the sun/moon tides (Fig. 7). The temperature followed both seasonal and daily variations, the former linked always to the hydrogeological cycle, the latter due to the closeness of the highway tunnel with the S13 hall, while the electrical conductivity did not show any significant periodical variations. Apart from these expected trends, at the end of 2015 the hydraulic pressure data showed coseismic effects related to two large earthquakes occurring during the monitored period, respectively located at global and regional distance: the 16 Sept 2015, Illapel, central Chile earthquake (M w 8.3, distance about of 12,000 km), and the 17 Nov 2015, Lefkada, western Greece earthquake (M w 6.5, distance about 700 km). Details of the waveforms recorded can be found in [80]; the arrival of Love and Rayleigh waves from the Chile earthquake and P-and S-waves from the Greece earthquake are very clear 80 . The signal recorded in the pressure channel during the earthquake represents the response induced by wave arrivals exciting the whole aquifer.
A few days after the 24 August 2016 (01:36:32 UT, doy 237) Amatrice earthquake (M w 6.0), we recovered our data to search for any changes possibly induced by the earthquake. We expected the signals to be large due the relative closeness of the Amatrice event (about 39 km to the monitored site) and its notable magnitude. As we observed an unexpected behaviour only in the hydraulic pressure and electrical conductivity channels, in the following we will refer only to these observations.   Fig. 9a,b). This was true until 19 August 2016 (five days before Amatrice earthquake) when large and asymmetric fluctuations appeared, lasting till the end of August 2016 (Figs 8 and 9b). In the period from May 2015 to July 2016 we did not observe any variations. In Fig. 8 and in a minor way also in Fig. 6, moreover, we observe a variety of steps. In most cases the steps are up or down with durations of less than 1 minute. These are probably related to pore-pressure phenomena. Coseismic steps are comparatively slower (about ten minutes) 9,13,16,18-21,23-25 . We next focus on the regime change that occurred on 19 August 2016 (doy 232) and on the resulting anomalies in the hydraulic pressure signal characteristics from S13 borehole. Figure 10a shows the onset of this anomalous behaviour on 19 August 2016. A transition from one regime to a different regime was recorded over a period of about 6 hours (from 14:00 to 20:00 UT). This change occurred not progressively but following a series of "flip-flop" cycles. This change is marked by two distinct features: an increase in the average hydraulic pressure (about 0.004 MPa), and a more dramatic change in the fine structure of the signal. The latter effect can be better appreciated by looking at the two panels (Fig. 10b,c) representing the behaviour of hydraulic pressure data before and after the transition. The presence of asymmetric fluctuations (negative micropulses) is evident following   the transition (Fig. 10b). This behaviour is shown more quantitatively in the histograms displayed in Fig. 10c. From 19 August 2016 (doy 232) onwards, this behaviour is observed systematically, although to a varying degree. Figure 11 shows the hourly number of negative micropulses (red line) from 19 August to 2 September 2016, while the blue line shows the integral trend on an arbitrary scale. It is of interest to note that, starting from about 36 hours before the 24 August mainshock there is a large increase of negative micropulses, culminating at the time of occurrence of the Amatrice earthquake. Moreover, in the period from May 2015 to July 2016 we did not observe any negative micropulses.
To explain the presence of these negative micropulses we go back to the details of the hydraulic setup. The flux of water from the S13 borehole is determined by the balance of two effects: the high impedance of the source associated with the high level of fracturing of the rocks crossed by the S13 borehole, and the small area of the valve orifice (O(1) mm 2 ). The density and viscosity of the water limit the flow rate to about 2 l/min with an inlet/outlet pressure ratio of about 25 (and flow speed of about 10 m/s). Any increase in flow rate implies a sudden decrease of hydraulic pressure as showed by the trend of Fig. 4. Then when a macroscopic bubble (O(1) ml at 1 bar) of geogas reaches the orifice of the valve an instantaneous increase of flow rate occurs due to the sudden decrease of density and viscosity. This results in a rapid decrease of the hydraulic pressure, as measured above the valve, producing a negative micropulse in our recording (on the order of 2 × 10 −3 MPa). Therefore, the counting rate of negative micropulses is a proxy for the abundance of gas bubbles in groundwater. Figure 12 focuses on 24 August 2016 from 01:36:30 to 01:38:00 UT, when the Amatrice earthquake occurred. Data from the 3-components GIGS broadband seismic station are shown (red lines) next to the hydraulic pressure signal (black line) from S13 borehole as a reference. Figure 12 shows a main feature: the remarkable variation of hydraulic pressure (0.2 MPa pp) following the arrival of the S-waves; afterwards, a slow increase of about 0.02 MPa (corresponding to a piezometric head change of about 2 m) were recorded over a period of about 5 hours following the end of the earthquake signal (Fig. 8). Although similar observations have already been described in the literature 9,84 , the results we obtained exhibit a remarkable magnitude and a great level of detail 80 . The novelty of our time series of hydraulic pressure lies also in the presence of unambiguous signals starting several days before the mainshock of 24 August 2016.
To investigate these rather unexpected phenomena in a quantitative fashion, we performed an elementary statistical analysis of the hydraulic pressure time series looking at the first three moments of the complete time-series. More specifically, after subtracting the (60 s) moving average we evaluated the standard deviation (m 2 ) 1/2 , the skewness m 3 /(m 2 ) 3/2 and the kurtosis (m 4 /m 2 2 -3) of the resulting distribution. Figures 13a,b show the results of the statistical analysis of hydraulic pressure data from S13 borehole. The standard deviation (Fig. 13a) is computed from 27 July (doy 209) to 4 September 2016 (doy 248). The Kurtosis and Skewness (Fig. 13b) values are   . The evidence of a signal change was very clear about five days before the main shock, although the change was somewhat detectable starting 40 days before the earthquake when a slow decrease in kurtosis is particularly evident. We checked the relevance of this signal by looking at the kurtosis of the full dataset (1 May 2015 to 5 September 2016) and we found no significant deviation from the zero value until 40 days before the main shock. It is worth to note that the 24 August mainshock was not preceded by any sizable foreshock; particularly, in the circular area of radius 10 km from the M w 6.0 epicenter we counted only one event with M w greater to 1.8 (15 Aug 2016, doy 228, at 19:00:41 UT -M w 1.9) in the period 1 May -24 August 2016 53 . We next address the electrical conductivity data. Figure 14a shows the electrical conductivity signal averaged over 1 hour, starting 140 days before the Amatrice earthquake (occurring at the 0 in Fig. 14a) and lasting 20 days after. A significant change was noted starting 60 days before the mainshock, which can be correlated with the kurtosis analysis of hydraulic pressure (blue line in Fig. 13b). In the period from May 2015 to May 2016 we did not observe any significant variation of electrical conductivity; note again that no significant foreshocks or swarms were present in the area. Figure 14b shows a zoom of three hours (from 00:00 to 03:00 UT of 24 August 2016) of electrical conductivity data averaged every minute; a change in conductivity occurred at 01:51 UT, and the Amatrice earthquake strikes at 01:36. There were 15 minutes of lag time, but considering that the horizontal S13 tube contains about 2,000 l of water flowing at 2 l per minute, this was traced back, to the upper thrust zone, to around 17 hours before the mainshock.
During the period of data analysis and interpretation (between September and the middle of October 2016), which was taken a few months, we made further modifications and tests of our hydraulic system. At the end of October, the seismic sequence restarted with two events of magnitude greater than 5.0 (26 October 2016: 5.4 M w at 17:10:36 UT and 5.9 M w at 19:18:07 UT) culminating on October 30, 2016 with a 6.5 M w (06:40:17 UT) 53 . During this period, only the pressure sensor was continuously recording, but with a different hydraulic port setup, so a direct comparison with older data is difficult. Figure 15a shows the 30 October earthquake (6.5 M w at The only realistic comparison that can be made between these two periods concerns the coseismic effect following the two mainshocks of 24 August (6.0 M w ) and 30 October 2016 (6.5 M w ). The upward coseismic offset related to the two earthquakes is clearly shown by the superimposed plots of the hydraulic pressure (in MPa), with moving average at 1 minute (the average smoothes out the elastic part of earthquake signals, principally P-, S-and Coda-waves) (Fig. 16). The positive coseismic offset is almost the same (the data are aligned and superimposed forcing at (0,0) the start of two earthquakes) (Fig. 16).

Conclusions
The presence of the S13 horizontal borehole near the underground laboratories of Gran Sasso (LNGS-INFN) in central Italy, gave us the unique opportunity to investigate the pristine groundwater reservoir in the core of the Gran Sasso aquifer located in the very seismically active area of central Apennines (Fig. 1). The S13 horizontal borehole was drilled at the end of the 1980s during excavation works for a tunnel across the Gran Sasso chain, located about 39 km south-eastward from the Amatrice earthquake epicenter (M w 6.0) (Fig. 3).
Starting in May 2015, we used the S13 borehole to measure continuously hydraulic pressure, groundwater temperature and electrical conductivity with a high sampling rate (10 Hz to 50 Hz) for each channel. Monitoring was still ongoing when the 24 August 2016 (01:36:32 UT, doy 237), Amatrice earthquake (M w 6.0) hit about 39 km from the study site (Figs 1-3). The acquired record showed distinct and interesting signals in the hydraulic pressure and electrical conductivity data of the Gran Sasso groundwater before, during and after this earthquake. We recognized anomalies that can be classified both as long-term seismic precursors (anomalies of hydraulic pressure signals and electrical conductivity about 40-60 days before the Amatrice earthquake -Figs 13 and 14) and as short-term seismic precursors (about 5 days in the hydraulic pressure data - Figs 8,9,11 and 13). For many years, scientists have observed changes in the amount, direction, and discharge of water flow in streams and in the water level of wells. Furthermore, changes in fluid pressure in the subsurface have been documented following large earthquakes by observing fluctuations as large as several meters in water levels resulting from mid-size regional earthquakes or even from very large earthquakes located thousands of kilometers away. The authors of ref. 9 for example, have supplied a review of water level variations in surface wells coincident with the arrival of teleseismic waves. Our data, however, represent the first observations of earthquake-related signals in hydraulic pressure measured with continuous high-frequency sample rate (20 Hz) monitoring of the S13 borehole placed deep within a Gran Sasso carbonate aquifer; S13 also intersects a thrust fault near its end (Figs 2 and  3).
Groundwater variations are believed to respond to crustal deformation processes well before the occurrence of significant earthquakes 20,21 . Several mechanisms have been proposed to explain such variations during the seismic cycle, including the dilatancy theory which assumes that, during the pre-seismic phase, the rock volume increase, prior to rupture, is favoured by the micro-fracturing within the earthquake source volume. Micro-fracturing induces also rock permeability increase favouring groundwater variations and mobilization or uprising of geogas 9,50,85 .
Groundwater monitoring at the S13 borehole revealed the existence of highly dynamic aquifer behaviour. The data collected during these twenty-one months turned out to be quite a useful tool to: (i) characterize the background in groundwater monitoring, and (ii) identify more clearly the observed anomalies and their possible link with phenomena of tectonic interest. In this respect, the anomalies observed in the hydraulic pressure 5 days before the 24 August 2016, Amatrice earthquake are obvious.
We interpret the presence of negative micropulses in hydraulic pressure data (Fig. 10b) as the result of a possible uprising of deep geogas, mainly CO 2 , consistent with recent studies suggesting a link between fault mechanics, seismicity, underground fluids dynamics and CO 2 uprising [86][87][88][89][90][91][92][93][94] . More precisely, during the L' Aquila earthquake (April 6 th , 2009, 6.3 M w ), the evolution of seismicity was being driven in part by the poro-elastic response of trapped reservoirs of high-pressure fluid, presumably CO 2 , and postseismic fluid flow initiated by the main shock 44 .
Further investigation of the relationship between earthquakes and changes in groundwater parameters near large seismogenic faults, are needed for a full understanding of preseismic, coseismic and postseismic processes.