Abstract
Pipe structures are considered as fluid conduits beneath cold seeps. These structures have been observed in many geological settings and are widely accepted as the most critical pathway for fluid migration. One of such pipe structures in the Haima cold seep region is investigated herein. The pipe structure extends from below the BSR and reaches the seafloor. It is characterized by a string of events with short and strong seismic amplitudes, similar to the string of bead reflections (SBRs) associated with small-scale caves in carbonate reservoirs. This leads to the hypothesis that multiple small-scale bodies exist within the pipe structure. We test this hypothesis by analysis of diffraction waves and numerical seismic modeling. Travel time pattern analysis indicates that the diffractors within the pipe structure caused the rich diffraction waves on the shot records, and the reversed polarity indicates that the diffractors have a lower impedance than the surrounding sediments. These low-impedance bodies are interpreted as gas pockets within the pipe structures. Based on these interpretations, a conceptual model is proposed to describe the fluid migration process within the pipe. Briefly, we propose that gas pockets within the pipe structure could be analogue to the magma chambers located beneath volcanoes and this may provide a new insight into how gases migrate through the pipe structure and reach the seafloor.
1. Introduction
Cold seeps are areas on the seafloor where cold hydrocarbon-rich fluid leaks into the ocean water from the seafloor sediments. Cold seeps are widespread globally [1] and have attracted considerable attention because of their potential role in the global methane budget [2, 3] and their association with gas hydrates and seafloor ecosystems [4, 5]. The occurrence of cold seeps can be inferred from the interpretation of acoustic gas flares on high-frequency acoustic data from the water column and from the presence of chemosynthetic biotics and authigenic carbonate. Such gas flares have been observed in most of the oceans around the world, such as the Hydrate Ridge [6], Barbados [7], Costa Rica [8], Hikurangi Margin [9], Gulf of Mexico [10], Mediterranean [11], Black Sea [12], U.S. Atlantic margin [13, 14], Sub-Antarctic island [15], and the Svalbard continental margin [16]. Chemosynthetic biotics and authigenic carbonate are also widespread around the world and are commonly interpreted as local bright amplitudes in the seafloor [5].
Cold seeps are often associated with a pipe structure located beneath the cold seep site. The pipe structure is considered the fluid migration conduit feeding the cold seep. Such fluid escape pipes have been discovered around the world mainly relying on the interpretation of seismic reflection data [17]. However, the current understanding on the role of the pipe structure in the cold seep is still limited. Therefore, it is important to investigate the geometry and the internal structure of the pipe structure, since the nature of the pipe might determine the activity and life cycle of the overlaying cold seep. Several geophysical methods have been used to reach this objective. High-frequency subbottom profilers provide a very detailed image of the near-seafloor pipe structure [18]. However, the high-frequency nature of the data implies a short penetration of the signal in the subsurface. High-resolution seismic data have been successfully used to study the pipe structure in the Nyyagg pockmark [19]. Multichannel seismic represents the most common data to study the pipe structure. The existing knowledge about this type of pipe structure is based on the interpretation of 2D and 3D seismic reflection data [17, 20–22]. Multichannel seismic method can provide not only the image of internal structure but also the image of the water column above, which could help to detect the fluid process associated with the pipe structure. For instance, according to numerical studies, seismic waves could be influenced by gas plumes and it is possible to detect gas plumes by using seismic waves [23, 24]. Recently, the controlled source electromagnetic (CSEM) method has been used to study the internal structure in several sites around the world, such as at the Vestnesa Ridge [25], in the western Svalbard continental slope [26], and off Norway [27]. CSEM has proved to be a very promising method to characterize these subsurface structures. However, and despite all efforts, the internal structure of the pipe structure is still not well understood, partially due to the small width of the conduit compared to the spatial resolution of most geophysical methods.
Seismic diffraction is a type of acoustic wave produced by the scattering of a seismic wave after it encounters a discontinuity. Diffraction often appears as hyperbolic, or umbrella-shaped, events on a seismic profile. Unlike the reflection waves associated with smooth variations of the subsurface properties, diffraction waves are generated from small-scale bodies and discontinuities, such as cracks, caves, fractures, and pinch-outs [28–30]. Modelling and interpreting seismic diffractions have proven their value in providing information about the nature and presence of the small-scale bodies in the subsurface [31–34].
In this work, one of such fluid pipe structures within the Haima cold seeps [35] is investigated by seismic diffraction analysis. First, we analyze the diffraction waves on the field shot records. This interpretation was complemented with numerical seismic modeling considering two different scenarios to investigate the nature of the internal structure of the pipe. Then, we infer and propose a genetic model for the fluid migration process occurring within the pipe.
Our work is mostly motivated by the presence of rich diffractions on the shot records and the presence of a series of very short strong and bright reflections on the migration profile, very similar to what has been observed in the carbonate reservoir [36, 37]. In the latter case, caves filled with low-density materials originate a series of short reflections on the seismic data, commonly designated as a string of bead reflections (SBRs) [36]. In this geological setting, numerical and physical modeling combined with the interpretation of seismic data has shown that diffractions are widespread and are the reason for the observed SBRs [38, 39]. Based on the similarity between recorded and modelled seismic signatures, we propose that the pipe structure within our study area is filled in with similar small-scale bodies.
2. Geological Background and Previous Works
Haima cold seeps were discovered in the Qiongdongnan Basin in 2016 [35]. Four methane acoustic flares were observed in the water column based on multibeam water column data acquired in 2016 [40]. Among them, plume B is located directly above the pipe structure. It has a lateral scale of about 30 m and a vertical scale of about 630 m. Though theoretically this plume may be detected by the seismic method [23, 24], we do not observe such an image of plume from a previous work [41]. This could be attributed to the different acquisition times of MB data (in 2016) and seismic data (in 2007) when the intermittent of the gas venting is considered.
Massive hydrates located between 4 and 8 m deep below the seafloor were recovered at the ROV1 and ROV2 sites (Figure 1(b)) [41]. The oxygen isotopic results indicate probable destabilization of gas hydrates in the past, and the carbon isotopic ages suggest a major episodic of carbonate precipitation between 6.1 Ka and 5.1 Ka BP [40]. Moreover, dead bivalves are common at both ROV1 and ROV2 sites, which suggests a decline in seepage activity over time. The variation of seepage intensities is also indicated by lipid biomarker analysis [42] and the compositions of sediments [43]. However, the discovery of gas bubble plumes nearby indicates that methane seepage is still very active in the region [40, 41]. Pore fluid analysis showed lateral migration of methane rich from the ROV1 site to the nearby site which contributes to the enhanced methane flux at the nearby site [44]. The lateral migration of fluid may result from the shallow gas hydrates that might have clogged the fluid channel in the ROV1 site. Chemical and structural analyses of gas hydrates indicate that gases are mixtures of biogenic and thermogenic gases [45]. Although extensive studies have been conducted in this area, the fluid migration process is still not fully understood. Seismic imaging of the subsurface structure of ROV1 and ROV2 showed that magmatic activities have contributed to the formation of the fluid flow and minor faulting may act as fluid conduits [46].

(a) Study area

(b) Seismic acquisition geometry
Yang et al. show seismic profiles passing through three of the plumes and revealed the pipe structures beneath the plumes [41]. These pipes presumably act as the gas migration pathway. However, the seismic expression of pipes is not like a typical pipe structure [17]. The pipe structures identified in this region are characterized by a series of very short strong reflections (Figures 2(a) and 2(b)). These observations resemble SBRs (Figure 2(c)) that have been observed in carbonate formations where caves filled with low-density materials generate this kind of seismic signature [39].

(a) Seismic profile passing through the pipe structure

(b) Zoomed image of the pipe

(c) String of bead reflections
3. Data and Methods
We use a two-dimensional multichannel seismic reflection profile passing through plume B, named 1236, to study the pipe structure below the cold seep location. The data were acquired by the Guangzhou Marine Geological Survey (GMGS) in 2007 aiming at investigating the nature and distribution of gas hydrates. The source used during the acquisition is a GI gun array with a total volume of 160 cu. in. The shot spacing is 25 m. A streamer of 240 channels spaced by 12.5 m was used to record the seismic data. A schematic representation of the seismic acquisition geometry is displayed in Figure 1(b). The dominate frequency of the data is about 60 Hz. The data used in this study have been previously processed using a Kirchhoff prestack time migration method [41]. The resulting migrated seismic profile showed the characteristics of a string of strong and short amplitudes associated with the pipe.
We use the final migrated profile as a reference to display the pipe structure and to determine some distinct layer boundaries (Figure 2(a)). It could be observed that the pipe extends from below the BSR and reaches the seafloor (Figure 2(b)). Several boundaries including the seafloor, H1 and the BSR are labeled. SBRs are distributed in three zones designated as zone A, zone B, and zone C (Figure 2(b)). We performed velocity analysis to the common middle point (CMP) gather close to the pipe structure. Velocity analysis is focused on three layers defined by the horizons of the seafloor, H1, and the BSR (Figure 3).

(a)

(b)

(c)
3.1. Diffraction Analysis
We analyzed the characteristics of diffraction waves on shot records from the perspectives of the travel time pattern and waveform polarity. To analyze the waveform, we mainly focus on the polarity of the diffraction. To analyze the travel time pattern, we picked the diffraction events on the real-shot record and plotted the travel times from different shots together. To extract more information from the travel time patterns, we plotted the travel times against the channels and the receiver coordinates. When travel times from different shot points are plotted against the channel numbers, by plotting the data in different domains, we gain insights into how the diffraction event changes as the shot points move away from the diffractor. When travel times from different shot points are plotted against the receiver coordinates, the locations of the diffraction apex could be determined and this helps us determine the possible subsurface source of the diffractions.
3.2. Seismic Modeling
Seismic modeling is a powerful technique that can be used to help seismic interpretation. In summary, seismic modeling means creating synthetic records for a given acquisition geometry and velocity model. We considered two different scenarios with different degrees of complexity.
The first scenario considers a relatively simple velocity model, which is aimed at studying the travel time patterns associated with the diffractor. The model and the source-receiver configuration are displayed in Figure 4. The velocity is set constant and equal to 1500 m/s above the seafloor and 1600 m/s below the seafloor. The source configuration and receiver configuration are the same as those in the real-data acquisition. Travel times were calculated for six shot points. Travel times are analytically calculated by dividing the distance by the velocity.

The second scenario considers a more realistic configuration. The complex model is aimed at (1) studying the waveform polarity characteristics, (2) investigating the effects of the size of diffractors, and (3) demonstrating how the diffraction waves change over the distance between the diffractor and the source. This scenario comprises alternative models with the presence of small-scale bodies at different depths (Figure 3). The background model has a width of 7 km and a depth of 2 km (Figure 5(a)). Models with one diffractor (Figure 5(b)) are used to study the waveform polarity characteristics and to investigate the effects of the size of diffractors. The nature of and the size of the diffractor can vary. Models with multiple diffractors (Figure 5(c)) are used to demonstrate how the diffraction waves change over the distance between the diffractor and the source.

(a)

(b)

(c)
To model the arrival times, we used a finite difference method [47] to numerically solve the seismic wave equation and create synthetic records. The source configuration and receiver configuration used in the modeling match those used in the field acquisition. To simulate the characteristics of the real data as closely as possible, we used a wavelet with a dominant frequency of 60 Hz which is the same as that of the real data. The grid size is 1.25 m in both the vertical and lateral directions.
4. Results
4.1. Characteristics of Diffracted Waves on the Real Shots
Diffractions are mostly observed between shot records 1760 and 1830 (Figure 6). Shot record 1750 does not show any diffraction, and shot 1830 shows very weak diffractions with the tendency to move out of the recording range. Some short diffractions can be observed between the H1 reflection and the BSR reflection on shot 1760 and shot 1770, such as those labeled D1 and D2 (Figure 6(c)). Though diffractions are abundant, most of them are difficult to map on successive shots, with an exception of the diffraction labeled D3. As the dashed lines show, the diffraction D3 is recorded at later times as the source moves away from the pipe structure.

Shot points 1754 and 1755 are the closest to the pipe structure (Figure 1). These shot records are displayed in Figure 7 and interpreted in detail. On the four shot records, two diffractions labeled D1 and D2 can be observed and traced on successive records.

The travel times of diffraction D3 have been picked for shot records between shots 1780 and 1830 with an interval of 10 shots, and the travel times are plotted in two domains: (1) two-way travel time versus the channel number (Figure 8(a)) and (2) two-way travel time versus the receiver coordinates (Figure 8(b)). As expected, as the source moves away from the pipe structure, the travel times of the diffraction increase (Figure 8(a)). However, the locations of the apex of the travel time are fixed at the site of the pipe structure, independent of the source position (Figure 8(b)).

(a)

(b)
The shot records shown in Figures 6 and 7 show that the seafloor reflection is characterized by blue-red-blue. Some of diffractions mentioned above are definitely characterized by red-blue-red, such as the diffraction D1 and D3 (Figure 7). Figure 9 shows a zoom of shot record 1810 where the diffraction can be easily interpreted as red-blue-red. Therefore, the diffractions have reversed polarity with respect to the seafloor reflections.

4.2. Seismic Modeling
We start by modelling the first simple scenario (Figure 4) by calculating the diffraction travel times for six shot records. The modeled diffraction travel times are plotted in the same domains as the recorded ones (Figure 10). In general, the modeled diffraction exhibits a travel time pattern similar to that of the real diffraction on the field records. Both for real data and modelled data, as the shot points go further away from the pipe, the diffraction travel times increase (Figures 8(a) and 10(a)). In addition, the apexes of the travel time curves are fixed at the site of the pipe structure (Figures 8(b) and 10(b)).

(a)

(b)
The second scenario considers a more complex model. To study the nature of the diffractors, different velocities are assigned to the diffractor in Figure 5(b). Figure 11 shows the shot records corresponding to different velocity anomalies. For the diffractor with a higher velocity than the surrounding sediments, the diffraction has the same polarity with respect to the seafloor reflection (Figure 11(a)), while for the diffractor with a lower velocity than the surrounding sediments, the diffraction has a reversed polarity with respect to the seafloor reflection (Figure 11(b)).

(a)

(b)
To demonstrate the effects of the size of the diffractor on the wave characters, diffractors of different heights are considered. Figure 12 shows shot records corresponding to different heights. When the height is 10 m, the diffraction is characterized by two very close events (Figure 12(a)), while as the height increases, several events can be distinguished (Figures 12(b)–12(d)).

(a)

(b)

(c)

(d)
A velocity model with multiple diffractors (Figure 4(c)) is also considered. For this model, two shot records (SP1810 and SP1768) are created and the shot records are displayed in Figure 13. SP1768 is closer to the diffractors. For SP1768, the apex of the diffraction associated with the shallow diffractor between H1 and BSR occurs between the H1 reflection and BSR reflection (Figure 13(a)), while for SP1810 which is further to the diffractors, the apex of the diffraction associated with the shallow diffractor between H1 and BSR goes below the BSR reflection (Figure 13(b)).

(a)

(b)
5. Discussion
5.1. Small-Scale Bodies as the Source of Diffraction
For the marine seismic survey, possible sources of diffraction include the following origins: (1) other nearby vessels, (2) in-line rugged seafloor morphology, (3) out of the plane source, and (4) subseafloor small-scale geological bodies. (1)The vessel itself does not transmit signals into the water, but the propeller can generate acoustic signals. These signals can be recorded by streamers and be a possible source of diffractions on the shot record. This type of diffraction exists from top to bottom on the field records, since the propeller continuously generates signal. Moreover, this type of diffraction is generally more prominent at arrival times before the seabed reflection. We did not observe either characteristics on our shot records. Therefore, it is unlikely that diffractions are from the nearby ship(2)From the migrated seismic profile (Figure 2), we cannot identify any in-line seafloor morphology. Previous studies have also shown that the seafloor is flat and smooth [40]. Therefore, the in-line irregular seafloor body is not the source(3)Seafloor mounds, canyon walls, slab type, and mud volcanoes that are not located in line with the seismic profiles may act as the source of diffraction. But previous studies have also shown that the seafloor is flat and smooth with no such type of structures observed [41]. Furthermore, these out of plane signals could produce a side sweep on the stack profile. However, such a type of side sweep is not observed (Figure 2). Therefore, the out-of-plane source is also unlikely(4)We believe that the diffractions on the shot record are from subsurface small-scale geological bodies. Fault planes and pinch-outs are common geological features that originate diffractions. We did not see any pinch-out structure in the migrated section, though a fault plane associated with small-scale faults may explain some of the diffractions observed on the shot record and on the stack profile. However, faults planes cannot explain the diffraction with reverse polarity and the string of bead reflections located in the otherwise blanking zone. Therefore, we conclude that some of the diffractions, particularly the strong diffractions on the shot record from 1750 to 1850, are associated with small-scale bodies within the fluid migration pipe
Based on the modeling results (Figure 13), when the source location is close to the diffractors, the relationship between the diffraction and the reflections could tell the spatial relationship between the diffractor and the layer boundaries. From this perspective, based on the shot records close to the pipe structure (e.g., 1754, 1758, 1764, and 1768), we could infer that diffraction D1 and D2 are responses to the diffractors between H1 and BSR and contribute to the SBRs in zone B in Figure 2. The diffraction D3 is a response to the diffractors below BSR and contributes to the SBRs in zone C.
5.2. The Nature of the Diffraction Bodies
Numerical modeling conducted in the case of carbonate caves has shown that even one single cave may produce several strings of bead reflections [39]. But it is unlikely that the long SBRs within the pipe are produced by only one single diffractor. First, at least three zones of SBRs could be identified on the real migration profiles (Figure 2(b)). Moreover, several distinguishable diffractions, such as D1, D2, and D3, are observed on real-shot records (Figures 6 and 7). It is reasonable to infer that at least three small-scale bodies are located within the pipe structure. Further, we could propose that at least one diffractor is located between the seafloor and H1, one between the H1 and the BSR, and one below the BSR.
As to the material nature of the bodies, impedance contrasts between the body and the surrounding sediments could provide vital clues. Free gas, gas hydrate, and carbonate are three common components, associated with the gas migration process in the gas hydrate and cold seep areas. The presence of gas hydrate or carbonate will increase the acoustic velocity and therefore the impedance, while the presence of even a small amount of gas will dramatically decrease the wave velocity and consequently decrease the impedance [48]. The reversed polarity of the diffraction (Figure 9) suggests that the corresponding body has a lower impedance compared to the surrounding sediments. We interpret this feature as a gas pocket or gas reservoir. The presence of gas pockets has been widely observed in the shallow sediments close to the seafloor [49–52]. Considering that free gas must migrate upward through the pipe to feed the active seep, the presence of gas pockets within the pipe is reasonable. In summary, we conclude that the bodies are low velocity anomaly and interpreted as gas pockets.
As to the size of a small-scale body, the lateral size can be inferred from the two-dimensional seismic profile. From the profile displayed in Figure 2(b), it could be estimated to be about 30 m. However, it cannot be as easy to determine the height, as this would require a detailed quantitative seismic interpretation based on the inversion of the seismic data for high-resolution elastic models [53]. According to previous physical modeling [39], a short SBR could be the response of a cave with a height of less than 60 m, while a long SBR could be the response of a single cave with larger height or two caves distributed with a distance less than 60 m. The physical modeling results may help to estimate the height of the diffractors in our work. However, estimation based on these results could be too rough. Considering that the real data has a dominant frequency of 60 Hz, we may obtain better estimation based on the observations from real data and from modeling. From the shot records, most of the diffractions are characterized by a single event. This may lead us to infer that the height of the small-scale body is less than the resolution of the data. Assuming that the velocity is 1500 m/s, the wavelength is about 25 m. Taking 1/4 wavelength as the criteria, we may infer that the height is less than 7.5 m. From the modeling results (Figure 13), when the height is set to 10 m, the diffraction is characterized by two very close events. As the height increases, the events become more separated. Hence, we could infer that the height is less than 10 m from the modeling side. Combining the two independent estimates, we may provide a good estimate of the height of the diffractor and the height is probably less than 10 m.
5.3. Gas Migration Pattern beneath the Cold Seep
Pipe structures have been observed in many geological settings and are widely accepted as the most important pathway of fluid migration. Commonly, pipes are characterized on the seismic data by columnar zones of disrupted reflection continuity and/or the complete absence of reflection due to the attenuation of gas or the scattering effect [17]. Pull-up or push-down reflections are also commonly associated with fluid escape pipes. Pull-up deformation is considered to be caused by the presence of hard materials (i.e., with high wave velocity), such as gas hydrates or carbonates [54, 55], while the push-down deformation is thought to be caused by the presence of free gas [52]. The diversity of seismic expressions of pipe structures may indicate a diverse pattern of the fluid migration process.
In the study area, the migration pipe shows a different type of seismic expression and is characterized by a series of short and strong reflections stacked vertically (Figure 2). We have interpreted that a series of gas pockets exist in the pipe. Based on this interpretation, we can infer the fluid migration process as displayed in Figure 14. In stage one, free gas migrates upward and accumulates beneath the BSR due to the sealing effect of the gas hydrate layer (Figure 14(a)). In stage two, gas pressure exceeds the overburden and keeps migrating upward and forms a series of gas pockets (Figure 14(b)). In stage three, gas migrates toward the seafloor producing gas bubble plumes (Figure 14(c)). In this process, the presence of gas pockets can be analogous to the presence of magma chambers beneath volcanoes [56, 57]. The gas pockets play the role in the cold seep as the magma chamber plays in the volcano process. The presence of gas pockets may be determined by the lithology.

A conceptual model is proposed to describe the gas venting process. First, gas leaks from the deep reservoir along the faults induced by the uplift. Gas accumulates below the BSRs due to the sealing effect of gas hydrates and forms the weak amplitude zone. The pipe structure forms owing to the overpressure induced by continuous accumulations of free gas. Gas migrating upward through the pipe enables the formation of shallow gas hydrates and produces cold seeps on the seafloor.
6. Conclusions
A pipe structure beneath an active cold seep is studied in this work by analysis of diffraction waves from real and modelled synthetic data. On the migration profiles, the pipe is characterized by a series of short vertically stacked strong amplitudes, very similar to the seismic response of carbonate caves. Diffractions are abundant on shot records and are caused by small-scale bodies (or diffractors) within the pipe structure. Reversed polarity with respect to the seafloor reflections confirms that the bodies have lower -impedance than the background material. These types of small-scale bodies are interpreted to be gas pockets within the pipe. The lateral size of the small-scale bodies is about 30 m inferred from the migration profile. The vertical size is estimated to be less than 10 m by combining the resolution analysis of the real data and the qualitative analysis of the modeled results..
A conceptual model is proposed to describe the gas venting process to interpret the active Haima cold seeps associated with the fluid escape pipes. Beneath the BSR, gas exceeds the overburden pressure and keeps migrating upward and gas pockets emerge, which may be determined by the lithology. The presence of gas pockets is the core of this interpretation. The presence of gas pockets within fluid escape pipes can be analogous to the magma chambers beneath volcanoes. From this point of view, we may gain a new insight into how free gas migrates upward through the pipe and reaches the seafloor.
Data Availability
The synthetic data and the matlab code used to support the findings of this study are available from the first author upon request ([email protected]). The field data is not available due to the confidential issue.
Conflicts of Interest
The authors declare that there is no conflict of interest regarding the publication of this paper.
Acknowledgments
Great thanks are due to the crew members and scientists who participated in the acquisition of the multibeam data in 2016 and seismic data in 2007 and 2013. We thank the GMGS for their permission to publish the results. This study was financially supported by the Dedicated Fund for Promoting High-Quality Economic Development in Guangdong Province (Marine Economic Development Project) (GDNRC [2020] 045), the Laboratory for Marine Mineral Resources, Qingdao National Laboratory for Marine Science and Technology (Grant no. MMRKF201810), and the Key Special Project for Introduced Talents Team of Southern Marine Science and Engineering Guangdong Laboratory (Guangzhou) (GML2019 ZD0207). Dr. Jiangxin Chen is funded by the Shandong Province “Taishan Scholar” Construction Project. LA acknowledges the support from the CERENA (strategic project FCT-UIDB/04028/2020).