Stress Constraints From Shear-Wave Analysis in Shallow Sediments at an Actively Seeping Pockmark on the W-Svalbard Margin

Mechanisms related to sub-seabed fluid flow processes are complex and inadequately understood. Petrophysical properties, availability of gases, topography, stress directions

. Several tens of pockmarks occur at the crest of Vestnesa Ridge, a sedimentary drift on the western-Svalbard Margin (Figures 1b, 1c, and 2;Bünz et al., 2012;Vogt et al., 1994).Many of the pockmarks on the eastern segment of Vestnesa Ridge actively seep gas into the water column and show gas-chimney features in subsurface data below them (Bünz et al., 2012;Plaza-Faverola et al., 2015;Vogt et al., 1994).Why seepage is sustained exclusively at few pockmarks on the eastern Vestnesa Ridge segment is the question that motivates our study of shallow stress and anisotropy.Previous studies have hypothesized that the eastern Vestnesa Ridge segment is under an extensional stress regime that facilitates fracture/fault dilation and favors gas release via fluid advection (Plaza-Faverola et al., 2015;Plaza-Faverola & Keiding, 2019).However, the western segment may be under a more compressive regime, where any potential gas release would occur at lower rates, for example, via diffusion (Plaza-Faverola & Keiding, 2019;Ramachandran et al., 2022).
Permeability, availability of gases along with various other factors (e.g., topography, stress directions, hydrate formation, pore pressure and other petrophysical properties like porosity, grain size, etc.) control the seepage through a gas chimney (J.Liu et al., 2019;Ramachandran et al., 2022).Chimneys have a higher secondary permeability (i.e., the one created by faults/fractures/cracks) compared to surrounding sediments, mainly due to a series of interconnected subvertical or radial fractures, which facilitate the vertical flow of gas as well as gas hydrate accumulations in the shallow subsurface (Cartwright et al., 2007;Cartwright & Santamarina, 2015;X. Liu & Flemings, 2007).Gas overpressure forms self-enhanced hydrofractures, which are a major cause of secondary permeability in gas chimneys (Daigle et al., 2011;Ramachandran et al., 2022).Permeability, especially one driven by faults and fractures, is dependent on stress orientations and often changes with a change in stress direction (Rutqvist, 2015;Zheng et al., 2015).The present-day effective stress fields might help with opening or closing of pre-existing faults and fractures, depending on orientations of stress fields and faults/fractures.Gravitational stress due to seafloor topography is a particularly important contributor in the overall stress field at shallow depths.It is also important to study at which depths tectonism-related stresses dominate over gravitational/slope related stresses (Haacke et al., 2009).The in-situ measurement of stress in marine sediments is expensive and instruments do not penetrate very deeply (<5 m) (Boggess & Robertson, 2011;Lunne et al., 1997).Therefore, remote methods are necessary to study stress fields in sediments.The most conventional approach relies on subsurface seismic imaging of faults and fractures, developing tectonic stress models and linking them together to get a comprehensive understanding of the subsurface stress field.
By implementing integrated S-wave and 3D seismic data analysis, our aim is to study subsurface anisotropy in shallow (<250-300 m) sediments below the seafloor at an active seep site on Vestnesa Ridge.Concretely, we are interested in understanding the relationship between anisotropy, the orientation of faults and fractures and seepage distribution through vertical pathways.We discuss our observations in relation to processes that exert a control on gas seepage activity in this deep marine Arctic setting (e.g., topography, sediment cohesion, gas saturation and distribution, pore pressure, lithology, etc.).

Study Area
Vestnesa Ridge is a contourite sedimentary drift in the western Svalbard margin and northeast to the Molloy Transform Fault (MTF) in ∼79°N (Eiken & Hinz, 1993;Vogt et al., 1994; Figures 1b and 1c).The ridge is at a water depth of ∼1000 to ∼1700 m and has up to ∼100 m elevation and ∼1° slope relative to the deeper side toward the continental shelf (Figure 1b).The sediment drift bends anticlockwise from SE-NW to ESE-WNW in the north (Figure 1b) where it terminates at the Molloy Ridge (MR).In previous studies, the SE-NW and ESE-WNW oriented segments are referred as the eastern and western segments, respectively; and we will also follow the same naming convention (Figure 1b).The occurrence of contourite mounds, moats and migrating wave features in the region indicate sedimentation under the influence of bottom water contour currents (Eiken & Hinz, 1993;Ottesen et al., 2005;Stein et al., 2005).Sedimentary deposits along the western Svalbard margin slope consist mostly of turbiditic, glaciomarine and hemipelagic sediments (Howe et al., 2008).Sediment core analyses from nearby Ocean Drilling Program wells reveal a high (∼105 cm/year) sedimentation rate, mostly of silty turbidites, from the mid Weischselian to the Last Glacial Maximum (LGM) and followed by a relatively slow (<10 cm/year) sedimentation rate, mostly of muddy-silty contourites with abundant ice-rafted debris, from the LGM to the early Holocene (Eiken & Hinz, 1993;Forsberg et al., 1999;Howe et al., 2008).Acoustic flares from six pockmarks on the eastern segment of the ridge are profusely evident (up to 400 m water depth from the seafloor) in sonar data (Smith et al., 2014).Lunde Pockmark is one of them (Panieri et al., 2017;Figures 1c and 2).It is a ∼400-500 m wide complex structure, with small-(∼0-10 m) to medium-scale (∼10-50 m) depression features, referred to as pits, within the overall perimeter of the pockmark (Figure 1c; Himmler et al., 2019;Panieri et al., 2017).Hydroacoustic and seismic data show an intense advective gas release occurring through these pits, whereas slow diffusive gas release is more likely to dominate elsewhere within the pockmark (Hong et al., 2016(Hong et al., , 2021;;Panieri et al., 2017).Present-day seepage is not evident from sonar data in the western segment of the ridge; however, the occurrence of pockmarks and gas chimney features indicate fluid seepage activity in the past (Plaza-Faverola et al., 2015).Even at the Lunde Pockmark site, seepage activity has changed over time scales of thousands of years; and several methane seepage episodes have been documented via sampling and dating of authigenic carbonate for the last ca.160,000 years (Himmler et al., 2019).

Structural Background
The continental break-up and seafloor spreading in the early Eocene along Reykjanes, Aegir, and Mohns Ridges, followed by a change in the plate motion (∼33 Ma) between Svalbard and Greenland from strike slip to oblique divergence, led to the initial formation of Fram Strait (Demenitskaya & Karasik, 1969;Eldholm et al., 1987;Myhre & Eldholm, 1988;Talwani & Eldholm, 1977;Vogt, 1986).The shearing along faults between Greenland and Svalbard resulted in the Western Spitsbergen Orogeny, which continued until the early Oligocene, when the spreading direction changed from NNW-SSE to NW-SE (Harland et al., 1974;Steel et al., 1985).
The opening formed Spitsbergen Shear Zone, which developed into an asymmetric, ultra-slow and obliquely spreading ridge system in the region.This rifting system continued as an asymmetric pure shear or a high angle simple shear mode, pivoted around a system of faults adjacent to the continental margins of Svalbard (G.L. Johnson et al., 1972;Kovacs & Vogt, 1982;Nunns, 1983;Nunns & Peacock, 1983;Vogt et al., 1982).By the end of the Oligocene (∼23 Ma), rifting shifted northwards along MR.The process of continental break-up continued, opening up Fram Strait by 10-15 Ma for deep water circulation, establishing bottom water-current driven sedimentary drifts (Ehlers & Jokat, 2009;Jakobsson et al., 2007;J. E. Johnson et al., 2015).
At present, shearing along Spitsbergen Transform Fault and rifting at Knipovich and Molloy Ridges influence overall stresses at a regional scale.The effect of these forces can also be seen in the formation of half graben structures with bounding faults rooted deep in the basement and dipping toward the ridge axis (Amundsen et al., 2011).These deep-rooted faults in the region provide pathways for fluid flow and bring deep crustal fluids into the sedimentary cover (Madrussani et al., 2010;Waghorn et al., 2018).Analytical tectonic stress modeling predicts an extensional stress, due to oblique spreading, in an area extending northward from Knipovich Ridge (KR), encompassing the eastern Vestnesa Ridge segment (Plaza-Faverola & Keiding, 2019).Glacial stress modeling predicts low magnitudes (<6 MPa) of maximum and minimum horizontal stresses in the region comprising the eastern Vestnesa Ridge segment due to isostatic adjustment of the ice-sheets forebulges (Vachon et al., 2022).Glacial stresses are in addition to the background tectonic stresses, whose magnitudes at Vestnesa Ridge are unknown.A northward progradation of the extensional stress regime from KR affects fault behavior in the eastern segment of Vestnesa Ridge (Crane et al., 2001;Plaza-Faverola & Keiding, 2019;Vanneste et al., 2005).NW-SE trending sedimentary faults are aligned nearly parallel to the modeled maximum tectonic stress direction in the eastern segment of Vestnesa Ridge, which favors the occurrence of seepage (Plaza-Faverola et al., 2015;Plaza-Faverola & Keiding, 2019).However, whether overpressured fluids create faults and fractures or pre-existing faults and fluids act as pathways for the upward migration of fluid is still an unanswered question in this area.An upward fluid flow through faults can be obstructed at shallower depths, because of the formation of gas hydrates or due to changes in stress regimes in shallow sediments, as both effects can reduce secondary permeability in the system (J.Liu et al., 2019).
The preferential alignment of minerals cannot cause anisotropy in shallow (<300 m) poorly consolidated sediments at Vestnesa Ridge, as this process becomes dominant in fully consolidated metamorphic rocks (Zinke & Zoback, 2000).In addition, the hypothesis of aligned minerals affecting shear-wave polarization is not a universally accepted, as Crampin and Peacock (2008) give strong evidence in favor of rejecting this hypothesis.
Near-vertical faults and fractures often expressed as vertical fluid migration pathways along Vestnesa Ridge (Plaza-Faverola et al., 2015;Singhroha et al., 2016Singhroha et al., , 2020) ) suggests that the orientation of faults and fractures may be an important cause of anisotropy in this setting.Tensional faults and fractures form perpendicular to the minimum horizontal stress direction (C.J.Evans & Brereton, 1990;Pastori et al., 2019;Zinke & Zoback, 2000).In shallow sediments, where vertical stress dominates over horizontal stresses, hydro-fractures, created by trapped overpressured fluids, are also typically preferentially oriented perpendicular to the direction of minimum principal stress (Hubbert & Willis, 1957).Crampin and Peacock (2008) argue that stress-aligned fluid-saturated microcracks should be the default interpretation of azimuthally aligned S-wave splitting, and they further argue that other interpretations are either fallacious or have few supporting arguments.Small-scale seismic wave induced stresses interact in a non-linear fashion with effective differential stresses in fluid saturated media, which results in the polarization of fast-S-waves aligned parallel to the maximum principal horizontal stress direction (Crampin & Peacock, 2008;Crampin & Zatsepin, 1997;Zatsepin & Crampin, 1997).It is therefore reasonable to associate S-wave splitting from sub-seabed sediments of Vestnesa Ridge with principal horizontal stresses, with the polarization of fast S-waves occurring parallel to the maximum horizontal principal stress direction (more arguments provided in Text S1 in Supporting Information S1).

Data
During a research expedition (CAGE-2019-1: 29 June 2019 to 8 July 2019) on-board R/V Helmer Hanssen, we acquired OBS data using an array of 22 ocean bottom seismometers (OBS) around Lunde Pockmark on the eastern segment of Vestnesa Ridge, W-Svalbard Margin (Figure 1).The OBS instruments were equipped with a three-component seismometer (4.5 to a few 100 Hz frequency range and the sensitivity of ∼28.8 V/m/s) and a hydrophone to record the seismic wavefield.
The OBS distribution and survey configuration were intended to maximize the chances of studying the anisotropy in and around the pockmark (Figures 1c and 1d).OBS stations were initially planned as a dense grid with 250-350 m spacing between stations covering an area close to and within Lunde Pockmark.However, it was not possible to deploy the OBS equipment with a high spatial precision on the seafloor, due to the failure of an ultra-short baseline release system.Consequently, the OBSs were freely dropped in the water, and they drifted sideways as they sank through the water column, landing on the seafloor at locations different to the planned positions (Figures 1c and 1d).We used the least-squares inversion of direct-arrival travel-times to precisely relocate the instruments.A large number (>33,000) of shots from different directions constrained OBS positions, with an expected error in OBS locations of less than 1 m (Figures 1c and 1d; Table S1 in Supporting Information S1).The position of OBSs with respect to the seafloor pockmark was considered suitable for the analyses of S-waves despite the divergence from the planned positions.
A total of 50 seismic profiles (i.e., air gun trajectory) were shot around these stations (Figure 1d).In addition, we recorded several shots along circular trajectories (Figure 1d).These inlines, crosslines and circular trajectories ensured a good azimuthal coverage (Figure 1d).Weather conditions were good (<1 m high waves) throughout the experiment.The seismic signal was generated using a mini generator-injector air gun (Sercel) configured as 30/30 in 3 .The air gun generated signal with frequencies up to 300 Hz at the shot interval of 5 s.OBS data were recorded at a sampling rate of 0.5 ms (OBS1-20) and 0.4 ms (OBS21-22).Data processing after relocation of the instruments included band pass filtering using frequency ranges between 5 and 100 Hz, depending on the signal to noise ratio characteristic of each instrument.Low frequency components of the signal were preserved, since a significant portion of the converted shear energy appear in the low frequency range (<20 Hz).

Data Quality Control (QC)
Down-going P-wave energy is reflected into upcoming P (PP)-wave and converted S (PS)-wave.PS-waves are recorded in horizontal seismometer components and are used in our S-wave analysis.We analyzed the quality of recorded S-wave data using three parameters/approaches: (a) signal/noise (S/N) ratio, (b) seismometer tilt and (c) vector fidelity analysis (Table 1).Any tilt from the vertical axis due to the placement of the seismometer will lead to a distribution of S-and P-wave energies in all three components of the seismometer.The S/N ratio decreases with increase in the noise, and the S-wave energy increases in the vertical component of the seismometer with increase in the seismometer tilt.The S/N ratio of S-waves in horizontal components and the seismometer tilt can be visually inspected by looking at all three components of the recorded OBS data set (Table 1 and Figure 3).The absence of good strength of S-wave energy in horizontal components is primarily due to the poor coupling between sediments and the seismometer.A poor coupling can occur due to bad placement of the seismometer on the seafloor.An ideally placed seismometer with no tilt will record S-and P-wave energies in horizontal and vertical components of the seismometer, respectively (Figure 3).By visual inspection we classify the OBS data sets as good/medium/bad based on the S/N ratio and seismometer tilt (Table 1).

Vector Fidelity Analysis
The vector fidelity analysis is carried out to ensure that the orientation of particle motion in recorded data matches with the orientation of the arriving wavefield (e.g., Exley et al., 2010;Figure 4 and Figure S1-8 in Supporting Information S1; Text S2 in Supporting Information S1).A bad vector fidelity can be primarily due to the bad quality of sensors inside a seismometer (e.g., differences in the sensitivity of horizontal components), instrumentation (e.g., electronics) related noises, and roughness of the surface on which a seismometer is placed (Exley et al., 2010).
We analyze the vector fidelity of data recorded in the different data sets using the particle motion of first S-wave arrivals from different shots (Figure 4).The consistent median difference observed in the particle motion and the shot-receiver directions is used to derive the orientations of the two horizontal components.After correcting for the horizontal-component orientations, the deviations observed in different shots between particle motion and shot-receiver directions are used to study the vector fidelity of the seismometer (Figure 4).In a seismometer with a good vector fidelity, the orientation of particle motion inferred from the first arrivals will roughly match the orientation between source and receiver (Figures 4a and 4c).On the contrary, for data sets with a bad vector fidelity, the particle motion direction cannot be related with sufficient confidence to the direction of the arriving wavefield (Figures 4b and 4d; more details on vector fidelity analysis provided in Text S2 in Supporting Information S1).This Quality Control (QC) process ensures that the particle motion analysis using S-waves converted at greater depths is reliable.Based on the vector fidelity analysis results, we classify different OBS stations as good/ medium/bad (Table 1).
Using the three approaches discussed above, we qualitatively classify data from all OBS stations (Table 1) and select nine OBS stations for further analysis, which satisfy all three conditions as either medium or good.Out of these nine stations, two stations (OBS13 and OBS15) can be classified as good stations and one station (OBS8) as a very good station, with a good/very good quality data assessment from all three approaches used for classification (Table 1).

Methods
When an incoming P-wave encounters a seismic reflector, a part of its energy gets reflected in the form of P-wave (PP) and another part gets reflected as S-wave (referred as converted/PS-wave) (Figure 5).The polarization of this PS-wave in an isotropic medium is in the radial direction (Figure 5).The radial direction is the direction from which the incident P-wave arrives to the receiver, and perpendicular to this is the transverse direction (Figure 5).The distribution of PS energy in radial and transverse components depends on the anisotropy in a medium (Exley et al., 2010;Haacke et al., 2009;Haacke & Westbrook, 2006;Robinson et al., 2022).The PS-wave energy splits along fast and slow S-wave polarization directions in an anisotropic medium and hence the transverse direction also receives energy (Exley et al., 2010;Haacke & Westbrook, 2006).When the radial direction (the direction of arriving P-wave energy) matches with a plane of symmetry in an anisotropic medium, the amount of energy in the transverse direction is smaller compared to other directions (Crampin, 1981).Hence, by looking at the PS-energy in the transverse component, we can reveal the splitting of energy in fast and slow S-wave components and can also identify planes of symmetry in an anisotropic medium.
We assume that the geological system we are investigating (i.e., characterized by near-vertical faults and fractures) can be best studied by assuming horizontal transverse isotropic (HTI) medium (i.e., a medium with vertical faults and fractures predominantly oriented along one orientation).The planes along and perpendicular to vertical faults and fractures are symmetric planes in a HTI medium.If the radial direction matches with one of the symmetry planes of a HTI medium, the energy in PS transverse component is less compared to other directions (Pastori et al., 2019;Schutt et al., 1998).This happens because the planes of HTI symmetry also match with fast and slow S-wave polarization directions, and S-wave splitting is minimal when radial direction matches with fast or slow S-wave polarization direction (Crampin, 1981).Hence, a relative decrease in energy along certain azimuths in the transverse direction indicates axis of HTI symmetry, fast or slow S-wave polarization direction, and direction parallel or perpendicular to predominant fault-and fracture-orientation (Li, 1998;Robinson et al., 2022).The fast-polarization direction is expected to be parallel to the maximum horizontal stress (Shmax) direction (Crampin & Peacock, 2008;Pastori et al., 2019).
The analysis of S-wave splitting in an anisotropic medium with HTI symmetry can be best performed by transforming the data into an alternative radial-transverse coordinate system (Exley et al., 2010;Haacke & Westbrook, 2006).The radial-transverse rotation of data can be achieved in two ways: (a) reorienting data by maximizing the first S-wave arrival energy in the radial component of each shot, or (b) finding an optimum solution (using the entire data set) for the orientation of horizontal components of a seismometer from the particle motion of first S-wave arrivals, and then using these orientations to align energy parallel to the arriving wavefield (radial) and perpendicular to the arriving wavefield (transverse) directions (Figure 5).Considering that a seismometer gets placed in the sediment at the beginning of the survey, the directions of horizontal components are fixed throughout the survey.Hence, we use the second approach and perform the radial-transverse rotation of data using estimates of orientations of horizontal components and shot-receiver directions.This approach does not always ensure maximization and minimization of the first S-wave arrival energy in radial and transverse components, respectively.However, it provides a better rotation of data from greater depths, as first S-wave arrivals from a single shot can be influenced by a local very-shallow depth anisotropy and may not be an optimum criterion for the rotation of the entire data from that shot.
The lack of energy in certain azimuths in a transverse component, hereafter referred to as energy nulls, indicates the axes/directions of symmetry planes in an anisotropic HTI medium (Figure 6) (i.e., based on the HTI assumption).In a homogeneous anisotropic medium, peaks and troughs of energy occur twice (i.e., 2 energy minima with 180° azimuth spacing) in the radial component and four times (i.e., 4 energy nulls with 90° azimuth spacing) in the transverse component (Exley et al., 2010).In a HTI medium, one of the energy null symmetry axes in the transverse component corresponds to the fast-polarization direction, whereas the other symmetry axis (i.e., perpendicular to the first axis) corresponds to the slow-polarization direction (Robinson et al., 2022).

Shot Selection, Binning and Stacking
The entire survey consists of 33,774 shots.We sort recorded shots using their azimuth and offset for each OBS station.Shots are grouped in 360 azimuths (with each shot grouped to the nearest azimuth).For example, the azimuth of 1° contains all shots with azimuths lying between 0.5° and 1.5°.Next, we also sort the  common azimuth in relation to offset, with offsets grouped in 250 m long radial bins.Then we combine shots from azimuth and offset groups in their corresponding bins.For example, all shots from azimuths of 0.5°-1.5°and offsets 750-1000 m are put together in one bin.After applying moveout correction (using velocity analyses from Singhroha et al., 2019) and rotating the data to the radial-transverse domain, the shots in each bin are stacked together.
We have stacked data to 2.4 s arrival-time to not include P-wave multiple reflections (Figure 3), which would negatively affect the S-wave analysis.This is especially important since we aim to find azimuths where seismic energy is low, and the presence of noise will degrade the overall quality of contrast between azimuths with low and high energy in the transverse component section.The primary PP signal will also contaminate the data in the 0.8-1.1 s arrival-time interval.We expect a good separation of low and high energy azimuths within the 1.1-2.4s arrival-time interval, where there is no primary or multiple PP energy (Figure 3).
Energy conversion from P-wave to S-wave at a seismic reflector (Figure 5) is a function of the incident angle, and by dividing data in offset groups of 250 m, we take the utmost care to avoid a biased interpretation due to differences in incident angles.However, we find no significant differences in amplitudes, neither in the radial nor the transverse component for 750-1000 m, 1000-1250 m and 1250-1500 m offset groups.Here, we present stacked results with 750-1500 m offset to ensure a complete azimuthal coverage and a high signal to noise ratio.The offset range (750-1500 m) is optimal in terms of improving the overall detectability of amplitude nulls, while ensuring the effect of angular dependence of PP-PS energy conversion remains insignificant.

Results
The occurrence of S-wave splitting is evident in different OBS stations (Figures 6 and 7).In most OBS stations, there are visible variations in the amplitude with azimuth in the transverse component (Figures 6 and 7).All OBSs commonly show energy along all azimuths in the radial and transverse sections within the 0.8-1.1 s arrival-time interval.Barring a few exceptions (Figures 6-8): in OBS6, OBS9 and OBS11, we see azimuths with a low energy (e.g., azimuths of: 0°-40°, 180°-220° in OBS6 radial, 115°-125° in OBS9 transverse, 350°-70° in OBS11 radial and 310°-350° in OBS11 transverse components), although most of them are associated with polarity changes (Figures 6a,6b,6g,6h,7a,and 7b).Energy anomalies in the radial and transverse sections have clear separation in OBS6, OBS8 and OBS9 compared to OBS7 or OBS13 (Figure 6).The vector fidelity analysis shows that OBS6, OBS8 and OBS9 have better quality than for example, OBS7.As expected, stations with the high-quality vector fidelity, separate transverse and radial energies better and amplitude nulls are sharper at these sites.
OBS15 is on the northeastern rim of Lunde Pockmark (Figure 9a).In OBS15, there are some zones of relatively low energy in the transverse component, but they are not as sharply defined as zones observed in most of the other OBS stations (Figure 7f).Azimuths of 125°-150°, 220-245°, 200°-315°, 345°-355° show relatively low energy (Figure 7f), however, there is some energy which appears in-between, and the azimuthal nulls are not distinct (Figure 7f).The locations of OBS17 and OBS22 are separated only by 15 m, and their seafloor positions lie northeast of Lunde Pockmark (Figures 1c, 1d, and 9d).Transverse components of both OBS show an absence of amplitude nulls (Figures 8b  and 8d).However, in both OBS, we observe a sharp change in amplitude strengths at azimuths of 133° and 314° across the entire arrival-time window, in both radial as well as in the transverse component (Figure 8).
We integrate results from shear-wave analysis with high-resolution 3D seismic data (Plaza-Faverola et al., 2015).Seismic attribute analysis ("variance" and "ant tracking") of this data allows us to statistically measure predominant fault and fracture orientations in the subsurface around and within the Lunde pockmark.Details on the methodology of this approach are included in Appendix A.

Uncertainties and Limitations
Before discussing the implications of the observed PS-wave energy distribution for understanding the fluid migration system and associated stress field at Lunde we discuss data uncertainties and the limitations of the approach implemented.

S-Wave Polarization Changes With Depth
One major assumption in our analyses is that fast and slow S-wave polarization directions does not change with depth.Amplitude null symmetries represent summation of effects of fast and slow S-wave polarization directions of different anisotropic layers through which S-wave travels (Robinson et al., 2022).Interpretation of fast and slow S-wave polarizations from amplitude null symmetries requires a layer-stripping approach, especially at greater depths (>150-200 m below the seafloor [bsf]), because of changes in S-wave polarization directions with depth (Haacke et al., 2009;Robinson et al., 2022;Slack et al., 1992).However, in shallow sediments (<150 m bsf), particularly in cases with no other visible symmetries, amplitude null symmetry directions correspond to fast and slow S-wave directions (Crampin, 1993a(Crampin, , 1993b;;Robinson et al., 2022;Slack et al., 1992).Since we primarily study shallow sediments and we do not see changes in amplitude null symmetry with depth, we avoid the application of relatively complicated layer-stripping approach.Because most amplitude null symmetries are unambiguous, and we analyze data only up to 2.4 s (∼220-250 m bsf, depth converted considering S-wave velocity model from Singhroha et al., 2019), the azimuths of amplitude null symmetry correspond to the planes of symmetry in the anisotropic medium (Crampin, 1993a(Crampin, , 1993b;;Pastori et al., 2019;Robinson et al., 2022).

Noise From PP Arrivals
The second assumption in our analyses is that the interference of PP signal is insignificant.Such an assumption here allows us to interpret the absence of amplitude nulls between the arrival-time 0.8-1.1 s (∼0-20 m bsf) as the absence of the planes of symmetry in the anisotropic shallow sub-surface and the potential lack of horizontal stress differences in sediments near the seafloor (Figures 6-8).Nevertheless, the upper 20 m in the sedimentary column is a low confidence interval and any interpretation should consider potential interference with PP energy.

Inhomogeneity
The third and last major assumption is that the stratum where the PS-wave energy is propagating is homogeneous.This assumption allows us to relate the amplitude nulls to symmetry planes, and therefore to slow and fast energy traveling directions.The presence of vertical cracks in a medium makes it homogenous and anisotropic, whereas the presence of inhomogeneity (e.g., the presence of randomly distributed methane-derived authigenic carbonates, gas hydrates or free gas in a gas chimney) can create amplitude contrasts due to differences in reflectivity and seismic attenuation.In OBS data, since the point of reflection varies with azimuth and offset, the presence of symmetries will not occur at regular 180°/90° intervals in an inhomogeneous medium.Our case study is not an ideal homogeneous medium (i.e., we expect inhomogeneity at small-scales [<10 m], especially at OBS stations which lie within the pockmark, due to changes in the sediment and pore-fill properties [Figures 1 and 2]) and this is considered when discussing energy nulls in radial and transverse components at different OBS stations.
The seafloor position of OBS11 lies 30 m away from a pit (at azimuth of 320° from OBS11; marked as pit1 in Figures 12a and 12c) and 70 m away from another pit (at azimuth of 10° from OBS11; marked as pit2 in Figures 12a  and 12b) within Lunde Pockmark.We see a significant decrease in S-wave energies arriving from azimuths of 320° and 10° (Figures 7a and 7b), however, we do not see 180° symmetry in amplitude nulls.The overall energy in azimuths of 270°-70° is significantly lower than the energy in azimuths of 70°-270° (Figures 7a and 7b).The rest of Lunde Pockmark and the gas chimney below, also lies along azimuths of 270°-70° from OBS11 site (Figure 12).The observed anomalies in OBS11 radial and transverse components are most likely due to inhomogeneities created by pits, carbonates, or a free gas, and the effect of differences in energy with azimuth are possibly related to energy scattering or seismic attenuation related issues (Figure 1c).The pits within the pockmark often have uneven structures which can scatter arriving seismic energy (Panieri et al., 2017; Figure 12a).The occurrence of pits on the seafloor and strong-reflections around 1.70-1.75s two-way travel time, indicates the presence of carbonates with high-reflectivity (which potentially can scatter most of the arriving energy) or free gas with high-attenuation (Figures 12b and 12c; Bünz et al., 2012;Himmler et al., 2019;Singhroha et al., 2016).The vertical anomalous changes in S-wave energy, in OBS17 and OBS22, also appear to be the result of inhomogeneity (Figure 8).Orientations of azimuths corresponding to this anomaly are aligned parallel to the ridge crest.On other OBS sites near Lunde Pockmark (OBS7, OBS9, and OBS15), we observe 180° symmetry in amplitude nulls, suggesting that inhomogeneity due to the pockmark and gas chimneys plays an insignificant role in our observed results.

Discussions
Keeping in mind that some of the PS-wave energy distribution trends in the transverse and radial components may reflect amplitude attenuations due to chimney related inhomogeneities and may be mixed with PP energy arrivals, we analyze the potential relationship between amplitude null symmetries, anisotropy in sub-seabed stress fields and fracture orientations.
First, the occurrence of energy in the transverse direction per se (Figures 6-8) confirms the existence of sub-surface anisotropy.In addition, the presence of amplitude null symmetry in the transverse direction (Figures 6, 7, and 9) can be correlated to the symmetry axis in the anisotropic shallow sediments below and southwest of Lunde Pockmark.Assuming a horizontal symmetry axis (HTI), the observed azimuths of amplitude null symmetries in the transverse component can be correlated with directions parallel and perpendicular to the predominant fault/ fracture strike direction (Li, 1998).
In OBS6-9 and OBS13, we see consistency in the azimuths of amplitude null symmetries in the 1.1-2.4s arrival-time interval (∼20-250 m bsf) (Figures 6 and 7).Amplitude null symmetries vary within the azimuth ranges of (40°-50°)-(130°-140°)-(220°-230°)-(310°-320°) (Figures 6 and 7).We observe a change in the S-wave amplitude patterns across the Vestnesa Ridge crest from the southwestern side (OBS6-9) to the northeastern side (OBS15, OBS17, and OBS22) of Lunde Pockmark (Figures 1 and 6-9).There are no distinguishable symmetry planes at OBS17 and OBS22 sites (Figure 8).Even OBS15, close to the ridge crest, shows a significantly different strength of amplitude nulls compared to OBS6-9 (Figures 6 and 7f).We interpret these significant differences in amplitude strengths as an indication of changes in the anisotropic conditions on the Vestnesa Ridge crest and across Lunde Pockmark.Southwest of Lunde Pockmark, the orientation of faults and fractures are dominantly NW-SE, presumably creating symmetry planes in the anisotropic medium and explaining the observed amplitude null symmetry (Figures 6,9e,and 9g).However, across the ridge crest (northeast of Lunde Pockmark), faults and fractures do not show a dominant orientation and the amplitude null symmetry planes are therefore not expected (Figures 8,9e,and 9g).Overall, the structural attribute analysis suggests that the likelihood of fracturing northeast of Lunde Pockmark is lower for all the azimuths in comparison to the predictions for the southwestern region (Figure 9).Therefore, the strata toward the northeast of Lunde Pockmark are characterized by the absence of amplitude nulls and by a low likelihood of fracturing (Figures 8 and 9; example, OBS17 and OBS22).
Assuming that hydrofracturing develops documented faults and fractures, southwest of Lunde Pockmark (i.e., the area within the blue rectangle in Figure 9a), our study suggests that there is a ∼14-20% higher probability of fault orientation following stress conditions (7.7% fault likelihood along stress directions as opposed to ∼6.4-6.8% of fault likelihood in other directions [Figure 9e]).In the area containing Lunde Pockmark (i.e., the area within the red rectangle in Figure 9a), there is a ∼32% higher probability of fault orientation following stress conditions (7.9% fault likelihood along stress directions as opposed to ∼6.0% of fault likelihood in other directions) (Figure 9e).In the area northeast of Lunde Pockmark (shown by the yellow rectangle in Figure 9a), stresses do not favor formation of any preferential fault orientation potentially due to undetectable differences in maximum and minimum horizontal stresses.The absence of predominant fault orientation in the structural attribute analysis and the absence of the occurrence of symmetry axis in the area support our interpretation.
The gas leakage distribution within the pockmark coincides with inferred changes in the subsurface anisotropic conditions (Figure 9).We speculate that leakage points may relieve pressure and form zones where the sub-seabed stress field is locally modified (i.e., explaining the contrasting fracturing and faulting patterns to the northeast and southwest of Lunde Pockmark; Figure 9).Ultimately, the relationship between stresses, faulting and leakage may be the other way around; potential changes in the differential stress conditions generated by regional processes (e.g., glacial isostasy, mid-ocean ridge spreading, differential compaction; Fejerskov & Lindholm, 2000;Lindholm et al., 2000;Olesen et al., 2013, Vachon et al., 2022) may trigger leakage at this location by modifying the permeability of already in place faults and fractures.
The fracture/crack-induced permeability is important, considering hemipelagic sediments with low primary permeability in the region (Eiken & Hinz, 1993;Vogt et al., 1994).Vertically connected fractures increase permeability of a medium, promote leakage, and create anisotropy as inferred from the results we document here.In shallow sediments with a typical normal-faulting regime, fracture and crack planes oriented perpendicular to maximum stress can be closed due to stress (Nur & Simmons, 1969).Previous studies document that some faults and fractures are permeable to fluid flow whereas other faults and fractures are closed to fluid flow in the region (Plaza-Faverola et al., 2015;Plaza-Faverola & Keiding, 2019;Singhroha et al., 2016Singhroha et al., , 2019Singhroha et al., , 2020)).The relative decrease in faults detected along azimuths of 43°-47° in the area southwest of Lunde Pockmark (shown by the blue rectangle in Figure 9a) could be due to the presence of closed faults oriented perpendicular to the present day maximum horizontal stress (blue curve in Figure 9e).There is a relatively high occurrence of fractures and faults, (>7% fault likelihood in azimuths of 33°-51° as shown by the red curve in Figure 9e) in the direction perpendicular to the maximum horizontal stress, in the active seep area (shown by the red rectangle in Figure 9a).In the direction perpendicular to the maximum horizontal stress, fault likelihood decreases in the area without observed gas leakage (shown by the blue rectangle in Figure 9a), whereas it increases in the area with the observed gas leakage (shown by the red rectangle in Figure 9a).This finding suggests that faults are open due to a high gas pressure in the system.
With the available data, we are not able to resolve potential sources of local stress fields.However, a similar OBS study conducted ∼40-45 km south from Lunde Pockmark (north of the KR and MTF intersection; Figures 1a  and 1b) concludes that gravitational forcing may be dominating the stress regime within the upper 100 m (Haacke et al., 2009;Haacke & Westbrook, 2006).They document a transition in the orientation of maximum and minimum horizontal stresses at ca. 200 m and argue that regional tectonic stress may be overcoming the effect of gravitational stress at these depths (Haacke et al., 2009;Haacke & Westbrook, 2006).Our approach on Vestnesa Ridge only allows us to infer stress orientations within the upper 200 m.In line with the study by Haacke and Westbrook (2006) we think it is likely that the inferred orientations of minimum and maximum principal stresses on Vestnesa Ridge reflect the effect of gravitational stress acting beneath the sub-surface given that the Vestnesa Ridge crest is NW-SE oriented and has a southwest-dipping slope at its steepest flank (i.e., its southern flank) (Figures 1a, 1b, 9a, and 9b).The inferred sub-seabed minimum and maximum stress orientations based on amplitude null symmetries are consistent with the Vestnesa Ridge geometry.

Conclusions
We analyze active-source PS-wave energy recorded by OBS recorders (OBSs) from a seabed pockmark on Vestnesa Ridge offshore west-Svalbard.We document amplitude null symmetries in the transverse component from different OBS stations and interpret them as an evidence of anisotropy.Simultaneously, we analyze the likelihood of forming faults and fractures at specific azimuths via implementation of structural attributes on high-resolution 3D P-cable seismic data.
The combined interpretation from both analyses suggests the presence of a minimum horizontal stress exerted along the ∼45°-225° (NE-SW) direction and a maximum horizontal stress exerted along the ∼135°-315° (NW-SE) direction within the upper 100 m of sedimentary column.A change in the fracture pattern together with a lack of amplitude nulls symmetries to the northeast of the investigated pockmark indicate a drastic change in the local stress regime coinciding with the zone where seepage is concentrated at present day.The change in the differential stress across the seep site could have triggered seepage or conversely, the occurrence of seepage could be responsible for the observed stress differences inside and outside the pockmark region.The change in stresses across the seep site establish a correlation between the occurrence of seepage and the stress regime in shallow sediments in this region which favors (∼32% higher probability) the occurrence of faults perpendicular to the minimum horizontal stress.
The sources of sub-seabed stress cannot be resolved with the implemented approach.However, a close correlation of inferred principal stress orientations with the geometry of the Vestnesa Ridge suggests a dominance by gravitational forcing in line with similar OBS studies farther south on the same contourite system.Future analyses of the propagation of PS-wave energy from deeper layers may help to identify vertical changes in the stress field orientation potentially indicating the effect of glacio-tectonics on sedimentary deformation and seepage evolution in this deep marine Arctic setting.

Appendix A: Fault Analysis
We use structural seismic attributes (i.e., Petrel's variance and ant tracking attributes) derived from high-resolution 3D seismic data (Plaza-Faverola et al., 2015) to study fault orientations in our study area.The seismic variance attribute measures the difference in the shape of waveforms of adjacent traces, over a given lateral and vertical window (Van Bemmel & Pepper, 2000).This attribute detects discontinuities, faults, and changes in lithology.The ant tracking method, an even more powerful automatic tool, by using an analog of swarm intelligence inspired by the behavior of ants, detects faults in the seismic data with a very high level of accuracy (Pedersen et al., 2002).We use these two methods to estimate the likelihood of the occurrence of faults in different orientations (details on parameters used for attribute analysis available in Text S3 in Supporting Information S1).In order to do this, we employ a statistical measure to identify predominant fault and fracture directions (Figure A1).Fault analysis is performed on three subsets of data within the 3D seismic data, published and interpreted in Plaza-Faverola et al., 2015 (blue, red and yellow rectangles [shown in Figure 9a]).The area bounded by the blue rectangle covers the undeformed strata to the southwest of Lunde Pockmark and has data within inlines: 50-110 and crosslines: 700-800.The area bounded by the red rectangle contains Lunde Pockmark and active venting sites and has data within inlines: 110-166 and crosslines: 700-800.The area bounded by the yellow rectangle comprises undeformed strata to the northeast of Lunde Pockmark and has data within inlines: 166-226 and crosslines: 700-800.
After estimating variance and ant tracking attributes, we normalize values in the attribute volume by assigning the maximum value to one and minimum value to zero.We calculate the mean of attribute values along different azimuths in subsets of the data, for example, for data bounded by inlines: 50-110 and crosslines: 700-800, we calculate the mean attribute value along all azimuths from a reference common depth point (CDP) bin corresponding to inline 80 and crossline 750 (Figure A1), that is, if 20 CDP bins occur along a certain azimuth, we calculate the mean value of the normalized attribute values, corresponding to traces from 20 CDP bins.where, M(A(X), y) is the mean attribute value along azimuth of X° (A(X)) from a reference bin y.NAV X (t, bin) is the normalized attribute value at time t in a bin b, where bin b has an azimuth of X° or (X° + 180°) from the reference bin.N is the total number of bins that have azimuth of X° or (X° + 180°) from the reference bin and NT is the total number of samples between t1 (seafloor) and t2 (BSR).
We repeat the same procedure while shifting the reference CDP bins to other places (e.g., shifting the reference CDP bin to inline 81 and crossline 750) and estimating mean attribute values there (Figure A1).Subsequently, we derive fault likelihood as a mean value taken from all reference bins (e.g., if 0.01, 0.02, and 0.06 are mean ant tracking attribute values along azimuth of 30° from three different reference bins, fault likelihood along azimuth of 30° is 0.03% or 3%).where, fault likelihood (X) (in %) is the likelihood of fault occurrence within a fault strike azimuth of X° (or fault strike orientation of X° − (X° + 180°)).K is the total number of reference bins.
The averaging ensures that we statistically capture fault orientations even with very small variations in likelihood.Azimuths with 180° spacing are considered the same for the purpose of fault analysis in this study (e.g., fault strikes of 10° and 190° have the same fault strike orientation [10°-190°] and are together referred to as a fault with 10° azimuth [or fault azimuth of 10°] in this study).

Figure 1 .
Figure 1.(a) Regional map showing the location of the study area along with Knipovich Ridge (KR), Molloy Transform Fault (MTF), Molloy Ridge (MR), Spitsbergen Transform Fault (STF).(b) Bathymetry map showing various pockmark features along Vestnesa Ridge (ESVR, eastern segment of Vestnesa Ridge and WSVR, western segment of Vestnesa Ridge) and the location of Lunde Pockmark on Vestnesa Ridge.(c) Relocated OBS positions (shown in white circles) are laid over the seafloor map picked from 3-D seismic data and high-resolution bathymetry map of Lunde Pockmark.Dotted lines show location of inline and crosslines seismic sections in Figure 2. (d) Location of different shots fired using seismic source and red points show the relocated OBS positions.The green line shows shot locations for data shown in Figure 3.

Figure 2 .
Figure 2. Inline and crossline seismic sections imaging the Lunde Pockmark (inside the rectangle).Inline and crossline locations are shown by black dotted lines in Figure 1c.Locations at the top (marked in green) show the intersection of the inline and crossline.

Figure 3 .
Figure 3. Seismic data recorded by vertical (a) and horizontal (b, c) components of OBS13 instrument from shot numbers 9757-10,190 (shot locations plotted in green in Figure 1d).

Figure 4 .
Figure 4. Vector fidelity plots for OBS8 (a) and OBS3 (b) show the first S-wave arrival polarization azimuth and OBS locations (blue dots).The center of plotted vectors (a, b) is the position of shot and dotted radial lines from the OBS locations show the radial direction from both OBS8 (a) and OBS3 (b).The histogram plots for OBS8 (c) and OBS3 (d) show the distribution of the difference angle (corrected first S-wave arrival polarization azimuth-radial azimuth; more details in Text S2 in Supporting Information S1) for different shots.

Figure 5 .
Figure 5. Schematic of the radial (shown using the red color) and transverse (shown using the green color) polarization directions for a near-vertically reflected converted PS-wave (shown using the blue color).

Figure 6 .
Figure 6.Radial and transverse components of recorded S-waves for OBS6-9.Black ellipsoids show regions of low energy in the transverse components from S-wave splitting.

Figure 7 .
Figure 7. Radial and transverse components of recorded S-waves for OBS11, OBS13, and OBS15.Black ellipsoids show regions of low energy in the transverse component from S-wave splitting.

Figure 8 .
Figure 8. Radial and transverse components of recorded S-waves for OBS17 and OBS22.

Figure 9 .
Figure 9. (a) Maximum (Sh max ) and minimum (Sh min ) horizontal stress orientations plotted over the seafloor map from 3D seismic and high-resolution bathymetry.(b) Conceptual diagram of a gas chimney, free gas flow and the base of gas hydrate stability zone (BGHSZ) laid over the seafloor path-profile along the red line (partly plotted in (8a) and full length plotted as red line in Figure 1b).(c) Seismic section along crossline 750 illustrating gas chimney features and OBS locations with respect to gas chimney.(d) Inline 70 and crossline 750 plotted along with ant tracking-attribute time-slice.(e) Likelihood of occurrence of fault in different orientations (10° fault azimuth corresponds to fault strike orientation of 10°-190°) based on the fault extraction method using ant tracking in the areas shown by colored rectangles.Blue, red and yellow rectangles have data within crosslines: 700-800 and inlines: 50-110, 110-166, and 166-226, respectively.(f) Inline 70 and crossline 750 plotted along with variance attribute time-slice.(g) Mean variance amplitude in different orientations (data sets shown with colored rectangles in (a)).

Figure 10 .
Figure 10.The ratio of radial and transverse energies in the 1.1-2.4s arrival time interval.

Figure 11 .
Figure 11.The ratio of radial and transverse energies in the 1.1-2.4s arrival time interval.

Figure A1 .
Figure A1.Ant tracking-attribute time-slice at 1.692 s two-way time.A(X) represents azimuth of X° from a given point or reference bin.A(X) and A(X + 180) are considered as A(X) for the purpose of fault analysis.Blue lines, red lines, green lines and black lines show data points in the attribute time-slice along azimuths of 0°, 46°, 90°, and 136°, respectively from two different reference bins.M(A(X), y) in equations in Appendix A represents mean of all normalized attribute values (NAV) along a A(X)-My line.The A(X)-My line has azimuth of X° from a reference bin y (RB y ) (e.g., value of A(46)-M1 will be mean of all values along the red line passing through RB 1 [the red line lies along azimuth of 46° from RB1]).

Table 1
Analysis of Recorded Data Set From Different OBS Stations Note.The pale green, pale orange and pale red color cells signify acceptable, medium and poor data quality, respectively.