# Methane Occurrence and Quantification in a Very Shallow Water Environment: A Multidisciplinary Approach

Giulia Matilde Ferrante, F. Donda, Valentina Volpi, Umberta Tinivella, Mathilde Adelinet

<figure><img src="image_1.png" /></figure>

## To cite this version:

Giulia Matilde Ferrante, F. Donda, Valentina Volpi, Umberta Tinivella, Mathilde Adelinet. Methane Occurrence and Quantification in a Very Shallow Water Environment: A Multidisciplinary Approach. Geochemistry, Geophysics, Geosystems, 2022, 23 (3), pp.e2021GC010276. 10.1029/2021GC010276. insu-03643038

HAL Id: insu-03643038

https://insu.hal.science/insu-03643038

Submitted on 25 Apr 2022

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Geochemistry, Geophysics, Geosystems®

## RESEARCH ARTICLE

10.1029/2021GC010276

**Key Points:**

* Integration of seismic and well-log data allows to evaluate the gas occurrence when the shallow P-wave seismic velocity field cannot be obtained

* Resistivity anomalies can be used to assess the gas occurrence and its concentration

* Methane occurrence and quantification in the sediments is a crucial step to evaluate the potential impact of its leakage in the marine environment

**Supporting Information:**

Supporting Information may be found in the online version of this article.

**Correspondence to:**

F. Donda, 
fdonda@inogs.it

**Citation:**

Ferrante, G. M., Donda, F., Volpi, V., Tinivella, U., & Adelinet, M. (2022). Methane occurrence and quantification in a very shallow water environment: A multidisciplinary approach. *Geochemistry, Geophysics, Geosystems*, 23, e2021GC010276. https://doi.org/10.1029/2021GC010276

Received 21 NOV 2021

Accepted 27 DEC 2021

**Author Contributions:**

**Conceptualization:** G. M. Ferrante

**Data curation:** G. M. Ferrante

**Formal analysis:** G. M. Ferrante, V. Volpi, U. Tinivella, M. Adelinet

**Funding acquisition:** F. Donda

**Investigation:** G. M. Ferrante, V. Volpi

**Methodology:** V. Volpi, U. Tinivella, M. Adelinet

**Project Administration:** F. Donda

**Resources:** F. Donda

**Software:** G. M. Ferrante, V. Volpi, U. Tinivella

## Methane Occurrence and Quantification in a Very Shallow Water Environment: A Multidisciplinary Approach

G. M. Ferrante¹ Ⓓ, F. Donda² Ⓓ, V. Volpi² Ⓓ, U. Tinivella² Ⓓ, and M. Adelinet³

¹IPGP, Paris, France, ²Istituto Nazionale di Oceanografia e di Geofisica Sperimentale – OGS, Trieste, Italy, ³IFP Energies nouvelles, Rueil-Malmaison, France

<figure><img src="image_3.png" /></figure>

**Abstract** In the last decades, geological degassing in shallow marine environments has been recognized as a significant contributor to atmospheric methane, hence to the global-warming potential. Especially in shallow water environments, a proper assessment of the distribution, quantification and migration pathways of methane within the sediments is fundamental to help forecasting the amount that could leak and eventually reach the atmosphere. Traditionally, velocity anomalies from seismic data are the ones used to assess the occurrence of gas and its concentration. However, in shallow-waters (<30 m), the post-critical conditions make the near-surface velocity estimation from P-wave reflections extremely challenging, requiring an integrated approach. Here, we propose an original joint analysis of seismic data and geophysical logs, together with information from drilling reports, with the aim of characterizing and quantifying the gas along two crossing multichannel seismic profiles in the Northern Adriatic Sea, a very shallow marine basin where methane occurrence within the sedimentary succession is widespread. We estimated the gas distribution from resistivity anomalies, which are correlated with the seismic response associated with the presence of gas through the signal frequency content. Our results show a different concentration pattern in the two seismic profiles, revealing that gas is both diffuse and concentrated in local accumulations, in agreement with the gas-related features already identified on the seismic data. Gas concentration appears to be locally associated to the tectonic features identified in the area, indicating that faults act as preferential conduits for gas migration, locally reaching the seafloor and seeping in the water column.

**Plain Language Summary** Methane is a powerful greenhouse gas, with a global warming potential that is ca. 23 times that of carbon dioxide. The most important source of uncertainty in the global methane budget is attributable to natural emissions, and geological degassing is today recognized as a major contributor to atmospheric methane. In particular, when methane is released from shallow-water marine sediments, it might easily bypass the water column and reach the atmosphere. Thus, a proper assessment of gas distribution and its quantification within the sediments in such contexts play a crucial role. In this paper, we propose a joint analysis of multichannel reflection seismic data and borehole information in the Northern Adriatic Sea, a very shallow basin, where the presence of gas in the subsurface is widespread. Due to the fact that gas is able to influence the electrical properties of the rock formations, we were able to identify its imprint and provide its amount hosted in the sediments. Our method represents an ad hoc procedure to study methane emissions in every similar marine environment and can be upscaled to help improving the global methane budget estimates.

## 1. Introduction

Methane ($CH_4$) is a powerful greenhouse gas and, over a 20-year period, the global-warming potential of 1 tonne of atmospheric $CH_4$ is similar to ca. 85 tonnes of $CO_2$ (Hartmann et al., 2013). Global $CH_4$ emissions have increased nearly 10% over the past two decades (Schiermeier, 2020) and its atmospheric concentration is now more than 2.5 times above pre-industrial levels, approaching 1.9 ppm. $CH_4$ has a short lifetime in the atmosphere, ca. 10 years (Prather et al., 2012), and hence a stabilization/reduction of emissions will lead in a few decades to a stabilization/reduction of its atmospheric concentration and therefore its radiative forcing (Saunois et al., 2020). The most important source of uncertainty in the global $CH_4$ budget is attributable to natural emissions (Saunois et al., 2020), and particularly to geological seepage. Indeed, once considered a minor natural $CH_4$ source globally, geological degassing is today recognized as a major contributor to atmospheric $CH_4$ (Etiope et al., 2019). A challenging aspect on the quantification of present-day emissions is that seepage sites may be more widespread and abundant than expected, especially in the marine environment (Thornton et al., 2020). It is estimated that 30–45 Tg of $CH_4$ per year are released from marine seeps to the atmosphere (Etiope & Ciccioli, 2009;

© 2022. The Authors.
This is an open access article under
the terms of the Creative Commons
Attribution License, which permits use,
distribution and reproduction in any
medium, provided the original work is
properly cited.
**Supervision:** F. Donda, V. Volpi, U.

Tinivella

**Visualization:** G. M. Ferrante

**Writing – original draft:** G. M. Ferrante, F. Donda

Schwietzke et al., 2016), contributing to the atmospheric CH₄-pool and acting as a powerful greenhouse gas (Mitnick et al., 2018). In very shallow waters (<30 m water depth), this aspect is particularly important since gas bubbles can quickly reach the atmosphere (Borges et al., 2016). A proper identification of CH₄ occurrence in gas-bearing sediments and its quantification is therefore a crucial step to help forecasting the potential contribution of CH₄ escape from the sedimentary succession in the water column and in the atmosphere. Seismic data is a powerful tool to identify fluid accumulation and migration in the subsurface, as evidenced by acoustic anomalies such as bright spots, blanking, chimneys, and pipes (Cartwright et al., 2007; Løseth et al., 2009), but it may also allow for a quantitative estimation of the gas saturation. Traditionally, velocity anomalies are used to estimate the amount of fluid and, in particular, gas in the sediment pore space. In shallow marine environments, however, the estimation of a reliable near-surface velocity distribution from P-waves may be tricky because of the strong interference with unwanted events, such as seabed multiples and guided-waves (Yilmaz, 2001). Reflections at the seafloor can occur at post-critical angles (Tinivella et al., 2017; Verschuur, 2006), preventing the recording of primary water-bottom reflections at the near-offset and further complicating the multiples attenuation. Recently, Giustiniani et al. (2020) show how the integration of different seismic wave types, such as P- and interface waves, can help in constraining shallow sediment characteristics in shallow water environments. However, a petrophysical characterization of the gas-bearing sediments is required to validate the seismic data analysis and provide the physical properties for the gas estimation. The aim of this work is to provide a multidisciplinary approach by means of seismic and well-log data to determine the gas distribution and quantification within the sedimentary succession in the Northern Adriatic Sea (Figure 1), a very shallow marine basin, where both deep- and shallow-source CH₄ are widespread distributed within the sedimentary succession, at the seafloor and in the water column (Donda et al., 2015, 2019; Gordini et al., 2012).

## 2. Gas Occurrence in the Northern Adriatic Sea

In the Northern Adriatic Sea, several gas fields were discovered and exploited during the 1970s. Fluid extraction has been then forbidden for environmental reasons, but the gas system is anyway considered mature (Cazzini et al., 2015). Microbial gases occur within Pliocene-to-Pleistocene turbiditic sands and they are commonly characterized by multiple “pools” within thin sand beds at approximately 1,200–1,500 mbsf (Bertello et al., 2008; Casero, 2004; Casero & Bigi, 2013; Donda et al., 2015). Geophysical data acquired in the last few years clearly revealed the occurrence of both shallow and deep gas in the subsurface. Donda et al. (2013); Donda et al. (2015) identified several acoustic gas-related features, such as signal amplitude anomalies (wipe-out zones and bright spots), or evidence for gas-charged fluid migration pathways, for example, gas chimneys. These features are mainly found within the Plio-Quaternary sequence, although chimneys are also recognized in the Eocene-Miocene succession, vertically extending for up to 2,000 m. Gas is also present in shallow, laterally persistent peat layers, which are widely distributed throughout the Northern Adriatic Sea and formed between 16,000 and 24,000 years BP in an alluvial plain environment (Correggiari et al., 1996; Donda et al., 2019; Gordini, 2009; Zecchin et al., 2011) during periods of relatively warmer climate in the last glacial phase (Donda et al., 2019; Zecchin et al., 2011). This gas is able to seep and produce 4–20 m-high gas flares and intermittent bubbling within the water column. Seeps appear to be locally associated with distinct rock outcrops, irregularly distributed on the seafloor, represented by bio-concretionned carbonate rocks, interpreted as authigenic carbonates and locally named “Trezze” (Gordini et al., 2002, 2012; Donda et al., 2013, 2015, 2019).

## 3. Material and Methods

### 3.1. Geophysical Data

For this work, we processed two perpendicular multi-channel seismic profiles, the STENAP 08 and the GANDI 09, acquired with the main aim to (a) constrain the relationship between the gas emissions and the regional geological setting of the area (Donda et al., 2013) and (b) characterize the gas-charged fluids occurring within the sedimentary succession (Donda et al., 2015, 2019). The acquisition parameters of the two related seismic cruises are reported in Table 1.

Our velocity estimation from the seismic data suffers from the shallow water conditions and it cannot be uniquely used to quantify the gas. However, drilling results show that resistivity logs anomalies correlate with the seismic response associated with the presence of gas and can be used to constrain the velocity information (Figure 2).
<figure><img src="image_4.png" /><figcaption>**Figure 1.** (a) Position map of the STENAP (black lines) and the GANDI (purple lines) multichannel seismic lines. The two analyzed profiles are highlighted in white. Yellow dots: wells used in this study, from west to east: Arlecchino 1, Arcobaleno 1, Tiglia Mare 1, Rachele 1, Amanda 1bis. In green: gas field (from ENI Fact Book Exploration and Production). Gas fields bordered in red are the ones crossed by the two analyzed seismic profiles. (b) Sketch of the regional structural setting, modified from Donda et al. (2015).</figcaption></figure>

To this aim, we derived the main petrophysical properties, namely P-impedance and porosity, from the geophysical log data of 5 hydrocarbon wells: Arlecchino 1, Arcobaleno 1, Triglia Mare 1, Rachele 1, and Amanda 1bis (Figure 1). We used sonic, resistivity, spontaneous potential, and gamma ray logs. We correlated available and calculated well-log curves with a family of seismic attributes to predict resistivity and porosity along the two seismic lines and, taking advantage of Archie's law, we estimated the spatial distribution and concentration of the gas along the two profiles. Figure 3 summarizes the entire procedure with a flow diagram, whose steps are explained below.

## 3.2. Seismic Data Processing

The seismic data are strongly affected by short and long period multiples and several issues arise from the very shallow water condition. Indeed, due to the shallow water column (25 m at maximum), the critical distance for re-

fractions is very small, ~90 m. This means that, apart from the first 4–5 channels, refracted arrivals are present and interfere with reflections (post-critical condition). Wavelets are distinguishable only in the very near offsets and are distorted both in phase and in amplitude. Furthermore, the occurrence of Scholte waves, that is, surface waves propagating along the seafloor, is able to mask deeper reflections (Giustiniani et al., 2020). Therefore, an estimation of a reliable velocity distribution from P-wave is not trivial. As expected, semblance velocity analysis (Yilmaz, 2001) does not provide satisfying results for the STENAP 08 and the GANDI 09 seismic lines. For this reason, we incorporated available sonic-log data in the velocity estimation. We calculated the velocity field of the GANDI 09 using the log information for the shallowest 1.6 s, laterally extrapolated along the seismic line, and we performed semblance analysis for the deepest part. Velocity field of the STENAP 08 is entirely built using sonic-log of Arcobaleno 1 well. These velocity fields are used in the processing flows. Note that no indications of gas-related anomaly are present in the velocity fields. We processed the seismic data by applying two different sequences: (a) for seismic imaging purpose and (b) to preserve the amplitude information. Details concerning the seismic data processing are available in the Text S1 in Supporting Information S1.

<table><caption>**Table 1**
Seismic Acquisition Parameters</caption><thead><tr><th></th><th>STENAP 08</th><th>GANDI 09</th></tr></thead><tbody><tr><td>Source</td><td>Sleeve gun</td><td>GI-gun</td></tr><tr><td>Volume (cu.in.)</td><td>1,180</td><td>840</td></tr><tr><td>Pressure (psi)</td><td>2,000</td><td>2,000</td></tr><tr><td>Source depth (m)</td><td>5 ± 0.5</td><td>5 ± 0.5</td></tr><tr><td>SP interval (m)</td><td>12.5</td><td>25</td></tr><tr><td>Streamer length (m)</td><td>1,200</td><td>1,500</td></tr><tr><td>Streamer depth (m)</td><td>3 ± 0.5</td><td>3 ± 0.5</td></tr><tr><td>N of channels</td><td>96</td><td>120</td></tr><tr><td>Ch. interval (m)</td><td>12.5</td><td>12.5</td></tr><tr><td>Minimum offset (m)</td><td>25</td><td>50</td></tr><tr><td>Fold</td><td>48</td><td>24</td></tr><tr><td>Sample rate (ms)</td><td>1</td><td>1</td></tr></tbody></table>
<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /><figcaption>**Figure 2.** Velocity-resistivity cross-plots for (a) Arcobaleno 1, (b) Arlecchino 1, and (c) Rachele 1 boreholes.</figcaption></figure>

### 3.3. Well-Log Data Analysis

Sonic, gamma-ray, spontaneous-potential and resistivity logs, only available in raster format, were digitized, interpolated with a spline function, smoothed and resampled every 10 m to be consistent with the seismic frequency content. The borehole composite logs do not contain every geophysical log mentioned above. For example, the sonic log was recorded only in 3 of the 5 analyzed wells. Moreover, none of the borehole lies exactly on the seismic profiles. To verify possible correlations between the measured physical properties, we produced cross-plots of each pair of logs. This confirmed the presence of lithological units with similar petrophysical characteristics that can validate the extrapolation of the log data along the seismic profiles. Our analysis revealed that velocity and resistivity have a good cross-correlation in all the wells where the sonic log is available, namely Arcobaleno 1, Arlecchino 1, and Rachele 1 (see Figure 2). Apart from the shallowest formations, resistivity always appears to be increasing with velocity. In the Arcobaleno 1 and Arlecchino 1 cross-plots, three correlated depth intervals can be distinguished: 0–300, 300–1,300, and 1,300–1,800 m. The shallowest interval is characterized by high resistivity, up to 9 Ω·m, quickly decreasing to values <2 Ω·m with increasing velocity, most likely because of compaction. At depth >300 m, resistivity starts to increase with velocity and shows a shift at ~1,300 m, reflecting a change in the petrophysical parameters of the sediments. Rachele 1 well shows resistivities clustered around 1 Ω·m in the shallowest 900 m and an increasing behavior in deeper strata. Velocity-resistivity cross-plots not only show a good correlation, but also suggest that resistivity is able to discern the geological formations. This can be explained by the different electrical properties of both the rock and the fluid saturating the pores. Therefore, even though velocity is the parameter traditionally used to assess gas presence, we decided to use resistivity to detect gas bearing zones. The gas content quantification requires an estimation of the sediments porosity, which can be in turn obtained from a P-impedance inversion procedure. Unfortunately, neither P-impedance (IP) nor porosity profiles are available amongst the selected boreholes logs. We provided an estimation of these properties at well locations (pseudo-logs) through an integrated approach by using all the other log data and the qualitative or sparse information reported in the borehole reports. It should be noted that these are the only available measurements in the study area.

#### 3.3.1. P-impedance Assessment

To estimate IP pseudo-logs, we carried out a lithostratigraphic interpretation of the sedimentary column analyzing the local trend of each log curve and their interrelations. Geological information reported in the well reports and headers guided our analysis and led us to distinguish five lithostratigraphic facies: sandy/peaty clay, clay, sand, marlstone, and limestone. We built a rough preliminary density profile $ρ_0$, increasing with depth, for each well. Values range between 1,700 $kg/m^3$ at ~60 m to 3,200 $kg/m^3$ at ~1,700 m. We then multiplied $ρ_0$ by the available velocity logs to obtain preliminary $IP_0$ curves for each well. We statistically extracted wavelets from the near offset stack sections in the vicinity of each well and we calculate reflectivity and synthetic seismograms at well locations. Then, comparing real and synthetic traces, we iteratively adjust $IP_0$. The lithostratigraphic interpretation is shown for Rachele 1, which is the most recent amongst the selected wells and the one containing the higher quality log curves (Figure 4). Figure 5a shows the final synthetic trace of Arcobaleno 1 well and its comparison with the real trace. Final IP profiles for Arcobaleno 1, Arlecchino 1, and Rachele 1 wells are shown in Figure 5b.
<figure><img src="image_8.png" /><figcaption>Figure 3. Flow diagram summarizing the procedure used to quantify the gas within the sediments.</figcaption></figure>

### 3.3.2. Porosity Assessment

To estimate porosity at well locations, we used the Effective Medium Theory (EMT) with an homogenization approach based on Eshelby's inclusion theory (Adelinet & Le Ravalec, 2015; Eshelby, 1957; Le Ravalec et al., 2014). This type of inversion estimates the microstructural parameter from the known seismic impedances minimizing a given cost function which quantifies the data mismatch in a least-squares sense. The process is repeated until the cost function is small enough. To better constrain the porosity estimation, we inverted separately each of the five facies recognized in the lithostratigraphic interpretation. We modeled every facies considering a solid matrix and a fully saturated porous part. We supposed a mineral composition for the solid matrix and we calculated the correspondent elastic moduli (bulk modulus $K_s$ and shear modulus $G_s$) using the Hill's average of the Voigt and Reuss bounds equations (Schon, 2015):

$$K_s = \frac{1}{2} \cdot \left[ \sum_{i=1}^{m} f_i K_i + \left( \sum_{i=1}^{m} f_i / K_i \right)^{-1} \right] \quad (1)$$

$$G_s = \frac{1}{2} \cdot \left[ \sum_{i=1}^{m} f_i G_i + \left( \sum_{i=1}^{m} f_i / G_i \right)^{-1} \right] \quad (2)$$
<figure><img src="image_9.png" /><figcaption>**Figure 4.** Well-logs and lithostratigraphic interpretation of Rachele 1 borehole.</figcaption></figure>

where *m* is the number of mineral constituents and *f*<sub>i</sub> is the volumetric fraction of the *i*th constituent having *K*<sub>i</sub> and *G*<sub>i</sub> in the solid phase.

We considered pores as isotropically distributed and ellipsoidal-shaped, with a constant aspect ratio for each facies (Bailly et al., 2019). Sand and limestone (facies 3 and facies 5) have spherical pores (aspect ratio = 1); marlstone (facies 4) has quasi-spherical pores (aspect ratio = 0.5) and sandy/peaty clay and clay (facies 1 and facies 2) have ellipsoidal pores (aspect ratio = 0.1). For every facies and every well, we assigned a mixed-fluid (i.e., salt water and gas) composition, under the reasonable assumption of total pore space saturation. We assumed a first-order mixture composition based on gas saturation values measured during drilling and reported in the well reports at sparse depths. In absence of direct information, we assumed that the gas is uniformly distributed in the pore space. In this case, the effective bulk modulus of the mixture is given by the isostress (Reuss) average (Lee, 2004):

$$ \frac{1}{K_f} = \frac{S_w}{K_{sw}} + \frac{S_g}{K_g} \qquad (3) $$

where $S_w$ and $S_g$ are the salt-water and gas partial saturations and $K_{sw}$ and $K_g$ are the respective bulk moduli. Table 3 shows the mineral composition and the elastic moduli assigned to the interpreted facies. Moduli are calculated with Equation 1 using mineral values reported in Table 2. “Common sand” values are revised by Vernik and Kachanov (2010). Table 4 reports the pore fluid bulk modulus $K_f$ of the considered facies. Values are calculated with Equation 4 supposing $K_{sw} = 2.3$ GPa for salt water and $K_{CH_4} = 0.024$ GPa for CH$_4$.

We iteratively inverted IP logs into porosity, treating each facies separately. We performed two different inversions: a first one, assuming a fully water-saturated medium and a second one, considering a mixed-fluid (Figure 5c).

<figure><img src="image_10.png" /><figcaption>**Figure 5.** (a) Comparison between synthetic and real trace at Arcobaleno 1 well location; (b) estimated IP pseudo-logs at wells locations; and (c) estimated porosity pseudo-logs at wells locations. Porosities of fully water-saturated sediments are reported in blue, porosities of mixed fluid-saturated sediments are reported in purple.</figcaption></figure>
<table><caption>**Table 2**
Pore Fluid Bulk Modulus $K_f$ of the Considered Facies</caption><thead><tr><th></th><th>Mineral composition</th><th>K<sub>s</sub> (GPa)</th><th>G<sub>s</sub> (GPa)</th></tr></thead><tbody><tr><td>Facies 1 (Sandy clay)</td><td>70% Clay1<br>30% Common sand</td><td>20.06</td><td>10.55</td></tr><tr><td>Facies 2 (Clay1)</td><td>40% Ilite<br>35% Smectite<br>25% Kaolinite</td><td>15.45</td><td>5.50</td></tr><tr><td>Facies 3 (Common sand)</td><td>100% Quartz<br>30% Common sand</td><td>35.06</td><td>33</td></tr><tr><td>Facies 4 (Marlstone)</td><td>60% Calcite<br>40% Clay1</td><td>40.35</td><td>15.84</td></tr><tr><td>Facies 5 (Limestone)</td><td>100% Calcite</td><td>75.50</td><td>40</td></tr></tbody></table>

It is important to underline that inclusion-based EMT models should be used to estimate high-frequency elastic moduli related to unrelaxed states (Adelinet et al., 2010). To compare these porosity results with seismic-derived data, it is necessary to transfer them from high to low frequency domain. Furthermore, saturated media moduli depend on frequency and the dispersion of the P- and S-waves should be taken into account (Adelinet et al., 2011). However, we did not take into account low frequency effects and our model should be considered as a first approximation. Further data on micro-structure of involved materials could allow a more precise analysis, considering the transition between high and low frequency domain.

### 3.4. Seismic-Well Logs Correlation

We imaged the spatial distribution of petrophysical properties, namely P-wave velocity, resistivity and porosity, along the two analyzed seismic profiles using a seismic-log correlation method called multi-attribute analysis (Emerge-Hampton Russell, CGG Software). The aim of the method is to use the seismic data to predict the property of interest (a new geophysical log) at every CDP location of the seismic profile (Coren et al., 2001). This is done, at each time sample, solving a generalized linear regression between the parameter to be predicted and a wide family of seismic attributes. We collected 10 sample data traces around the well locations and find a relationship between the seismic at those locations and the measured logs. To find the best set of attributes which can optimally describe the target log, we used a step-wise regression (Hampson et al., 2000), which add one attribute per time and perform the optimization using all the wells that constrain the seismic line. A mean squared prediction error was calculated with each attribute addition and stabilized below a threshold of ~5%. Furthermore, we evaluated the correlation coefficient between the original and the modeled log. Then, we assumed that the derived relationship is valid for the entire seismic section and applied that relationship to the entire data set.

For each seismic section, we obtained P-wave velocity, resistivity and porosity distributions. We predicted resistivity in two ways: including (Resistivity 1) or excluding (Resistivity 2) the frequency-related seismic attributes, a class of attributes that is likely to be affected by the gas presence. Furthermore, considering the uncertainties related to our estimated porosity profiles, we applied a further smoothing to the final 2D porosity distribution. Table 5 reports the prediction parameters for every target log of the STENAP 08 and the GANDI 09 seismic lines. P-wave velocity and porosity results are reported in the Text S2 in Supporting Information S1.

### 3.5. Gas Content Quantification

Electrical resistivity data can be extremely valuable for analyzing the subsoil gas behavior and its relation to surface emissions. This is valid for hydrocarbon exploration and reservoir monitoring (Jiang et al., 2021; Senger et al., 2021), for gas hydrates studies (Coren et al., 2001; Kannberg & Constable, 2020; Weitemeyer et al., 2011) but also for evaluating other types of gas emissions, as, for example, from landfills (Johansson et al., 2011). Normally, these gas assessment studies are based on well-logs or CSEM (marine controlled-source electromagnetic) data. Here, we estimated gas saturation using Archie's laws and resistivity theory (Archie, 1942; Glover, 2016). Following the theory, the ratio of the resistivity of a gas-bearing formation to the resistivity of the water-bearing formation is related to the water saturation in the formation when it contains both fluids, that is, gas and water. Therefore, assuming that the pores are fully saturated, the gas saturation $S_g$ is estimated by the following formula:

$$S_g = 1 - S_w = 1 - \left( \frac{R_{bg}}{R_{tot}} \right)^{\frac{1}{n}} \quad (4)$$

where $R_{tot}$ is the gas-bearing sediment resistivity, $R_{bg}$ is the fully water-saturated (gas free) sediment resistivity, $n$ is the saturation exponent and $S_w$ is the fractional water saturation. The gas concentration $C_g$ is the product of gas saturation and porosity:

$$C_g = S_g \cdot \phi \quad (5)$$
<table><caption>**Table 3**

Bulk Modulus ($K_s$) and Shear Modulus ($G_s$) of Some Key Minerals</caption><thead><tr><th></th><th>K<sub>s</sub> (GPa)</th><th>G<sub>s</sub> (GPa)</th><th>Reference</th></tr></thead><tbody><tr><td>Illite</td><td>27</td><td>17</td><td>Greenberg and Castagna (1992)</td></tr><tr><td>Smectite</td><td>9</td><td>1.5</td><td>Mondol et al. (2008)</td></tr><tr><td>Kaolinite</td><td>13</td><td>3</td><td>Mondol et al. (2008)</td></tr><tr><td>Quartz</td><td>36</td><td>44</td><td>Greenberg and Castagna (1992)</td></tr><tr><td>Feldspar</td><td>48</td><td>24</td><td>Wang (2001)</td></tr><tr><td>Calcite</td><td>75</td><td>31</td><td>Greenberg and Castagna (1992)</td></tr></tbody></table>

An estimation of the saturation exponent *n* is needed to solve Equation 5. Generally, *n* is set on the basis of literature values. Here, to get a more constrained estimate of *n*, we exploited the same scattered measurements of gas saturation $S_g$ from the wells as those used in the porosity estimation. The procedure used to estimate *n* for the STENAP 08 seismic line is the following:

1. We extracted depth traces of the STENAP 08 $R_{bg}(z)$ and $R_{tot}(z)$ sections at Arlecchino 1 and Arcobaleno 1 locations.

2. Knowing $S_{z_i}$ at sparse depth $z_i$ from the well reports and using Equation 4, we calculated the saturation exponent $n_{z_i}$, at each $z_i$.

3. From our previous IP estimation, we calculated density profiles $\rho(z)$ at Arlecchino 1 and Arcobaleno 1 locations, in order to extract $\rho(z_i)$, that is, at the same sparse depth $z_i$.

4. We fit the cross-plot between $\rho(z_i)$ and $n(z_i)$ and find an empirical polynomial relationship between $\rho$ and $n$.

5. We used such a relationship to calculate *n* along the whole STENAP 08 seismic line (Figure S3A in Supporting Information S1).

We applied the same procedure to the GANDI 09 seismic line, but we didn't find any satisfying empirical relationship between $\rho$ and $n$. This is probably due to the fact that local measurements of $S_g$ are available only for very shallow depths for Rachele 1 borehole, and are only partly representative of the gas distribution in the whole sedimentary succession. So, because of the lack of constraints, we decided to assign *n* values to the main geological units crossed by the seismic line. We extracted a profile of the STENAP 08 *n* section at the crossing point with the GANDI 09 seismic line. We smoothed this profile and calibrated *n* for the main lithological units, which have been interpreted in both seismic lines (Figure S3B in Supporting Information S1).

## 4. Results

For each seismic section, we obtained P-wave velocity, porosity and resistivity distributions. P-wave velocity values (Figure S1 in Supporting Information S1) show an overall common trend of increasing velocity with depth and, as well as the velocity fields estimated from semblance analysis and well-logs, they do not show significant anomalies. The porosity distributions also have the expected increasing trend with depth (Figure S2 in Supporting Information S1).

<table><caption>**Table 4**

Pore Fluid Bulk Modulus ($K_f$) of the Considered Facies</caption><thead><tr><th></th><th>Arcobaleno 1</th><th>Arlecchino 1</th><th>Rachele 1</th></tr></thead><tbody><tr><td>Facies 1</td><td>97% S<sub>w</sub><br>3% CH<sub>4</sub><br>K<sub>f</sub> = 0.6 GPa</td><td>91% S<sub>w</sub><br>9% CH<sub>4</sub><br>K<sub>f</sub> = 0.24 GPa</td><td>98% S<sub>w</sub><br>2% CH<sub>4</sub><br>K<sub>f</sub> = 0.79 GPa</td></tr><tr><td>Facies 2</td><td>96% S<sub>w</sub><br>4% CH<sub>4</sub><br>K<sub>f</sub> = 0.48 GPa</td><td>90% S<sub>w</sub><br>10% CH<sub>4</sub><br>K<sub>f</sub> = 0.22 GPa</td><td>97% S<sub>w</sub><br>3% CH<sub>4</sub><br>K<sub>f</sub> = 0.6 GPa</td></tr><tr><td>Facies 3</td><td>92% S<sub>w</sub><br>2% CH<sub>4</sub><br>K<sub>f</sub> = 0.79 GPa</td><td>92% S<sub>w</sub><br>8% CH<sub>4</sub><br>K<sub>f</sub> = 0.27 GPa</td><td>99% S<sub>w</sub><br>1% CH<sub>4</sub><br>K<sub>f</sub> = 1.18 GPa</td></tr><tr><td>Facies 4</td><td>99% S<sub>w</sub><br>1% CH<sub>4</sub><br>K<sub>f</sub> = 1.18 GPa</td><td>99% S<sub>w</sub><br>1% CH<sub>4</sub><br>K<sub>f</sub> = 1.18 GPa</td><td>99.5% S<sub>w</sub><br>0.5% CH<sub>4</sub><br>K<sub>f</sub> = 1.56 GPa</td></tr><tr><td>Facies 5</td><td>99.5% S<sub>w</sub><br>0.5% CH<sub>4</sub><br>K<sub>f</sub> = 1.56 GPa</td><td>99.5% S<sub>w</sub><br>0.5% CH<sub>4</sub><br>K<sub>f</sub> = 1.56 GPa</td><td>99.5% S<sub>w</sub><br>0.5% CH<sub>4</sub><br>K<sub>f</sub> = 1.56 GPa</td></tr></tbody></table>

### 4.1. Resistivity Anomalies

Resistivity sections (Figures 6a, 6b, 6d and 6e) show strong anomalies in those areas where seismic data have a pronounced variation in frequency content. Being less conductive than salty water, gas is able to influence the resistivity of formations and to affect the frequency of the seismic data. We suggest that these anomalies are induced by the presence of gas. Therefore, we predicted resistivity in two ways: excluding or including all the frequency-related seismic attributes (Figure 6). In the first case, the estimated resistivity section is thought to be a background resistivity. In the other case, the resistivity section is considered as total resistivity, being representative of the gas-induced anomalies over the background resistivity field.

The STENAP 08 total resistivity panel (Figure 6b) represents the resistivity pattern along the seismic line and it is characterized by strong positive anomalies coincident with the high frequency zones of the seismic line. These positive anomalies occur in the shallowest 0.4 s of the whole seismic profile and are distributed along vertical bands affecting, in particular, a large zone between CDPs 1 and CDP 2500. The vertical distribution pattern can be recognized also around CDP 3200, CDP 6400, and CDP 7500. Considering that the resistivity anomalies are as strong as the main reflectors resistivity, the
<table><caption>**Table 5**

*Prediction Parameters for Well-Log Correlation*</caption><thead><tr><th></th><th></th><th>P-wave (m/s)</th><th>Porosity</th><th>Resistivity (Ω·m)</th><th>Background resistivity (Ω·m)</th></tr></thead><tbody><tr><td rowspan="2">STENAP 08</td><td>Number of attributes</td><td>24</td><td>28</td><td>28</td><td>24</td></tr><tr><td>Correlation</td><td>0.9691</td><td>0.8902</td><td>0.9321</td><td>0.8648</td></tr><tr><td rowspan="2">GANDI 09</td><td>Number of attributes</td><td>24</td><td>6</td><td>28</td><td>10</td></tr><tr><td>Correlation</td><td>0.8865</td><td>0.8789</td><td>0.9519</td><td>0.7097</td></tr></tbody></table>

seismostratigraphy of the line is not visible. Instead, the background resistivity (Figure 6a) correlates much better with the main reflectors and with the seismic amplitude distribution.

As for the STENAP 08 seismic line, the GANDI 09 total resistivity panel (Figure 6e) represents the effective resistivity distribution and the background resistivity (Figure 6d) represents the gas-free case. In this case, the total resistivity shows higher correlation with the stratigraphy, but scattered high values up to 10 Ω·m and more characterize the whole section. In some zones, as at CDP 2500 and between CDP 6500 and CDP 7500, a mainly vertical distribution of positive anomalies can be recognized. Furthermore, and unlike the STENAP 08 seismic line, the GANDI 09 seismic section does not show any significant difference in the frequency content. This is partly related to the widespread poor amplitude content of the GANDI 09 seismic line, which had to be corrected during the processing, but it is also likely due to a different gas distribution within the two seismic lines. As for the STENAP 08, the resistivity anomalies appear to be related to the frequency-derived seismic attributes, even if with a minor impact.

<figure><img src="image_11.png" /><figcaption>**Figure 6.** STENAP 08 (a) and GANDI 09 (d) background resistivity. STENAP 08 (b) and GANDI 09 (e) total resistivity. STENAP 08 (c) and GANDI 09 (f) gas content. Arcobaleno 1, Arlecchino 1, and Rachele 1 well locations are shown. The black box indicates a zone of no coverage.</figcaption></figure>
## 4.2. Gas Content

Figures 6c and 6f show the estimated gas distribution from resistivity anomalies and highlight different concentration patterns in the two seismic profiles. In particular, the GANDI 09 line has a more diffuse gas presence, while the STENAP 09 seismic line reveals zones of gas accumulation.

Along the STENAP 08 (Figure 6c), values up to 0.3% can be reached in the shallowest (<0.25 s TWTT) part in CDPs 1–2000 and CDPs 5000–6500. Gas content of 0.15%–0.3% characterizes the SE part of the seismic line, up to 1.2 s of depth; it appears mostly horizontally distributed, consistent with the seismostratigraphic setting. Gas is also present along sub-vertical paths affecting the whole sedimentary column, showing values up to 0.2%.

The GANDI 09 seismic line is affected by zones of poor or no seismic coverage and is, furthermore, very poorly balanced. In fact, the amplitude preserving stack section (Figure 9b) appears to be whitened except for a 1550 CDP-wide zone, between CDP 6000 and CDP 7500, where amplitudes are higher and the stratigraphy is better imaged. Overall, the amplitude information of the GANDI seismic stack section does not allow a reliable correlation with log curves. This is the reason why we used the migrated and gained version of the line for correlation and gas content assessment. In doing so, we affected the absolute amplitude information, which should be preserved for a trustworthy gas quantification. Nevertheless, we obtained a picture of the distribution of the gas along the GANDI 09 seismic line. Figure 6f shows a diffuse gas concentration, with no consistent accumulation at specific depths or stratigraphic layers. Note that local increase in gas concentration, up to 0.2%, is present along vertical paths (i.e., around CDP 2200, 6500, 7250, and 8500) which can extend up to 1 s of depth.

# 5. Discussion

## 5.1. Gas Distribution Versus Gas Content

The two seismic lines reveal peculiar characteristics, which are discussed separately below.

The STENAP 08 (Figure 7) is characterized by a heterogeneous gas distribution. In the SE portion of the seismic line, gas appears to be horizontally distributed following the seismostratigraphic setting. The sedimentary succession is constituted by fine turbidites (Malvic, 2016), which likely provide the pore space for accumulation in the coarser-grained sand beds. The highest concentration values, up to 0.3% of the total space, are distributed horizontally in the shallower levels (TWTT < 0.25 s) between CDPs 0–2000 and CDPs 5000–6500. In the rest of the seismic line, several vertical paths are recognizable, possibly related to the main tectonic structures identified in the study area (Brancolini et al., 2019), implying that gas preferably migrates along such sub-vertical paths, as also suggested by Donda et al. (2013); Donda et al. (2015). A 700 CDP-wide, 2s-long window characterized by pull-downs is recognizable in the middle of the seismic line (Figure 7). The amplitude preserving version of the stack section (Figure 9a) shows how this area is dimmed in amplitude. This entire pull-down zone is characterized by a degradation in the stack response. We interpret it as a giant chimney, rooted below the Messinian Unconformity and joining the sub-vertical faults in focusing the gas migration. Other smaller-scale pull-downs are present along the seismic line and are often associated to areas of amplitude loss, which results in a poor coherence if the non-amplitude preserving processing is adopted. They can be interpreted as gas pipes (Cartwright & Santamarina, 2015; Heggland, 2005; Løseth et al., 2009; Maestrelli et al., 2017). Note that, in the Baltic Sea, acoustic turbidity is associated with less than 0.5% of gas concentration (García-García et al., 2007). Bright spots also occur at various depths in the STENAP 08 seismic line; the two most evident are highlighted in Figures 7 and 9a. However, we are not able to provide a good estimate of gas concentration at these high amplitude anomaly zones. This is probably because of their extremely local nature, which cannot be resolved by our analysis. Petrophysical properties are laterally extrapolated and such small and local anomalies cannot be properly resolved. Figure 9a also shows that the shallow reflectors located near the sea bottom, at CDPs 0–2000 and CDPs 5000–6500, have low amplitude. In fact, as already explained, due to post-critical conditions, the whole shallow part (<0.25 s) of the stack section is very poor in reflected energy. Furthermore, an apparent low-frequency effect induced by the applied deconvolution is unavoidable (deconvolution often fails in low amplitude zones). However, small lateral variations in amplitude in these shallow zones can be recognized, suggesting the presence of gas leakage areas just below and at the seabed. Figure 7 shows its highest gas concentration just below these shallow low-amplitude zones.
<figure><img src="image_12.png" /><figcaption>**Figure 7.** (a) STENAP 08 time migrated stack section. Main faults are drawn and gas-related features are highlighted. The Messinian Unconformity is represented by the yellow horizon, whereas the top of Carbonate Succession is the orange seismic reflector. (b) STENAP 08 gas concentration superimposed to the true-amplitude stack section. Well locations are shown.</figcaption></figure>

In the GANDI 09 seismic line (Figure 8), gas is more scattered and is locally distributed along vertical paths, suggesting vertical migration. In some cases (e.g., at CDP 2500), faults seem to act as preferential conduits for gas migration. Figure 8a reports the main faults interpreted by Brancolini et al. (2019) as part of regional tectonic lineaments. However, sometimes gas appears to be dissociated from the faults and its presence causes seismic response degradation. Possibly, in these circumstances, gas is able to migrate up to the surface through local heterogeneities in the sedimentary succession, as commonly occurring in turbidite-dominated depositional environments. Because of lack of reflections, the shallowest tens of ms in the GANDI 09 seismic line show the same loss in amplitude found along some shallow parts of the STENAP 08. In the GANDI 09, however, these sediments do not show any change in amplitude response along the seismic line (Figure 9), suggesting a different shallow gas distribution between the two lines.

The two seismic profiles cross an area where the hydrocarbon industry identified several biogenic gas fields. In Figure 1, we report their location, according to the Fact Book Exploration and Production report ENI published in 2013 (http://report.eni.com/factbook-2013/it/i-settori/exploration-production/i-paesi-di-attivita/italia.html). The GANDI 09 crosses two of the identified gas fields: one in the very SW end of the seismic line and a second one around CDP 5000, in the vicinity of Triglia Mare 1 well. From our results, at the SW of the seismic line, in the vicinity of Arlecchino 1, gas is present along the whole sedimentary column. However, we want to underline one more time that, in treating this seismic line, we could not preserve the amplitude information, so that the absolute gas estimation could be compromised. At CDP 5000, the seismic data show response degradation and reveal a shallow pull-down effect. The STENAP 08 also crosses two gas fields, at CDP 4700 and CDP 7500. In these zones, gas concentration is higher with respect to the surrounding area and distributed along deeply rooted vertical paths. Here, the seismic data show gas-related pull down effects, which we interpreted as gas chimneys.

Not much information has been made publicly available on the gas fields exploited in the study area (Figure 1), but both Arcobaleno 1 and Arlecchino 1 wells, on the NW end of the STENAP 08 seismic line, were drilled to investigate the Plio-Pleistocene sandy series. These sands were found to contain only water and both wells were
<figure><img src="image_13.png" /><figcaption>**Figure 8.** (a) GANDI 09 time migrated stack section. Red arrows indicate CDPs of poor coverage. The black box indicates a zone of no coverage. Main faults are drawn and gas related features are highlighted. The Messinian Unconformity is represented by the yellow horizon, whereas the top of Carbonate Succession is the orange seismic reflector. (b) GANDI gas concentration superimposed to the true-amplitude stack section. Well locations are shown.</figcaption></figure>

classified as sterile (Relazione al comitato tecnico per gli idrocarburi, permesso A.R54.MI). Further to the north, NE of Triglia Mare 1 well, Pleistocene turbidites revealed to be a good reservoir for gas and have been exploited. Moreover, a strong gas-related bright spot was found in the Pleistocene sequence at Triglia Mare 1 well (Relazione tecnica allegata all'istanza di rinuncia del permesso A.R54.MI).

## 5.2. Comparison With Ground Truth Data

The only available ground truth data to test our results are sparse measurements of total-gas content in Arcobalenno 1 and Arlecchino 1 wells and a total-gas log in Rachele 1 well. In the case of Arcobalenno 1 and Arlecchino 1, these are punctual values scattered along the whole wells depth. They are just noted, often by hand, on the side of the other log curves and there is no information about how these values were obtained. Because of the lack of details about the actual gas measurement, we had to do several assumptions. For simplicity, and considering that the vertical resolution of these logs is around 50 cm, we supposed that the gas concentration was measured in 1 m-thick sediment layer. To make a comparison with our results, we needed to spread this amount of gas over a thickness equal to the vertical seismic resolution of the seismic sections, which is changing with depth. In this regard, we divided both wells in three depth intervals: 0–600, 600–1,200, and 1,200–1,600 m, which correspond to the main trend changes in the sonic logs. The two deepest intervals also coincide with the ones reported in the velocity-resistivity cross-plots in Figure 2. For each interval, we observed the dominant seismic frequency *f*. Then, we calculated the theoretical vertical seismic resolution, which is defined as $v/4f$, where *v* is the average seismic velocity at that depth interval, from the sonic log. We found a vertical resolution of 10, 35, and 46 m for each interval of Arcobalenno 1 and of 10, 38, and 50 m for each interval of Arlecchino 1. Finally, we spread the 1 m-gas measurements over the vertical resolution-derived sediment thickness and, using the sonic logs, we converted their depth in TWTT. In this way, these direct measurements were made comparable with our gas content results, which we obtained analyzing totally independent data. Both the original gas values (Gas) and the spread values (Gas<sub>sp</sub>) are reported in Table 6, in depth and in TWTT, for Arcobalenno 1 and Arlecchino 1 wells. We can observe that, Gas<sub>sp</sub> is of the same order of magnitude as our gas content results, which are mainly below 0.3%.
<figure><img src="image_14.png" /><figcaption>**Figure 9.** (a) STENAP 08, (b) GANDI 09 amplitude preserving stack section. Zones of anomalous amplitude are highlighted. Translucent box indicates an area of no coverage. Well locations are shown.</figcaption></figure>

These low values are coherent with the fact that none of these wells has ever been under production, and they are, instead, considered sterile. Obviously, we cannot expect a perfect correlation between the two sets of values. This has several reasons, amongst others the fact that the wells do not exactly lie on the STENAP 08 seismic line. Furthermore, we know that, at least within the Plio-Pleistocene formations, the gas is hosted in some meters-thick “pool” or “lenses” whose shape we totally ignore and that can cause significant 3D effects which cannot be resolved by 2D data. Indeed, out-of-plane events are present in several shot gathers of both seismic lines, possibly testifying for 3D structures.

In the case of Rachele 1 well, an almost continuous total-gas log is present, spanning from ~250 to ~1,100 m. This is not the same depth interval where the sonic log was measured, so that both the conversion in TWTT and the estimation of the vertical seismic resolution cannot be accurately made. However, extrapolating the velocity values at depths where they are missing, we could provide an approximated conversion in TWTT of the gas log. In the seismic data, we observed a unique dominant frequency of 45 Hz at all TWTT and we calculated a vertical resolution of 11 m. However, a proper comparison with our results is prevented by the fact that the true-amplitude version of the seismic data didn't allow a satisfying correlation with the log data. Rachele 1 is almost 10 km away from the GANDI 09 seismic line and, as we already explained, its projection on the line falls in the vicinity of a zone of very poor seismic signal coverage. The total-gas log in depth, its conversion in TWTT and the values spread on 11 m-sediment thickness are plotted in Figure 10.
<table><caption>**Table 6**

Sparse Total-Gas Measurements From the Well-Logs Reports</caption><thead><tr><th colspan="4">Arcobaleno 1</th><th colspan="4">Arlecchino 1</th></tr><tr><th>Depth (m)</th><th>TWTT (ms)</th><th>Gas (%)</th><th>Gas<sub>sp</sub> (%)</th><th>Depth (m)</th><th>TWTT (ms)</th><th>Gas (%)</th><th>Gas<sub>sp</sub> (%)</th></tr></thead><tbody><tr><td>245</td><td>266</td><td>1.5</td><td>0.15</td><td>485</td><td>485</td><td>6</td><td>0.58</td></tr><tr><td>325</td><td>373</td><td>2</td><td>0.2</td><td>720</td><td>663</td><td>7</td><td>0.18</td></tr><tr><td>460</td><td>484</td><td>4.5</td><td>0.45</td><td>850</td><td>763</td><td>7</td><td>0.18</td></tr><tr><td>470</td><td>491</td><td>2.5</td><td>0.25</td><td>900</td><td>822</td><td>7</td><td>0.18</td></tr><tr><td>490</td><td>500</td><td>4.5</td><td>0.45</td><td>1,030</td><td>851</td><td>12</td><td>0.32</td></tr><tr><td>520</td><td>530</td><td>4.5</td><td>0.45</td><td>1,100</td><td>932</td><td>9</td><td>0.24</td></tr><tr><td>605</td><td>666</td><td>2.8</td><td>0.08</td><td>1,220</td><td>935</td><td>10</td><td>0.26</td></tr><tr><td>765</td><td>686</td><td>3.2</td><td>0.09</td><td>1,280</td><td>1,105</td><td>9</td><td>0.18</td></tr><tr><td>870</td><td>834</td><td>3.5</td><td>0.1</td><td></td><td></td><td></td><td></td></tr><tr><td>910</td><td>872</td><td>2.8</td><td>0.08</td><td></td><td></td><td></td><td></td></tr><tr><td>930</td><td>794</td><td>2.5</td><td>0.07</td><td></td><td></td><td></td><td></td></tr><tr><td>1,090</td><td>953</td><td>2.5</td><td>0.07</td><td></td><td></td><td></td><td></td></tr><tr><td>1,150</td><td>933</td><td>3.2</td><td>0.09</td><td></td><td></td><td></td><td></td></tr><tr><td>1,260</td><td>1,042</td><td>4</td><td>0.09</td><td></td><td></td><td></td><td></td></tr><tr><td>1,350</td><td>1,147</td><td>2</td><td>0.04</td><td></td><td></td><td></td><td></td></tr><tr><td>1,445</td><td>1,123</td><td>1</td><td>0.02</td><td></td><td></td><td></td><td></td></tr><tr><td>1,490</td><td>1,063</td><td>0.8</td><td>0.02</td><td></td><td></td><td></td><td></td></tr></tbody></table>

In summary, our gas-content results are consistent with the few direct measurements available in the area. Our estimation also agrees with the coring results obtained by García-García et al. (2007) offshore the Po Delta. Of the 20 sites analyzed, 5 sites yielded anomalous gas concentrations, with more than 90 ppm (equal to 0.009%) of CH₄ and, amongst these, three sites yielded 4·10⁴ ppm (equal to 4%; García-García et al., 2007). Further south of the study area, the Barbara Nord 007 well is gas-charged and is part of the well-known ENI Barbara gas field. This field includes tens of wells under production and the well-logs from Barbara Nord 007 are the only ones made available by ENI (https://www.videpi.com/videpi/videpi.asp). Here, gas was found in levels of ~2–30 m thickness. The maximum registered value of total-gas was 1.2% in a ~30 m-thick sand layer, which is equivalent to ~0.04% for a 1 m-thick layer.

### 5.3. Environmental Forcing in the Gas Distribution

The gas distribution appears to be strongly correlated with the seismostratigraphy and tectonic setting of the study area. In both seismic lines, the major faults act as preferential paths for gas migration. However, gas conduits could also originate from other mechanisms, such as hydrofracturing. Overpressure conditions may have originated during Pleistocene times, when sedimentation rates were high (Amadori et al., 2018), as a consequence of rapid burial. Opening of conduits would have occurred during abrupt changes in the vertical stress field due to big sea level fluctuations. This is the case of the Messinian salinity crisis, which has been recognized as a triggering mechanism of fluid expulsion by several authors (Bertoni & Cartwright, 2015; Dale et al., 2021; Iadanza et al., 2015). In the Northern Adriatic Sea, ~800 m of the Messinian sea-level drop was estimated by several authors (Amadori

et al., 2018; Ghielmi et al., 2013; Mancin et al., 2016), suggesting a relevant change in the vertical stress field during this widespread event. Therefore, we suggest that paleobathymetric variations could have had a role in the

<figure><img src="image_15.png" /><figcaption>**Figure 10.** Total-gas measurements for Rachele 1 well. (a) Values in depth; (b) values in TWTT (after conversion); and (c) values in TWTT, spread over 10 m-thick sediments.</figcaption></figure>
opening of gas migration conduits, which possibly formed along fault zones. Moreover, the deepest gas chimneys are rooted in the Carbonate Succession, suggesting a further deep gas source. The gas migrates along the identified faults throughout the Plio-Pleistocene sedimentary sequence, up to the uppermost stratigraphic levels, where it possibly mixes with the shallow gas identified in the high resolution geophysical data (Donda et al., 2019).

## 6. Conclusions

Our results demonstrate that, despite a proper P-wave velocity field is not available from seismic data analysis, the joint use of seismic and well data can be used to assess gas saturation by using the estimated resistivity anomalies. This has a potential yet very practical outcome: in very shallow environments, where the water column is limited and water masses well mixed, the potential transfer of methane-predominant gas from sediment into the atmosphere could be significant and introduces a non-negligible doubt on the global budget of atmospheric methane (Sultan, 2020, and references therein). It has been recently proven that even a moderate sea-level rise (<1 m) may significantly impact gas emissions at deep-water depths and partially counterbalance future temperature effects on the global marine gas systems (Sultan, 2020, and references therein). In fact, while the emissions from seeps should be considered as natural sources in the global CH$_4$ budget, further warming of surface waters could increase CH$_4$ emissions and provide a positive feedback on warming climate (Borges et al., 2016). This feedback will be expected to be acute in shallow gassy areas such as those in the Northern Adriatic Sea, where the increase of temperature will stimulate the biogenic CH$_4$ production, as well as promote gas leakage from sediments. A proper, comprehensive assessment of the distribution, quantification and migration pathways of methane-predominant gas is thus fundamental to help forecasting the amount of methane that could leak from the sedimentary successions in shallow marine environments and eventually reach the atmosphere.

## Data Availability Statement

Data are stored in the repository “Seismic data Network Access Point (SNAP)” (https://snap.ogs.trieste.it/cache/index.jsp). Well-log charts and technical reports from the oil and gas exploration have been made available by the “Visibility of petroleum exploration data in Italy (ViDEPI)” project.

## Acknowledgments

The authors would like to acknowledge Prof. Angelo Camerlenghi for invaluable discussions during the data interpretation and Dr. Flavio Accaino for the insights on seismic processing. The authors thank IHS Kingdom Educational Grant Program for providing OGS with a free license for the use of KingdomSuite TKS complete package interpretation software. The authors are grateful to IFP-Energie Nouvelles for the use of the EasyTrace tool, which has been particularly useful in the log analysis. The authors thank Emerson EandP for a free license of Paradigm Epos, used to process the seismic data and CGG for a free license of Emerge, used to perform the seismic-log correlation. The seismic data analyzed in this work were collected by the Istituto Nazionale di Oceanografia e di Geofisica Sperimentale (OGS) with the R/V OGS Explora in the framework of the “Seismostratigraphic and Tectonic Evolution of the Northern Adriatic Sea in the Plio-Quaternary (STENAP)” and the “Gas Emissions in the Northern Adriatic Sea (GANDI)” projects funded by the OGS during which more than 1,200 km of 2D multi-channel seismic lines, sub-bottom CHIRP profiles and morpho-bathymetric data were collected. TRETAMARA “Trezze, Tegnue and marine environments of the Northern Adriatic: management proposals” Interreg Italy Slovenia project (2020–2022).

## References

* Adelinet, M., Fortin, J., & Gueguen, Y. (2011). Dispersion of elastic moduli in a porous-cracked rock: Theoretical predictions for squirt flow. *Tectonophysics*, **503**(1–2). https://doi.org/10.1016/j.tecto.2010.10.012

* Adelinet, M., Fortin, J., Guéguen, Y., Schubnel, A., & Geoffroy, L. (2010). Frequency and fluid effects on elastic properties of basalt: Experimental investigations. *Geophysical Research Letters*, **37**. https://doi.org/10.1029/2009GL041660

* Adelinet, M., & Le Ravalec, M. (2015). Effective medium modeling: How to efficiently infer porosity from seismic data? *Interpretation*, **3**(4). https://doi.org/10.1190/int-2015-0065.1

* Amadori, C., Garcia-Castellanos, D., Toscani, G., Sternai, P., Fantoni, R., Ghielmi, M., & Di Giulio, A. (2018). Restored topography of the Po Plain-Northern Adriatic region during the Messinian base-level drop—Implications for the physiography and compartmentalization of the palaeo-Mediterranean basin. *Basin Research*, **30**(6). https://doi.org/10.1111/bre.12302

* Archie, G. E. (1942). The electrical resistivity log as an aid in determining some reservoir characteristics. *Transactions of the AIME*, **146**, 54–62. https://doi.org/10.2118/942054-G

* Bailly, C., Fortin, J., Adelinet, M., & Hamon, Y. (2019). Upscaling of elastic properties in carbonates: A modeling approach based on a multiscale geophysical data set. *Journal of Geophysical Research: Solid Earth*, **124**, 13021–13038. https://doi.org/10.1029/2019JB018391

* Bertello, F., Fantoni, R., & Franciosi, R. (2008). Overview on Italy’s petroleum systems and related oil and gas occurrences. In *70th EAGE conference & exhibition, Rome, Italy, extended abstract & exhibitor’s catalogue*.

* Bertoni, C., & Cartwright, J. (2015). Messinian evaporites and fluid flow. *Marine and Petroleum Geology*, **66**, 165–176. (The Messinian events and hydrocarbon exploration in the Mediterranean). https://doi.org/10.1016/j.marpetgeo.2015.02.003

* Borges, A., Champenois, W., Gypens, N., Delille, B., & Harlay, J. (2016). Massive marine methane emissions from near-shore shallow coastal areas. *Scientific Reports*, **6**, 27908. https://doi.org/10.1038/srep27908

* Brancolini, G., Civile, D., Donda, F., Tosi, L., Zecchin, M., Volpi, V., et al. (2019). New insights on the Adria plate geodynamics from the northern Adriatic perspective. *Marine and Petroleum Geology*, **109**, 687–697. https://doi.org/10.1016/j.marpetgeo.2019.06.049

* Cartwright, J., Huuse, M., & Aplin, A. (2007). Seal bypass systems. *AAPG Bulletin*, **91**. https://doi.org/10.1306/04090705181

* Cartwright, J., & Santamarina, C. (2015). Seismic characteristics of fluid escape pipes in sedimentary basins: Implications for pipe genesis. *Marine and Petroleum Geology*, **65**, 126–140. https://doi.org/10.1016/j.marpetgeo.2015.03.023

* Casero, P. (2004). Structural setting of petroleum exploration plays in Italy. In *Special Volume of the Italian Geological Society for the IGC 32 Florence-2004*.

* Casero, P., & Bigi, S. (2013). Structural setting of the Adriatic basin and the main related petroleum exploration plays. *Marine and Petroleum Geology*, **42**, 135–147. https://doi.org/10.1016/j.marpetgeo.2012.07.006

* Cazzini, F., Dal Zotto, O., Fantoni, R., Ghilemi, M., Ronchi, P., & Scotti, P. (2015). Oil and gas in the Adriatic foreland, Italy. *Journal of Petroleum Geology*, **38**(3), 255–279. https://doi.org/10.1111/jpg.12610
Coren, F., Volpi, V., & Tinivella, U. (2001). Gas hydrate physical properties imaging by multi-attribute analysis—Blake Ridge BSR case history. *Marine Geology*, **178**(1–4). https://doi.org/10.1016/s0025-3227(01)00156-6

Correggiari, A., Roveri, M., & Trincardi, F. (1996). Late Pleistocene and Holocene evolution on the North Adriatic Sea. *Italian Journal of Quaternary Sciences*, **9**, 697–704.

Dale, M., Marín-Moreno, H., Falcon-Suarez, I., Grattoni, C., Bull, J., & McNeill, L. (2021). The Messinian salinity crisis as a trigger for high pore pressure development in the Western Mediterranean. *Basin Research*. https://doi.org/10.1111/bre.12554

Donda, F., Civile, D., Forlin, E., Volpi, V., Zecchin, M., Merson, B., & De Santis, L. (2013). The northernmost Adriatic Sea: A potential location for CO₂ geological storage? *Marine and Petroleum Geology*, **42**, 148–159. https://doi.org/10.1016/j.marpetgeo.2012.10.006

Donda, F., Forlin, E., Gordini, E., Panieri, G., Buenz, S., Volpi, V., et al. (2015). Deep-sourced gas seepage and methane-derived carbonates in the Northern Adriatic Sea. *Basin Research*, **27**, 531–545. https://doi.org/10.1111/bre.12087

Donda, F., Tinivella, U., Gordini, E., Panieri, G., Volpi, V., Civile, D., et al. (2019). The origin of gas seeps in the Northern Adriatic Sea. *Italian Journal of Geoscience*, **138**(2). https://doi.org/10.3301/IJG.2018.34

Eshelby, J. D. (1957). The determination of the elastic field of an ellipsoidal inclusion, and the related problems. *Proceedings of the Royal Society of London. Series A. Mathematical and Physical Sciences*, **241**(1226), 376–396.

Etiope, G., & Ciccioli, P. (2009). Earth's degassing: A missing ethane and propane source. *Science (New York, N.Y.)*, **323**, 478. https://doi.org/10.1126/science.1165904

Etiope, G., Ciotoli, G., Schwietzke, S., & Schoell, M. (2019). Gridded maps of geological methane emissions and their isotopic signature. *Earth System Science Data*, **11**, 1–22. https://doi.org/10.5194/essd-11-1-2019

García-García, A., Orange, D., Miserocchi, S., Correggiari, A., Langone, L., Lorenson, T., et al. (2007). What controls the distribution of shallow gas in the Western Adriatic Sea? *Continental Shelf Research*, **27**, 359–374. https://doi.org/10.1016/j.csr.2006.11.003

Ghielmi, M., Minervini, M., Nini, C., Rogledi, S., & Rossi, M. (2013). Late Miocene-Middle Pleistocene sequences in the Po Plain – Northern Adriatic Sea (Italy): The stratigraphic record of modification phases affecting a complex foreland basin. *Marine and Petroleum Geology*, **42**. https://doi.org/10.1016/j.marpetgeo.2012.11.007

Giustiniani, M., Tinivella, U., Parolai, S., Donda, F., Brancolini, G., & Volpi, V. (2020). Integrated geophysical analyses of shallow-water seismic imaging with Scholte wave inversion: The Northern Adriatic Sea case study. *Frontiers of Earth Science*, **8**, 532. https://doi.org/10.3389/feart.2020.587898

Glover, P. W. J. (2016). Archie's law – A reappraisal. *Solid Earth*, **7**(4), 1157–1169. https://doi.org/10.5194/se-7-1157-2016

Gordini, E. (2009). *Integrazione di metodologie geofisiche, geomorfologiche, sedimentologiche e geochimiche, per la definizione della genesi e dell'eta' degli affioramenti rocciosi presenti sul fondale marino dell'Adriatico settentrionale* (Unpublished doctoral dissertation). Universita' di Trieste.

Gordini, E., Falace, A., Kaleb, S., Donda, F., Marocco, R., & Tunis, G. (2012). *Methane-related carbonate cementation of marine sediments and related macroalgal coralligenous assemblages in the Northern Adriatic Sea in seafloor geomorphology as benthic habitat-GeoHAB Atlas of seafloor geomorphic features and benthic habitats* (1st ed.). Elsevier.

Gordini, E., Marocco, R., & Vio, E. (2002). Stratigrafia del sottosuolo della “trezza grande” (golfo di trieste, adriatico settentrionale). *Gortania*, **24**, 31–63.

Greenberg, M. L., & Castagna, J. P. (1992). Shear-wave velocity estimation in porous rocks: Theoretical formulation, preliminary verification and application. *Geophysical Prospecting*, **40**(2). https://doi.org/10.1111/j.1365-2478.1992.tb00371.x

Hampson, D., Todorov, T., & Russell, B. (2000). Using multi-attribute transforms to predict log properties from seismic data. *Exploration Geophysics*, **31**. https://doi.org/10.1071/EG00481

Hartmann, D. L., Klein Tank, A. M. G., Rusticucci, M., Alexander, L., Brönnimann, S., Charabi, Y. A. R., et al. (2013). Observations: Atmosphere and surface. In *Climate change 2013 the physical science basis* (Vol. 9781107057999, pp. 159–254). Cambridge University Press. https://doi.org/10.1017/CBO9781107415324.008

Heggland, R. (2005). Using gas chimneys in seal integrity analysis: A discussion based on case histories. *AAPG Hedberg Series*, no. 2, 237–245.

Iadanza, A., Sampalmieri, G., & Cipollari, P. (2015). Deep-seated hydrocarbons in the seep “Brecciated Limestones” of the Maiella area (Adriatic foreland basin): Evaporitic sealing and oil re-mobilization effects linked to the drawdown of the Messinian Salinity Crisis. *Marine and Petroleum Geology*, **66**, https://doi.org/10.1016/j.marpetgeo.2015.03.006

Jiang, Y., Zhou, J., Fu, X., Cui, L., Fang, C., & Cui, J. (2021). Analyzing the origin of low resistivity in gas-bearing tight sandstone reservoir. *Geofluids*. https://doi.org/10.1155/2021/4341804

Johansson, S., Rosqvist, H., Svensson, M., Dahlin, T., & Leroux, V. (2011). An alternative methodology for the analysis of electrical resistivity data from a soil gas study. *Geophysical Journal International*, **186**(2), 632–640. https://doi.org/10.1111/j.1365-246X.2011.05080.x

Kannberg, P. K., & Constable, S. (2020). Characterization and quantification of gas hydrates in the California borderlands. *Geophysical Research Letters*, **47**(6), e2019GL084703. https://doi.org/10.1029/2019GL084703

Lee, M. (2004). Elastic velocities of partially gas-saturated unconsolidated sediments. *Marine and Petroleum Geology*, **21**, 641–650. https://doi.org/10.1016/j.marpetgeo.2003.12.004

Le Ravalec, M., Doligez, B., & Lerat, O. (2014). *Integrated reservoir characterization and modeling*. IFPEN. https://books.ifpenergiesnouvelles.fr/ebooks/integrated_reservoir_characterization_and_modeling/

Loseth, H., Gading, M., & Wensaas, L. (2009). Hydrocarbon leakage interpreted on seismic data. *Marine and Petroleum Geology*, **26**, 1304–1319. https://doi.org/10.1016/j.marpetgeo.2008.09.008

Maestrelli, D., Iacopini, D., Jihad, A. A., Bond, C. E., & Bonini, M. (2017). Seismic and structural characterization of fluid escape pipes using 3D and partial stack seismic from the Loyal Field (Scotland, UK): A multiphase and repeated intrusive mechanism. *Marine and Petroleum Geology*, **88**, https://doi.org/10.1016/j.marpetgeo.2017.08.016

Malvic, T. (2016). Regional turbidites and turbiditic environments developed during Neogene and Quaternary in Croatia. *Materials and Geoenvironment*, **63**, 39–54.

Mancin, N., Barbieri, C., Di Giulio, A., Fantoni, R., Marchesini, A., Toscani, G., & Zanferrari, A. (2016). The Friulian-Venetian Basin II: Paleogeographic evolution and subsidence analysis from micropaleontological constraints. *Italian Journal of Geosciences*, **135**. https://doi.org/10.3301/ijg.2015.34

Mitnick, E., Lammers, L., Zhang, S., Zaretskiy, Y., & DePaolo, D. (2018). Authigenic carbonate formation rates in marine sediments and implications for the marine $\delta^{13}$C record. *Earth and Planetary Science Letters*, **495**, 135–145. https://doi.org/10.1016/j.epsl.2018.05.018

Mondol, N. H., Jahren, J., Bjorlykke, K., & Brevik, I. (2008). Elastic properties of clay minerals. *The Leading Edge*, **27**(6). https://doi.org/10.1190/1.2944161

Prather, M., Holmes, C., & Hsu, J. (2012). Reactive greenhouse gas scenarios: Systematic exploration of uncertainties and the role of atmospheric chemistry. *Geophysical Research Letters*, **39**, 9803. https://doi.org/10.1029/2012GL051440
Saunois, M., Stavert, A., Poulter, B., Bousquet, P., & Canadell, J. E. A. (2020). The global methane budget 2000–2017. *Earth System Science Data, 12*, 1561–1623. https://doi.org/10.5194/essd-12-1561-2020

Schiermeier, Q. (2020). Global methane levels soar to record high. *Nature*. https://doi.org/10.1038/d41586-020-02116-8

Schon, J. (2015). *Physical properties of rocks – Fundamentals and principles of petrophysics* (2nd ed.), Elsevier. https://www.elsevier.com/books/physical-properties-of-rocks/schon/978-0-08-100404-3

Schwietzke, S., Sherwood, O., Bruhwiler, L., Miller, J., Etiope, G., Dlugokencky, E., et al. (2016). Upward revision of global fossil fuel methane emissions based on isotope database. *Nature*, 538, 88–91. https://doi.org/10.1038/nature19797

Senger, K., Birchall, T., Betlem, P., Ogata, K., Ohm, S., Olaussen, S., & Paulsen, R. S. (2021). Resistivity of reservoir sandstones and organic rich shales on the Barents Shelf: Implications for interpreting CSEM data. *Geoscience Frontiers*, 12(6), 101063. https://doi.org/10.1016/j.gsf.2020.08.007

Sultan, N., Plaza-Faverola, A., Vadakkepuliyambatta, S., Buenz, S., & Knies, J. (2020). Impact of tides and sea-level on deep-sea Arctic methane emissions. *Nature Communications*, 11(1), https://doi.org/10.1038/s41467-020-18899-3

Thornton, B., Prytherch, J., Andersson, K., Brooks, I., Salisbury, D., Tjernström, M., & Crill, P. (2020). Shipborne eddy covariance observations of methane fluxes constrain Arctic sea emissions. *Science Advances*, 6, eaay7934. https://doi.org/10.1126/sciadv.aay7934

Tinivella, U., Giustiniani, M., & Vargas-Cordero, I. (2017). Wave equation datuming applied to seismic data in shallow water environment and post-critical water bottom reflection. *Energies*, 10(9). https://doi.org/10.3390/en10091414

Vernik, L., & Kachanov, M. (2010). Modeling elastic properties of siliciclastic rocks. *Geophysics*, 75(6). https://doi.org/10.1190/1.3494031

Verschuur, D. J. (2006). *Seismic multiple removal techniques; past, present en future* (Vol. Education Tour Series), 2006 EAGE Education Tour. EAGE. https://www.eageseg.org/en/2006/10/seismic-multiple-removal-techniques-past-present-and-future/

Wang, Z. (2001). Fundamentals of seismic rock physics. *Geophysics*, 66(2). https://doi.org/10.1190/1.1444931

Weitemeyer, K. A., Constable, S., & Tréhu, A. M. (2011). A marine electromagnetic survey to detect gas hydrate at Hydrate Ridge, Oregon. *Geophysical Journal International*, 187(1), 45–62. https://doi.org/10.1111/j.1365-246X.2011.05105.x

Yilmaz, O. (2001). *Seismic data analysis* (2nd ed.). SEG. https://library.seg.org/doi/10.1190/1.9781560801580

Zecchin, M., Tosi, L., & Caffau, M. (2011). Relationship between peat bed formation and climate changes during the last glacial in the Venice area. *Sedimentary Geology*, 238. https://doi.org/10.1016/j.sedgeo.2011.04.011