Seismic Ambient Noise Imaging of a Quasi-Amagmatic Ultra-Slow Spreading Ridge

: Passive seismic interferometry has become very popular in recent years in exploration geophysics. However, it has not been widely applied in marine exploration. The purpose of this study is to investigate the internal structure of a quasi-amagmatic portion of the Southwest Indian Ridge by interferometry and to examine the performance and reliability of interferometry in marine explorations. To reach this goal, continuous vertical component recordings from 43 ocean bottom seismometers were analyzed. The recorded signals from 200 station pairs were cross-correlated in the frequency domain. The Bessel function method was applied to extract phase–velocity dispersion curves from the zero crossings of the cross-correlations. An average of all the dispersion curves was estimated in a period band 1–10 s and inverted through a conditional neighborhood algorithm which led to the ﬁnal 1D S-wave velocity model of the crust and upper mantle. The obtained S-wave velocity model is in good agreement with previous geological and geophysical studies in the region and also in similar areas. We ﬁnd an average crustal thickness of 7 km with a shallow layer of low shear velocities and high Vp/Vs ratio. We infer that the uppermost 2 km are highly porous and may be strongly serpentinized.


Study Area and Motivation
The Southwest Indian Ridge (SWIR) is an ultra-slow spreading ridge that separates Africa and Antarctica. It reaches from the Bouvet triple junction (BTJ) in the southern Atlantic Ocean to the Rodrigues triple junction (RTJ) in the Indian Ocean ( Figure 1A). The mean spreading rate of this ridge is about 14 mm/year [1]. Our study area is located in the easternmost part of the SWIR, centered at 64°35 E and 28°S. This easternmost portion of the ridge comprises its deepest parts reaching more than 5500 m below sea level. Compared to other mid-oceanic ridges, it is relatively melt-poor [2].
Previous surveys revealed the apparent absence of volcanic activities both at segment ends and along the ridge axis for over 100 km [3,4]. A large amount of mantle peridotite was detected in the axial valley [3,5]. The considerable quantity of outcrops of serpentinized peridotite on the seafloor is ascribed to exhumation in the footwall of detachment faults known as Oceanic Core Complexes (OCC) [6]. Corrugated OCCs appear in more magmatically active regions of the SWIR and can be identified by spreading-parallel corrugations. OCCs also occur in amagmatic and smooth terrains of the SWIR, which are named non-corrugated OCCs [7,8]. Seafloor sampling in corrugated OCCs shows a variety of petrologic compositions including peridotites, gabbros, and volcanic types [9]. The dredging of the non-corrugated OCCs exhibited 90% serpentinized peridotites with a minor contribution of gabbros [8]. The exposed part of the detachment fault system related to these OCCs is divided into emergence and breakaway ( Figure 1B).
Sampling and direct geological observations were limited to the near-surface seafloor and could not define the extension and depth distribution of serpentinization, magmatic rocks, and detachment faults. Therefore, geophysical methods were recruited to better understand the processes involved in shaping this melt-poor, ultra-slow spreading ridge. Limited off-axis geophysical data in the late 1990s indicated significant crustal thickness variations across and along the ridge axis in the eastern part of the SWIR. These significant variations suggested that the melt supply is spatially more concentrated and short-lived, in comparison to faster spreading ridges [10]. The first comprehensive off-axis data set in the easternmost part of the SWIR was collected in 2003 by bathymetric, magnetic, and gravimetric methods, covering an area of up to 250 km away from the ridge axis on both sides [7]. It was found that 37% of the surveyed seafloor showed only little or no evidence of volcanic activity. This non-volcanic seafloor is termed smooth seafloor because of its rounded topography. It has no equivalent at faster spreading ridges [7]. Four percent of the mapped region was identified as corrugated OCCs [11] and the remaining 59% was shown to be formed by volcanic activity [1]. The melt supply appears to be focused beneath the volcanic centers. This focusing has been explained with melt originating near the base of the lithosphere and rapid extraction via dikes [3,10,[12][13][14]. Off-axis studies showed higher values of magnetic anomaly over volcanic seafloor areas and lower magnetization over smooth seafloor or non-volcanic regions. This weakly induced magnetization in nonvolcanic areas of the seafloor is supposed to be produced by serpentinized peridotites [1,15]. A similar separation can be seen from maps of residual mantle Bouguer gravity anomaly where lower values (negative) were found for volcanic seafloor and higher values (positive) for deep sections or non-volcanic regions [10,14].
Seismic reflection profiles near 64°E suggested reflectors dipping 45-55°down to 5 km below the seafloor beneath the emergence of the axial detachment fault. These reflectors are interpreted as tectonic damage zone in the basement [6,12]. Next to these damage zones, the thickness of the crust was greater than further away from it. The existence of a low Vp seismic crust of a total thickness of up to 5 km , defined as Vp < 7.5 km/s [6], is explained by a downward gradient of serpentinization, with small magmatic intrusions. Several dipping reflectors were imaged in the hanging wall of the axial detachment fault and interpreted as small offset faults and trapped volcanic rocks emplaced in the ultramafic basement [12]. East of our study area, at 66°E, Minshull et al. [16] report crustal thicknesses between 2.2-5.4 km (avg. 4.2).
Imaging the crustal structure is necessary to understand past and ongoing tectonic processes and provides an important prior to mantle tomographic studies. In this study, we investigate the robustness of the ambient-noise method in a marine exploration setting. The obtained subsurface shear-velocity model will be compared to previous studies and interpreted in the context of the regional tectonic evolution. . The thin black line shows SWIR axis. Variation of fracture zones orientation represents the variation of spreading direction with time [17]. The blue rectangle marks our study area. (B) Bathymetry map of the study area. The locations of ocean bottom seismometers are shown by reverse triangles. The final stations used for inversion are indicated by black color and the removed stations based on data analysis are shown by gray color. Most of the stations have been placed on smooth seafloor. The red thick line is the axis and the purple line is indicating the emergence of the axial detachment fault system. The breakaway of the active fault is located at the top of the northern axial valley and suggested breakaways of inactive faults are seen in the southern Antarctic plate.

Ambient Noise Interferometry
Seismic ambient noise is a relatively continuous signal generated by the coupling between ocean waves and solid earth [18,19]. It is dominated by surface waves, with a period-dependent ratio of Rayleigh-to Love-wave energy [20]. Cross-correlating and stacking continuous ambient-noise records of arbitrarily chosen station pairs yields the empirical inter-station Green's function which can be used to extract phase-velocity information [21]. Often, only the Rayleigh component of the surface waves is used because it is more energetic and can be measured on the vertical receiver component [22]. Rayleigh wave phase-velocity measurements are most sensitive to the subsurface Vs structure with a minor dependence on Vp and density. The resulting ambient noise models yield complementary information to more classical teleseismic earthquake tomography and active source experiments, with each of the methods being sensitive to a different frequency range and having a different ray path coverage [23]. Additionally, active source methods have become less popular due to the potentially detrimental effects on the marine fauna.
Passive seismic interferometry has been widely utilized to image the crust and upper mantle in continental areas [24][25][26][27]. There are significantly fewer studies in oceanic regions, mainly because of limited availability of seismic networks on the ocean bottom and high levels of local noise, e.g., instrumental noise which is unusable as opposed to the ambientnoise signal [28,29]. In the next steps, we will explain theoretical concepts of measurement, processing, and imaging with ambient noise data. Furthermore, we illustrate how we have applied the concepts to get the final S-wave velocity model.

Measurement of Ambient-Noise Data
The marine seismic data used in this study were acquired by 43 Ocean Bottom Seismometers (OBS), each having three components geophones and one hydrophone. To analyze the ambient-noise signals, we used only the vertical component of the geophones [29,30]. The data acquisition took place between 25 September and 30 October 2014 during the SISMOSMOOTH cruise [31,32]. As the instruments were from three different countries, we have labeled the stations Canadian (C), Taiwanese (T), and French (F) ( Figure 1B). Of the 43 OBS stations used in this study, 29 were Canadian, 7 Taiwanese, and 7 French. The sampling rate was 250 Hz. Based on the furthest inter-station distances, the east-west and north-south aperture of the array was 56 km and 19 km, respectively. The spacing between the OBS stations was 2.5-5 km, and their depth ranged between 3.9 and 5.01 km. The raw records have been split into their components and corrected for clock drift and relocation errors.

Probability Density Function (PDF)
We did a power spectral density (PSD) estimation of the recorded signals and their probability density function (PDF) at each station [33]. While this analysis is expected to be dominated by ambient-noise sources, it includes also earthquakes and undesired signals, i.e., instrumental glitches, body waves or signals that do not propagate between stations pairs, that may deteriorate the ambient-noise measurements. We use this information to (1) estimate the intensity of the ambient-noise signals at different frequencies, (2) estimate the portion of undesired signals, and (3) check the overall health of the instrument.
The PSD curves were estimated in terms of decibel at each frequency via direct Fourier transform introduced in [34]. PDFs were constructed at each frequency distribution bin by gathering PSDs in 1/8 octave of period intervals and 1 decibel of power intervals. We have shown 4 exemplary PDF plots in Figure 2. Each figure is for the vertical component of one station and during the whole time of recording. The PSD curves with higher probability distribution were interpreted as the ambient noise whereas the scattered PSD curves can be attributed to undesired signals [35]. The analysis indicates that the recorded energy is low at frequencies below~0.1 Hz. Between 0.1 and 3 Hz, the curves reach a flat level showing low amplitude variations and a high probability level. At frequencies above 3 Hz, the variance between PSD curves increases significantly, which represent likely a result of random, local noises [33,36]. We used the PDF analysis in the junction with an inspection of the cross-correlation curves shown in the next section to limit the range of frequencies used in this study. This reduction of the frequency range speeds up the processing.

Cross-Correlation
Of the 43 available stations, we made sure to only correlate station data with identical instrumentation to avoid any potential phase bias, as we did not have access to the instrument response data. This left us with 200 possible station pairs. The signals recorded at each station were cut into time windows of 10 min with an overlap of 60 percent between successive windows. More complete explanations about the effects of overlapping time windows can be found in [37]. We did not explicitly remove earthquake signals which are expected to have only a minor effect on the obtained phase velocities [38]. However, we applied spectral whitening [39], to equalize spectral amplitudes and down-weight the effects of earthquake signals and monochromatic sources. The whitened 10-min segments were cross-correlated between station pairs and stacked. The cross-correlation traces show forward and backward propagating signals (causal/acausal) from opposing directions along the inter-station path, visible at positive and negative lag times in Figure 3. At inter-station distances greater 28 km, we did not observe any clear Rayleigh-wave signals whereas at distances below 5 km causal and acausal signals are not well-separated. This limitation is similar to the criteria of inter-station threshold recruited in the literature (see, e.g., in [36,39]) and resulted in 176 acceptable cross-correlations used in this study ( Figure 3). The timedomain cross-correlations can be used to inspect the distribution of ambient-seismic-noise sources. Perfect symmetry between the cross-correlations at positive and negative lag times is an indicator of a homogeneous azimuthal noise-source distribution [28,29]. In theory, the waves moving in opposite directions should be identical between station pairs [40]. An asymmetric cross-correlation thus indicates a low signal-to-noise ratio.  Figure 3 with regard to the inter-station distances on the vertical axis. We can discern the symmetry between positive and negative time lags in the frequency range of less than 0.8 Hz. The arrival of the Rayleigh wave signal cannot clearly be detected at higher frequencies (0.8-1.2 Hz). This loss means that the seismic noise sources are too weak and/or the random signals level is too high for a recording time of one month. Based on this analysis and the PDFs shown in the previous section, we limited our study to the frequency range between 0.05 and 1 Hz. This frequency range is in agreement with oceanic microseism utilized also in other seismic noise tomographic studies and typically generated by the non-linear interaction between the ocean and the solid earth (primary microseisms 0.07 Hz, secondary microseisms 0.14 Hz [18]). The amplitude differences between causal and acausal parts of the cross-correlations in Figure 4 are low, indicating a relatively homogeneous source distribution. This is supported by the more detailed comparison provided in Figure 4.

Amplitude Ratio
The amplitude ratios can be used to estimate the homogeneity of the noise-source distribution [41]. The amplitude ratio is defined as the maximum amplitude at the positive lag time of the cross-correlations over the maximum amplitude at negative lag time.
The rose diagram in Figure 5 shows amplitude ratios of all station pairs in our array. This diagram is a proxy for the azimuthal distribution of noise sources. If two microseisms propagate along the inter-station path in opposite directions and with the same amplitude, the amplitude ratio will be one. If the ambient-noise energy coming from one direction is larger, we will observe a dominant direction in the rose diagram. A comparison of Figure 5A-D suggests that the noise-source distribution is more uniform at relatively high ( Figure 5D) compared to relatively low ( Figure 5A) frequencies. However, apart from some peaks in the rose diagrams, we find no clear preferential noise-source direction. In other words, the ambient noise sources are approximately equi-distributed within the studied frequency range which is another endorsement for the quality of time-domain cross-correlations.

Phase Velocity Determination
The cross-correlations can be used to obtain the empirical Green's functions and were analyzed to gain information on the subsurface structure between the receiver pairs [42]. Most authors either picked group velocities by identifying the maximum in the envelope of surface-wave packets [43][44][45] or measured the phase velocity dispersion [46,47]. Fewer authors preferred utilizing the full waveform in their analysis [48]. In this study, we extracted phase velocities, because of the advantages in terms of accuracy, depth of investigation, and less contamination by interfering phases [42,49,50].
Several methods of phase-velocity measurements have been introduced such as the image transformation technique in the time-domain [51], the Frequency-wave number (F-K) method [52,53], multiple filter analysis [54], and Bessel function methods based on the theory of Aki [55] such as direct fitting in the frequency domain [56] or the herein applied zero-crossing fitting [38,47,49]. The intrinsic limitation of these methods is poor performance at low levels of signal-to-noise ratio. Plus, if the inter-station distance is too small, only few zero-crossings can be picked along the frequency axis, thus degrading the resolution of this method [56]. To mitigate these limitations, as explained above, we have only cross-correlated the station pairs distant by more than 5 km and also we checked the homogeneous distribution of sources in the data analysis ( Figure 5D).
We used the previously computed cross-correlations and extract their zero crossings in the frequency domain ( Figure 6). The zero crossings can be related to the phase velocity by comparison with the known zero crossings of a Bessel function [47,55]. We apply a smoothing filter to the cross-correlations which removes spurious zero crossings and thus stabilizes the obtained phase-velocity curves [49]. The phase-velocity picks from the 176 station-station cross-correlations were then taken manually. The picking procedure was performed in two iterations: after picking the phase-velocities for all station pairs in the first round, we calculated an average dispersion curve that was used to guide the manual picking in the second round. We applied several quality criteria: as explained above, we only accepted station pairs with an inter-station distance between 5 and 28 km. Because of their low signal-to-noise ratio, the French station data were not used. In addition, cross-correlations with large imaginary parts were discarded.  [49]. According to the data analysis, the maximum range of frequency set to 1 Hz.
The average dispersion curve for 42 station pairs is shown in Figure 7. As expected from the noise analysis, only a few measurements at periods higher than 10 s could be taken. Therefore, we did not take any measurements above 10 s period into account in the depth inversion. To see the variations of phase velocities, the dispersion measurements were plotted at four periods (1.2, 2.5, 3.1, and 6 s) in Figure 8.  An increasing trend of the measured phase velocity from shorter to longer periods is observable.

Inversion of Phase-Velocity Dispersion Curves
We made use of the dispersive nature of Rayleigh waves to obtain a 1D depthdependent shear-velocity profile. We used the average phase-velocity curve shown in Figure 7, thus neglecting any lateral variations of the subsurface structure. The dispersion curve can be inverted with linearized methods [57][58][59]. However, we preferred applying a Bayesian approach because of the highly nonlinear relation between phase velocity and subsurface structure and the subsequent non-uniqueness of the problem. We used the conditional neighborhood algorithm implemented in the Dinver software [60]. As in other Monte Carlo methods, we had to pre-define a range of priors for the pseudo-random sampling approach ( Table 1). The range of statistical speculation and the role of prior information in this method are discussed in articles such as Scales and Tenorio [61] and Mosegaard and Tarantola [62]. Compared to other Monte Carlo approaches, the conditional neighborhood algorithm is more self-adaptive when searching through the parameter space. The search starts with a pre-defined number (n s0 ) of models distributed with a uniform probability in parameter space. The misfit function (L 2 norm) is calculated for every model and the n r models with the lowest misfits are selected for the next iterations. This is iterated n t times, and at each iteration the selected lowest-misfit models are decomposed into n s /n r new models. The data misfit is given in relative units, normalized by the input data. In the conditional neighborhood algorithm, the dependency between parameters is taken into account at the beginning of the inversion by using a variable transformation [60].
The model is parameterized in terms of thickness, shear wave velocity (Vs), compressional wave velocity (Vp), and density (ρ) for each layer. Among them, Vs has the most important role and ρ exercises the least influence on the dispersion curve. The Poisson's ratio, which links Vs and Vp, was allowed to vary between 0.2 and 0.5 as indicated in Table 1, according to the value of Poisson's ratio for geological materials [60]. Vp and Vs are free parameters of the inversion. The Poisson's ratio is not a free parameter, but if a solution is generated with Vp and Vs such that the value of the Poisson's ratio is outside of the allowed range, that solution is discarded. We chose the following model search parameters: n s0 = 50, n r = 50 and generated 10,000 models in total. We chose to use four layers whose prior ranges (Table 1) were defined based on the Preliminary Reference Earth Model (PREM) [63] and on several test runs in which we made sure that the search algorithm would generate a reasonably wide range of models centered around the measured data ( Figure 9C). The obtained P-wave and S-wave velocity models are shown in Figure 9A,B, respectively. The black lines represent the average of the 5000 models with the lowest misfit and will be used for interpretation. The shown Vs profiles are more reliable, because of the higher sensitivity of Rayleigh waves to Vs compared to Vp.

Results and Discussion
We limit the discussion to the velocity models obtained from phase velocities up to a period of 10 s. This constraint allows investigating the subsurface structures down to~15 km. In general, we consider our average Rayleigh-wave dispersion curve to be robust. Because Rayleigh waves are mostly sensitive to Vs, we are also confident that our Vs models are reasonably well constrained. Rayleigh waves are well known to be much less sensitive to Vp than Vs (see Figure 10), so our inferences are Vp and the Poisson's ratio should be taken with care. We shall show, however, that the values we obtained are reasonable based on independent knowledge of those parameters. We also provide in Figure 9A,B a rough estimate of model uncertainty, defining as our confidence interval the range of values from three standard deviations below to three above the average. This definition of confidence is complicated near interfaces, where uncertainties in seismic velocity and interface depth are entangled with each other. Figure 10. Sensitivity kernels relating Rayleigh-wave phase velocity estimates to Vs (left) and Vp (right) at depth. Note the difference in horizontal scale between the left and right panels: sensitivity of Rayleigh waves to Vp is one-order-of-magnitude smaller than their sensitivity to Vs, as expected. The reference assumed for the computation of the kernels is our final Vp, Vs profile (solid black lines in Figure 9A,B). Kernels were calculated via Keurfon Luu's "disba" Python package.
We first compare our average dispersion curve ( Figure 9C) with results from ambientnoise tomography obtained in the south central Pacific (Figure 6 in [28]) and in the transform fault region of East Pacific Rise (Figure 4 in [29]). These two areas were selected as two of the rare examples where ambient-noise interferometry has been used to image the oceanic lithosphere. While our study areas are different, this qualitative comparison provided us with some insights from similar tectonic settings. In the average dispersion curve in Figures 7 and 9C, we can identify an increasing trend of phase velocity up to a period of 20 s. This trend is similar to that found in the two other oceanic studies. Up to 16 s, both other regions have a phase velocity range of approximately 1.5-3.7 km/s while this range in our case is approximately 1-3.5 km/s. Lack of data in periods less than 2 s in Harmon et al. [28] probably account for the different minimum phase velocity.
Similar to our upper period bound of 10 s, Yao et al. [29] and Harmon et al. [28] mentioned that their phase velocity estimations for the fundamental mode of Rayleigh wave were subject to high uncertainty at periods greater than 8 s and 9 s, respectively. Nevertheless, longer recording time (200 days in [28] and over 1 year in [29]) enabled them to benefit from ambient noise data up to periods of 16 s and 30 s. A longer recording time will usually lead to a smoother distribution of incoming ambient-noise energy and suppresses errors in the dispersion-curve measurements [64].
The P-wave velocity in the crustal layer is expected to be lower than the velocity for olivine-rich gabbros, which is about 7 km/s [65] which we use to approximate the the Moho depth. Therefore, we infer an approximate crustal thickness of 7.2 km from Figure 9A. A similar conclusion can be drawn from the S-wave velocity model in Figure 9B, where the maximum expected crustal S-wave velocity is approx. 3.9 km/s [65]. This thickness is larger than the 4.2-5 km estimated from the P-wave velocity model in Momoh et al. [12] in the central part of the study area. The larger estimate of crustal thickness resulting from our approach could be due to the location of most used stations in the hanging wall compartment of the south facing active axial detachment fault (Figure 1). Both P-wave and S-wave velocity models represent almost linear and monotonic velocity to depth variations in Figure 9A,B.
As shown in Figure 9, the crustal range of S-wave velocities is approximately 1.1-3.9 km/s. It varies between 2.5 and 7 km/s for the P-wave velocity. We calculated the Vp/Vs ratio at different depths from our best Vp and Vs models (Figure 9). The comparison between the minimum and maximum values of S and P velocities gives a Vp/Vs ratio equal to 2.23 and 1.81, respectively ( Figure 11). The ratio of P and S wave velocities provides an additional constraint on the the lithology [66] and is sensitive to the presence of fluids [67]. The high ratio of Vp/Vs = 2.23 at a shallower depth of the crust is indicative of high porosity and some amount of serpentinized peridotite. The minimum ratio of 1.81 is predicting lower porosity at deeper parts and is in agreement with Vp/Vs ratio of SWIR gabbroic rocks analyzed in Miller and Christensen [65]. Fast variations of P-wave and S-wave velocities at depths below 6.5-7.2 km and slow changes at depths >7.2 km can be another sign of permeability and porosity regime. This regime can affect the serpentinization process in the crust and upper mantle. The variable range of Vp/Vs in the crust and considerable serpentinization can be explained by the fracturing of the crust associated with detachment faults systems [68]. Another evidence of the fracturing of the crust was discovered in south dipping reflectors by Momoh et al. [12] which was interpreted as the damage zone. This zone is located within 2-2.8 km horizontally from the emergence of the active axial detachment fault system (deep regions on the map in Figure 1). The suggested serpentinization expressed in the velocity models is in agreement with seafloor sampling in non-corrugated oceanic core complexes collected from smooth (non-volcanic) seafloor [8].  Figure 9A,B), as a function of depth. Dashed lines identify the confidence interval, defined as in Figure 10. Full distribution of sampled models with the lowest misfit are shown in color.
We also compare our P-wave model with other models of ultramafic crusts. The Pwave velocity model of ultramafic crust at the Mid-Atlantic ridge at 23°N [69] is more linear and the observed velocity of the uppermost crust is higher than in our study, reaching the mantle faster. This higher range of velocity may be caused by gabbroic intrusions or fewer deformations of the exhumed, ultramafic crust in that region. The magma-poor rifted margins in Iberia-Newfoundland [70] and in western Iberia [71] also demonstrated a higher range of P-wave velocity in the upper ultramafic crust. Sedimentation and older ages of these margins likely explain why we observe this higher range of velocity.
An S-wave velocity higher than 3.9-4 km/s is typical for mantle rocks and our model is therefore indicating that unserpentinized and unfractured mantle is found at depths greater than 7 km. The average shear wave velocity between the base of the crust (defined by a 4 km/s Vs threshold), and the maximum depth of our model (15 km) is 4.2 km/s which is less than the global reference value of 4.5 km/s [72], that concerns old, off-axis oceanic lithosphere, but close to the 4.29 km/s of young oceanic lithosphere in the Pacific [29]. This relatively low value suggests that pristine and unfractured uppermost mantle may in fact be found deeper than 15 km below the hanging wall seafloor. This is consistent with results from a microseismicity study that shows earthquakes down to that depth in the investigated area (Jie et al., 2020 AGU abstract and paper soon to be submitted). The absence of a seismically well-defined and distinct Moho is consistent with reflectivity studies of Momoh et al. [12] in our study area. Moreover, other reflectivity studies predicted that more visible and recognizable Moho should be formed at faster spreading ridges, where mafic oceanic crust emplaced on the ultra-mafic mantle [73,74]. This prediction was also confirmed by another magma-poor ultra-slow spreading ridge by [75].

Conclusions and Outlook
1. We analyzed and cross-correlated ambient-noise recordings from 43 OBS stations. The symmetry of the cross-correlations in the time-domain and the amplitude ratio analysis showed that the noise source distribution is very homogeneous, despite the short operation time of the recorders.
2. The phase velocities extracted from the inter-station correlation functions ranged between 1 km/s to 3.5 km/s for periods between 1 s to 10 s and was in good agreement with other oceanic studies. Outside this period range, the low signal-to-noise ratio prohibited further analysis. Joint long time recording of ambient noise and earthquake data is proposed to cover also longer periods.
3. While the general P-wave velocity structure of our study is similar to previous works, we find a significantly thicker crust of~7 km, averaged over the study region. At uppermost mantle depths, our results suggest an average shear velocity around 4.2 km/s which is less than the global average. This can be explained by the younger age of the seafloor in our area. 4. Our model suggests a very low S-wave velocity of around 1-1.5 km/s at shallow depths above 1 km in junction with a high Vp/Vs ratio of approx. 2.2. We infer that this may be due to the high porosity of the shallow, water-rich uppermost layer. The Vp/Vs ratio gradually decreases to ∼1.8 at a depth of 2.5 km. We propose that this transition reflects the decrease in porosity and assume a higher degree of serpentinization of the crustal material at shallow depths.
5. The ambient seismic noise method can be used as a complementary technique in geophysical studies of the oceanic lithosphere. The S-wave velocity model obtained provides new constraints to evaluate the velocity structure of the crust and uppermost mantle at a melt-poor portion of the ultra-slow spreading SWIR.

Conflicts of Interest:
The authors declare no conflict of interest.