Research ArticleGEOCHEMISTRY

Correlation between tectonic CO2 Earth degassing and seismicity is revealed by a 10-year record in the Apennines, Italy

See allHide authors and affiliations

Science Advances  26 Aug 2020:
Vol. 6, no. 35, eabc2938
DOI: 10.1126/sciadv.abc2938


Deep CO2 emissions characterize many nonvolcanic, seismically active regions worldwide, and the involvement of deep CO2 in the earthquake cycle is now generally recognized. However, no long-time records of such emissions have been published, and the temporal relations between earthquake occurrence and tectonic CO2 release remain enigmatic. Here, we report a 10-year record (2009–2018) of tectonic CO2 flux in the Apennines (Italy) during intense seismicity. The gas emission correlates with the evolution of the seismic sequences: Peaks in the deep CO2 flux are observed in periods of high seismicity and decays as the energy and number of earthquakes decrease. We propose that the evolution of seismicity is modulated by the ascent of CO2 accumulated in crustal reservoirs and originating from the melting of subducted carbonates. This large-scale, continuous process of CO2 production favors the formation of overpressurized CO2-rich reservoirs potentially able to trigger earthquakes at crustal depth.


The relation between seismicity and fluid circulation has been demonstrated considering both the potential role of deep fluids in earthquake generation and the close spatial correlation between seismic zones and areas of deep fluid discharge (1, 2). The role of deep fluids in earthquake triggering has been recognized in different tectonic settings, including compressional (35), strike-slip (6, 7), and extensional (812) regimes. Besides, the involvement of CO2-rich fluids in seismic sequences has been proven (6, 914). A worldwide spatial correspondence between seismically active areas and CO2 emissions has been demonstrated (15, 16) as well as the primary role of extensional tectonics on Earth degassing of CO2 (2). The few measured cases, Central-Southern Italy (17, 18) and Eastern African rift (19, 20), indicate the global relevance of tectonic CO2 fluxes, which are thought to have controlled the CO2 concentration of the atmosphere and thus of the climate (21).

The temporal variation of the rate of tectonic CO2 discharge with seismicity is, however, still unconstrained because long time series of CO2 discharge in seismically active areas have not been collected. Here, we show a 10-year-long (2009–2018) record of tectonic CO2 emissions compared with seismicity in Central Apennine (Italy), an area where devastating historical earthquakes have occurred [e.g., the 1461 Me (equivalent magnitude) 6.4 and 1703 Me 6.7 events (22);]. The 2009–2018 observation period is a period of repeated seismic sequences (Fig. 1) that are considered here as a unique 10-year-long sequence (Central Apennine Seismic Sequence, hereinafter referred to as CASS). CASS includes the devastating main shocks of 6 April 2009 [Mw (moment magnitude) = 6.3, L’Aquila earthquake], 24 August 2016 (Mw = 6.0, Amatrice earthquake), and 30 October 2016 (Mw = 6.5, Norcia earthquake). These earthquakes and their stronger aftershocks destroyed over a large area of Central Italy with more than 600 deaths, 2000 injured, and evacuation of about 120,000 persons. The three major earthquakes occurred at similar depths (8 to 12 km) with thousands of aftershocks concentrated in the upper 10 to 12 km of the crust. These shocks were characterized by normal focal mechanisms with a northwest-southeast (NW-SE) rupture strike consistent with that of the major Apennine faults (11, 23) that move in response to a northeast-southwest (NE-SW) extension. In the Apennine, crustal extension and devastating earthquakes occur in the overriding plate in response to regional uplift related to an uprising and eastward-moving mantle wedge (17, 2427). Our 2009–2018 degassing record of tectonic CO2 flux shows that more than 1800 kt of deeply derived CO2 were released from a relatively small area (~700 km2) with a rate varying in time accordingly to the evolution of the seismicity.

Fig. 1 Map of the seismicity and of CO2 Earth degassing and location of investigated aquifers.

(A) Map of the M > 2 1986 to 2020 Italian seismicity recorded by the INGV (Istituto Nazionale di Geofisica e Vulcanologia) seismic network (; last accessed January 2020). (B) Studied aquifers and springs and the post-2007 seismicity. The main flow lines of the groundwaters are reported as arrows. Map of the CO2 flux, modified from (17, 29), is shown in the background of both (A) and (B). Coordinates are reported in UTM-WGS84-32N.


Record of tectonic CO2 emission compared with seismicity

Abundant, diffuse CO2 degassing affects the western side of the Apennine chain (17) and confirms the occurrence of melting of subducting plate carbonates in the mantle (25). The upward migration of CO2-rich fluids generates two large degassing structures at the surface: the Tuscan Roman degassing structure (TRDS) and Campanian degassing structure (CDS) that extend from the Tyrrhenian Sea to the Apennines (Fig. 1) (17). The seismic activity ( mainly concentrates in the axial sector of the chain, at the eastern boundary of TRDS and CDS (Fig. 1). This intriguing spatial distribution, in light of the degassing structures, suggests that the gas, accumulated in buried highly pressurized structural traps, may favor seismicity (17, 27), because high-pressure CO2 at depths can promote fault weakening mechanism (9, 28). CASS began at the eastern border of the TRDS (Fig. 1), where a deep source supplying endogenous carbon dioxide to the groundwater is detected.

The geochemical investigations aimed to detect and quantify the possible earthquake-related emission of deep CO2 started immediately after the April 2009 mainshock (29); in particular, we focused on the CO2 dissolved in two large carbonate aquifers close to the epicentral areas (Gran Sasso and Velino aquifers; Fig. 1B). To date, the geochemical dataset includes 270 chemical and isotopic compositions of spring waters sampled from April 2009 to December 2018 (see Materials and Methods). These samples were collected at the main 36 springs (fig. S1) whose flow rate amounts to ~70% of the total discharge of the two aquifers (~80 m3 s−1).

Concentration of deeply derived CO2 (Cdeep) is computed from the total dissolved inorganic carbon (TDIC) content of the water samples by applying a carbon mass balance (30) based on the carbon isotopes (see Materials and Methods). Remarkable variations are observed for three springs of Velino aquifer that were monitored for the entire period (Peschiera, Canetra, and San Vittorino), while the five springs of the Gran Sasso aquifer, monitored until 2015, show either a weak signal or no substantial variation in Cdeep when compared to the Velino springs (fig. S2). Cdeep of these latter springs follows the temporal evolution of the seismicity (Fig. 2, A to C). Peak values of Cdeep occur concurrently with the main shocks [i.e., the devastating events of Mw 6.3 April 2009 and Mw 6 and Mw 6.5 August to October 2016; Fig. 2], and then CO2 emissions decrease following the seismicity decay in terms of magnitude and rate. The robustness of the correlation between CO2 and seismicity was investigated by testing the null hypothesis of no correlation (Pearson’s correlation coefficient equal to 0) through the Student’s t distribution (31). The Cdeep is compared with the earthquake number, and the seismic energy release occurred at variable distance from the springs in a period centered at the sampling date ± a time lag. The correlation is statistically significant (significance level, 0.01) for a wide range of values for time lag (>10 days) and distance (20 to 80 km) (see Materials and Methods and fig. S3). This correlation is also well highlighted in the binary plots of Fig. 2 (A to C), where the Cdeep content of the three springs is compared with the number of earthquakes for a time lag of 40 days, at a distance <45 km. It is worth noting the high correlations of the earthquake number with the Cdeep of San Vittorino spring (the water that showed the strongest variation, R2 = 0.80) and Peschiera spring (the spring of highest flow rate, R2 = 0.62).

Fig. 2 Temporal evolution of CO2 degassing and seismicity.

(A to C) Chronograms of the earthquake magnitudes compared with the concentration of deeply derived carbon (Cdeep) in the monitored springs. (D) Binary plots of Cdeep of the monitored springs against the log number of the earthquakes occurred at distances <45 km in a period of 80 days centered at any sampling date. (E) Chronogram of the earthquake magnitudes compared with the daily amount of deeply derived CO2 dissolved by the groundwaters of the Velino aquifer (FCO2).

Although geochemical anomalies in groundwaters associated with CASS were already reported in many studies [e.g., (3234) and references therein], we produce an unprecedented record of tectonic CO2 emission (Fig. 2E) by applying a carbon mass balance to the Velino aquifer springs of known flow rate (see Materials and Methods). The daily flux of deeply derived CO2 (FCO2), similarly to Cdeep of the springs, shows a statistically significant correlation with seismicity (significance level, 0.01; see Materials and Methods and fig. S3). The total cumulative amount of the CO2 released during the CASS period (10 years) is of ~1800 kt, i.e., of the same order of magnitude of the CO2 involved in volcanic eruptions [see Table 2 in (35)]. It is worth noting that this amount is the minimum estimate of the total CO2 involved, because it does not include neither the degassed fraction of CO2 (i.e., part of the CO2 is directly emitted in the atmosphere by gas emissions; see movie S1) nor the gas emitted in other areas.

Origin of the CO2

Decarbonation triggered by frictional heating in the slip zone of carbonate-hosted faults has been suggested as a possible CO2 release mechanism during earthquakes [e.g., (36, 37)]. However, this is not the case for CO2 at the Velino aquifer, because the thermometamorphic CO2 would have a very positive carbon isotope signature, as computable by theoretical fractionation and measured in laboratory experiments [δ13C > +3 per mil (‰); (36)]. Instead, the Velino waters show the input of CO2 with a slightly negative carbon isotope signature (δ13C ~ −1.5‰) during the entire observation period (Fig. 3), therefore excluding a major contribution of CO2 formed by frictional heating during the earthquakes.

Fig. 3 CO2 sources and dissolved carbon in the monitored springs.

δ13Cext versus Cext diagram of the monitored springs from the Velino aquifer. Measurements are compared with the theoretical compositions obtained adding a deep inorganic CO213C = −1.48‰) to normal groundwaters. Samples taken in periods close to the main earthquakes (red symbols) show evident shifts toward the same deep CO2 component (i.e., of similar isotopic carbon signature) that enter the aquifer during the quiescence periods (white symbols).

The gas source must be sought in the westward subduction of carbonates of the Adriatic plate beneath the Tyrrhenian mantle. According to (25), melting of such carbonates, due to the high temperatures at the interface between the mantle wedge and the subducting plate, generates carbonate-rich melts. These low density and viscosity melts can migrate upward through the mantle causing zones of partial melting (and CO2-rich compositions) whose presence at depths from 130 to 60 km is revealed by shear-wave velocity images (25). The lower pressures at shallower depths (<60 km) may induce massive CO2 degassing and accumulation beneath the Moho and within the lower crust (25). The presence of a crustal CO2 reservoir at 10- to 15-km depth in the L’Aquila basin, where the CASS started in 2009, is indicated by a sill-like negative anomaly of P-wave seismic velocity (Vp) (Fig. 4) and low ratio of P and S-wave seismic velocities (Vp/Vs) zone [see Fig. 4 in (38)]. This zone was interpreted as a CO2 storage zone (39) fed by fluids from the underlying mantle wedge. The feeding process explains the high geothermal advective heat fluxes (~200 to 300 mW m−2) transmitted through the aquifers located in the area (39). The deep CO2 storage zone is underneath a great portion of the Velino aquifer (Fig. 4), where CO2 emission has been documented. It is worth noting that the total thermal energy associated with the CO2 emission during 2009–2018 is two orders of magnitude larger than the total seismic energy released by CASS (~6 × 1023 versus ~6 × 1021 erg; see Materials and Methods).

Fig. 4 VP anomalies and source of the CO2.

The image highlights a negative anomaly in the seismic velocity variations (DVp %) at 16-km depth [P-wave model (38)]. We interpret this anomaly as a deep gas storage zone that provides deeply derived CO2 to the above located Velino aquifer. The top of the low VP zone is located at a depth of 10 to 15 km [see Fig. 7 in (39)]. Coordinates are reported in UTM-WGS84-32N.

The deep origin of the tectonic CO2 detected at Velino aquifer is also in agreement with previous seismological studies that highlighted the primary role of deep fluids in triggering and controlling the temporal evolution of CASS and of previous seismic sequences (10, 11, 40, 41). For example, gas pulses from deep, highly pressurized reservoirs are recognized as the main source of the 1997 Umbria-Marche aftershocks (9) and of the 2013 Sannio-Matese earthquakes, these latter related to gas release from a magma intrusion in the lower crust (12).


In the Apennine chain, the crustal CO2 reservoirs such as that located below the Velino aquifer (Fig. 4) are fed by the ascent of CO2-enriched slab-derived fluids (25), a process implying a continuous mechanism of CO2 production and accumulation in the crust. This ceaseless movement of CO2 from depth causes pressurization in the reservoir and consequent transfer of gas to the uppermost crustal layers and, ultimately, CO2 saturation-oversaturation of the overlying aquifer(s). At regional scale, this mechanism explains the formation of the large CO2 degassing structures (TRDS and CDS in Fig. 1) and the formation of numerous emissions of free gas (Fig. 1). The occurrence of CO2 saturated-oversaturated waters at the Velino basin is testified by the Terme di Cotilia springs, which discharge hundreds of liters per seconds of CO2 oversaturated waters, and by the direct emissions of a CO2-rich gas phase (movie S1). Different mechanisms have been proposed for the transfer of CO2-rich fluids from depth to surface. These include the dynamic, long-term fluid pumping of fluids related to grain-size reduction through the nucleation phase in the ductile mantle-crust shear zone (42), the ascent of fluids along deep-seated faults (43), the fluid pressurization by self-sealing processes and release (fault-valve) during earthquakes at the base of the seismogenic layer of the crust (44), and the emplacement of magma within the crust (12). Here, we propose that, in addition to the above processes, seismic shaking causes sudden ascent of already separated gas bubbles and, possibly, new gas exsolution from the CO2-rich solutions (mechanical bubbling such as a shaken bottle of a carbonated drink). This process of fluid release is consistent with (i) the observed, almost simultaneous, increase of deep CO2 with earthquakes at the surface (Fig. 2) and (ii) the occurrence at depth of earthquakes triggered by a fluid-related pressure front as suggested by the diffusive spatiotemporal evolution of the earthquakes (11). The increase of tectonic CO2 flux at the surface can be caused by the gas released at shallow depths, while the formation of fluid-related pressure fronts would be explained by the gas separated at depth and by the ascent of increasing amounts of CO2 from the deep reservoir. It is worth noting that this process implies a feedback because seismicity causes gas separation with the formation of ascending pressurized gas fronts that, in turn, favor seismicity.

The ascent of fluids in tectonically active zones may be a passive process favored by fracturing of the crust during earthquakes and/or an active process in which pressurized fluids rise to the surface along preexisting faults, thus triggering earthquakes. Taking into account (i) the close relationship between CO2 discharge and rate/magnitude of the 2009–2019 CASS earthquakes (Fig. 2), (ii) the diffusive processes of seismicity migration and rate (11), (iii) VP and VP/Vs tomographic models (45), (iv) variations in fluid pore pressure during the preparatory phase of the 2009 and 2016 mainshocks (46), and (v) the fluid pressure larger than the regional minimum horizontal stress in the upper crust during the Apennine seismic sequences (10, 11), we propose that the influx of CO2 modulates the evolution of CASS.

We conclude that the ascent of huge amount of slab-derived CO2 that continually accumulates at depth may significantly contribute to earthquake occurrence in the Apennines and, potentially, to the observed uplift (47) of the chain. The analysis of long-term time series of CO2 discharge as presented here should be extended to other seismically active areas to better constrain the role of fluids in modulating the seismic activity and to estimate the overall tectonic CO2 flux on Earth.


Field and laboratory measurements of groundwater

The complete dataset of the waters used in this work is reported in data file S1, while their location is reported in fig. S1. The first survey was performed in 1997, while the other surveys were carried out after the 2009 L’Aquila earthquake. Overall, 36 springs were sampled. Among these, 14 springs were monitored from 2009 to 2015: 4 springs belong to the Velino carbonate aquifer (total discharge Qtot = 30.0 m3 s−1; sampled discharge Qs = 24.75 m3 s−1), and 10 springs belong to the Gran Sasso carbonate aquifer (Qtot = 49 m3 s−1; Qs = 25.6 m3 s−1). The Velino aquifer springs, which showed remarkable variations, were sampled other six times after the earthquake of August 2016. Data collected in 1997 and from April 2009 to February 2010 were previously published (29, 30) (see data file S1).

For each spring, temperature, pH, Eh, electrical conductivity, and total alkalinity were determined directly in the field. The alkalinity determination was performed by acid titration with 0.01 N HCl. Water samples for chemical analyses were filtered with 0.45-μm filters and collected in three 100-ml polyethylene bottles. One aliquot was immediately acidified with HCl 1:1 diluted.

For the determination of δ13C of TDIC (δ13CTDIC), the dissolved carbon species were precipitated in the field as SrCO3 by adding SrCl2 and NaOH in solid state to the water sample. Carbonate precipitates were then filtered and washed with distilled water in a CO2-free atmosphere in the laboratory.

Dissolved ions were analyzed at the laboratory of Earth Science Department of Perugia University. Dissolved anions (Cl, F, SO4, and NO3) were determined by ion chromatography using a Dionex DX-120, Ca and Mg were determined by atomic absorption flame spectroscopy on the acidified sample, while Na and K were determined by atomic emission flame spectroscopy. All the laboratory analytical methods have an accuracy better than 2%.

Analyses of water hydrogen and oxygen stable isotopes (δD and δ18O) and of δ13CTDIC were performed at INGV (Istituto Nazionale di Geofisica e Vulcanologia, Sezione di Napoli), using a Finnigan Delta plus XP continuous flow mass spectrometer coupled with a GasBench II device. Each sample was analyzed in replicate. Carbon samples were analyzed versus the Working Standard Marmo Acqua Bianca (MAB) (marble, δ13C = 2.45‰ and δ18O = −2.43‰ versus Vienna Pee Dee Belemnite, V-PDB) (analytical errors: 1σ < 0.06‰). Oxygen and hydrogen isotope were analyzed versus five Working Standard [standardized versus Vienna Standard Mean Ocean Water (V-SMOW), Greenland Ice Sheet Precipitation (GSIP) and Standard Light Antarctic Precipitation (SLAP), International Atomic Energy Agency (IAEA) standards] using equilibration techniques (analytical errors: 1σ < 0.08‰ and <1.0‰ for oxygen and hydrogen, respectively). The TDIC and CO2 partial pressure (PCO2) were computed using the code PHREEQC (48) starting from the field determinations of T, pH, and alkalinity and the major ion concentrations.

Classification of the groundwaters

All samples show δD and δ18O values typical of local meteoric waters and a Ca(Mg)-HCO3 chemical composition (data file S1). The samples are classified into three groups (groups A, B, and C; Fig. 5A) based on their different TDIC and δ13CTDIC. According to (30), the TDIC of the groundwater of carbonate aquifers is given by the sum of the carbon from dissolution of the carbonate minerals (Ccarb) and the carbon from sources external to the aquifers (Cext). The different origin of Cext allows us to differentiate the normal groundwaters (low TDIC, group A) from group B and C samples (higher TDIC). In detail, Cext is given only by the dissolution of biogenic carbon during the infiltration of the recharging water in group A (Cinf, the carbon concentration of the infiltrating water), whereas Cext also originates from the addition of deeply derived CO2 (Cdeep) to normal groundwater in groups B and C (Fig. 5). The group C waters can be differentiated from group B given its higher amounts of added Cdeep that causes CO2 degassing (PCO2 ~1 bar in group C waters), calcite precipitation, and isotopic fractionation during the carbon removal.

Fig. 5 Carbon dissolved in Velino and Gran Sasso groundwaters.

(A) TDIC versus δ13CTDIC diagrams of the sampled water. Group A groundwaters (normal groundwaters) follow a trend expected for an infiltrating water dissolving an isotopically “light” biogenic CO2. Group B and C groundwaters follow a trend compatible with the input of inorganic CO2, characterized by a heavier δ13C value. The high TDIC and the δ13C values of group C suggest that these groundwaters are affected by CO2 degassing before the sampling. (B) 1/Cext versus δ13Ccext diagram. The diagram allows the estimation of the isotopic composition of the deeply derived CO2 and the infiltrating water end members.

Carbon mass balance

To quantify the contribution of the different carbon sources and to characterize their isotopic composition, we compute for each sample Cext and δ13Cext using the following carbon mass balance equationsCext+Ccarb=TDIC(1)δ13Cext×Cext+δ13Ccarb×Ccarb=δ13CTDIC×TDIC(2)where (i) TDIC and δ13CTDIC are analytically determined; (ii) Ccarb is computed as the sum of Ca + Mg − SO4 considering the dissolution of carbonate minerals (i.e., calcite and dolomite) and the possible presence of gypsum/anhydrite; and (iii) δ13Ccarb is assumed to be constant and equal to the average δ13C of numerous samples of carbonate rocks from the investigated aquifers [+1.8‰; (29) and references therein].

The two equations can be used to determine Cext and δ13Cext only for waters that have not experienced CO2 degassing and calcite precipitation, a condition that, for the studied aquifers, is verified for TDIC <40 mM (29), i.e., only for the samples of group A and group B (Fig. 5A).

The concentrations (Cinf and Cdeep) and isotopic compositions (δ13Cinf and δ13Cdeep) of the carbon from the different sources contributing to Cext of group B samples are computed by considering the following carbon balanceCinf+Cdeep=Cext(3)δ13Cinf×Cinf+δ13Cdeep×Cdeep=δ13Cext×Cext(4)

The system of two equations and four unknown variables is solved for each sample using the binary plot δ13 Cext versus 1/Cext (Fig. 5B) where mixtures among different sources show a linear trend. In particular, in the case of the springs of Velino aquifer:

(i) The isotopic compositions of the deeply derived CO213Cdeep = −1.48‰) is considered unique and derived in Fig. 5B as the δ13Cext intercept at 1/Cext = 0 (i.e., pure CO2 end member) of the best-fit linear regression (δ13Cext = −1.48 − 28.1/Cext) of the samples;

(ii) Cinf is determined at the interception of the line connecting each sample to the deep source (1/Cext = 0, δ13Cext = −1.48‰) with the infiltrating water line (δ13Cext = −25.0 + 13.2/Cext) computed as the best-fit linear regression of group A (normal groundwater). This sample by sample computations gives a narrow range of Cinf values, from 1.5 to 2.18 mM with a mean Cinf = 1.76 ± 0.18 mM;

(iii) the carbon concentration of the deeply derived CO2 (Cdeep) of each sample (Fig. 2, A to C) is given by inserting the computed mean Cinf in Eq. 3.

Note that, to have a minimum estimate of the Cdeep, a simplified computation is applied also to the group C Terme di Cotilia spring, where part of the carbon is removed during degassing and calcite precipitation. In details, we compute Cext from Eq. 1 and Cdeep from Eq. 3, assuming the mean Cinf value derived for the group B waters.

Total CO2 and energy output from Velino aquifer

The total deeply derived CO2 entering the Velino (FCO2) is computed by summing the contributions of the monitored springs and applying a correction to account for the portion of the nonsampled groundwaters. The different contributions are computed by multiplying the Cdeep of the four monitored springs for a corresponding water discharge. For this computation, we consider the Peschiera, Canetra, Terme di Cotilia, and San Vittorino springs as representative of the discharge from groups of springs of similar compositions sampled in the detailed survey of January–February 2010. We consider a flow rate of 18 m3 s−1 for Peschiera, 6 m3 s−1 for the Canetra group, 0.25 m3 s−1 for Terme di Cotilia, and 0.5 m3 s−1 for the San Vittorino group according to data from detailed hydrogeological studies (49, 50). Last, we scale the computed total CO2 of the monitored springs (24.75 m3 s−1) to the total discharge of the aquifer (i.e., 30 m3 s−1).

The total energy associate with the deep CO2 influx in the Velino aquifer is estimated starting from the good correlation between the deep CO2 and the advective geothermal heat computed by (39) and (18) for 11 aquifers of the Central Apennine [which includes also the Velino aquifer, named Marsica N in (39)]. The good correlation suggests that the gas and the heat are advectively transported into the aquifers by hot CO2-rich fluids. The estimated CO2 flux from the data of each survey was used to estimate the thermal energy by applying the relation between the FCO2 (in t day−1) and the geothermal heat flow (QH in erg day−1), QH = 3.3 × 1017 FCO2, computed by the data reported in (18).

Correlation between deep CO2 and seismicity

To track the deep CO2, we consider four time series: Cdeep from Peschiera, Canetra, and San Vittorino springs and FCO2. The earthquake catalog includes seismic events with Mw >2.0 from September 2007 to August 2019. We verified the completeness of the catalog for this period and magnitude range. For each Cdeep value, we extract from the catalog earthquakes within a given time lag before and after the sampling date and within a given distance from the corresponding spring (from Peschiera spring in the case of FCO2), evaluating the number and the energy release (51) of the extracted events. In this way, for each Cdeep and FCO2 time series and each couple time lag–distance, we build two seismic time series considering the logarithm of the earthquake number and of the total energy release.

We test the correlation between each Cdeep-FCO2 series and each corresponding seismic time series (log number and energy) for a wide range of time lags (from 7 to 70 days) and distance (from 15 to 100 km). In each case, we test the null hypothesis of no correlation (Pearson’s correlation coefficient equal to 0) through the Student’s t distribution (31). As shown in fig. S3, the null hypothesis of no correlation can be rejected at significance levels of 0.05 and 0.01 for all time lags (>15 days) and distances (>30 km), with the exception of the case of Cdeep at Canetra. In this specific case, the null hypothesis can be rejected at this significance level of 0.01 only for a limited set of time lags and distances for the log seismic energy and cannot be rejected for the number of earthquakes. Also considering the multiple tests performed, the overall correlation between deep CO2 and the seismicity can be considered statistically significant.


Supplementary material for this article is available at

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.


Acknowledgments: We thank the editor G. Gaetani and three anonymous reviewers for their comments and suggestions that improved the manuscript. Funding: This work was supported by the MIUR, project no. PRIN2017-2017LMNLAW “Connect4Carbon.” Author contributions: G.C. originally conceived the study and wrote the manuscript with the help of C.C., G.V., F.D.L., F.F., and J.S. G.C., C.C., F.F., and S.C. elaborated the geochemical data. F.D.L. and G.V. interpreted the seismological data. J.S. performed the statistical analyses. C.C., A.R., and G.B. performed the fieldwork. S.C. performed the isotopic analyses. All the authors reviewed the manuscript. Competing interests: The authors declare that they have no competing interests. Data and materials availability: All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Additional data related to this paper may be requested from the authors.

Stay Connected to Science Advances

Navigate This Article