Abstract
This study assesses the variability of shear-wave (
Keywords
Introduction
Background
Shear-wave velocity (
The leading quantitative measure of a site is
The
According to Boore (2006), there are differences between the methods, and no single method is the best for all circumstances. Thus, the use of several different tests is worthwhile for providing constrained
Xia et al. (2002), compared MASW
The estimation of the
In this study, I focus on six sites (Figures 1 and 2). I estimate the uncertainties associated with measuring and interpreting the near-surface

Location map with the main active faults of the area (faults modified after Hamiel et al., 2018; Hofstetter et al., 1996; Kagan et al., 2011; Sharon et al., 2020).

Schematic columnar section of the six sites (H for Holocene).
Site conditions and geological setting
The BET site is located near the town of Kiryat Malakhi inside the industrial area of Be’er Tuvia (Figure 1 and Supplemental Figure A1). The data were acquired near the Local Council of Be’er Tuvia, over a thin asphalt layer on top of Holocene alluvium. The Holocene Dune Sands and the Pliocene Calcareous Sandstone Units are assumed to comprise the upper tens of meters beneath the surface (Sneh and Rosensaft, 2004) (Figure 2). The BSH site is located near the city of Beit Shean within the Shan factories zone in the Valley of Springs (Figure 1 and Supplemental Figure A4). The data were acquired near a cooling facility where a thin asphalt layer overlies several tens of meters of the Holocene alluvium, underlain by the Pliocene Cover Basalt Unit (Hatzor, 2000) (Figure 2). The CAS site is located near the community of Pardes Hana-Karkur inside the Industries Park Caesarea (Figure 1 and Supplemental Figure A7). The data were acquired over the Holocene Dune Sand Unit underlain by the Pliocene Red Sands & Loam Unit (Sneh and Rosensaft, 2004) (Figure 2). The HOD site is located within the town of Hod Hasharon (Figure 1 and Supplemental Figure A10). The data were acquired near the community’s Local Council over rich clay soil. The Pleistocene Red Sand & Loam Unit is assumed to comprise the uppermost 10 m beneath the surface (Sneh and Rosensaft, 2008) (Figure 2). The MGH site is located near the town of Migdal HaEmek inside the industrial zone (Figure 1 and Supplemental Figure A13). The data were acquired inside the Sunfrost factory over a thin asphalt layer on top of up to a meter of artificial fill. Weathered basalt is exposed under the filling material (Figure 2). The RAM site is located in the western part of the city of Ramle, adjacent to a water reservoir belonging to Mekorot (Israel’s national water company) (Figure 1 and Supplemental Figure A16). The data were acquired near the reservoir over a thin layer of artificial gravel on top of the Pliocene Red Sand and Loam Unit (Sneh and Rosensaft, 2004) (Figure 2).
Seismic methods
Surface-wave propagation occurs at free boundaries, such as the surface of the earth. In a layered medium, surface-wave propagation is ruled by geometric dispersion, in which longer wavelengths penetrate to greater depths. Therefore, at each frequency, the phase velocity depends on the properties of the examined interval below the subsurface (Garofalo et al., 2016a). The wave equation has multiple solutions at given frequencies or wavenumbers, which explains the multi-modal character of surface wave dispersion. Therefore, in a layered medium, the surface-wave propagation can be a multi-modal phenomenon, where different modes can exist at the same frequency but with different propagation velocities (Foti, 2000). The phase velocity depends primarily on the shear-wave velocity and is relatively insensitive to realistic variations in density and compressional wave velocity with depth (Nolet, 1981). Therefore, the surface wave velocity is a reasonable indicator of
Most surface-wave methods are based on three main steps: (1) acquisition of experimental data, (2) signal processing to obtain the experimental dispersion curve, and (3) one-dimensional (1D) inversion to estimate
In this research, I used MASW (Park et al., 1999) to acquire both the vertical and the radial components of Rayleigh waves and the transverse component of Love waves (by Dal Moro’s notation (2014): ZVF, RVF, and THF, respectively (Figure 3). To provide better estimations, I attempted to resolve the different modes. I also acquired data from the ESAC (Ohori et al., 2002) and the ReMi (Louie, 2001) methods. These surface-wave methods are commonly divided into two main groups (Table 1): active methods, where the source of the seismic waves can be a sledgehammer or an accelerated weight, versus passive methods in which the source is uncontrolled. Thus the waves are the result of spatially random sources.

Schematic plan of sources and geophones orientation to obtain different surface-wave components (modified after Dal Moro, 2014).
Acquisition parameters
ZVF: Rayleigh vertical component; RVF: Rayleigh radial component; THF: love transverse component; MASW: multi-channel analysis of surface waves; ReMi: refraction microtremor; ESAC: extended spatial autocorrelation.
To calculate which frequencies bias each site, and if there are preferential directions to the noise in the wavefield, and how these noise properties vary with time, I used the HVSR. HVSR is based on measuring two horizontal and one vertical component at a single station. By comparing the Fourier spectrum of the average horizontal components to its vertical component, the fundamental frequency of the site can be defined. Over the last decades, methods using a single station to record environmental noise have provided promising results in estimating HVSR (e.g. Arai and Tokimatsu, 2004; Bard, 1998; Mucciarelli and Gallipoli, 2001).
By acquiring HVSR data, we can also produce: (1) a time-dependent spectrogram; (2) power spectral density curves, which are the signal’s power content versus frequency; and (3) HVSR directivity.
Acquisition
In this research, both active and passive surface wave data were collected by adopting an “L”-shaped array (Dal Moro, 2014) (Figure 4). This shape allows us to acquire five different types of data sets using vertical geophones: ZVF and ReMi for each array line segment (a total of four data sets) and ESAC from both line segments. By changing to horizontal geophones (identical cable layout) adequately oriented for RVF (and then reoriented for THF), I acquire four more data sets for both line segments. Ideally, this “L” array provides nine different types of data sets with minimum logistics. In practice, not all nine data sets were recorded due to limitations in time, logistics, and equipment. I used R.T. Clark’s 4.5 Hz horizontal and vertical geophones, and the data were recorded on a Geometrics Geode seismograph (Table 1). HVSR measurements were acquired for periods of 30–60 min by Gemini 2 Hz, three-component geophones.

“L”-shaped array (the Ramle site).
Active survey
The active surveys were designed to acquire both Rayleigh (ZVF and RVF) and Love (THF) wave data. Measurements were collected by a linear array of geophones (Table 1), uniformly spaced according to the resolution required, and the desired depth of investigation (Socco and Strobbia, 2004). For a seismic source, I used an 8-kg sledgehammer which is commonly used and might contribute to uncertainties in an urban environment (due to the stiffness of the top layers and the urban seismic activity, the data is less clear at greater depths). For the collection of Rayleigh wave data, I used a linear array of alternating vertical and horizontal geophones. The vertical component of Rayleigh waves (ZVF) was recorded with the vertical geophones, and the radial component (RVF) with horizontal geophones oriented parallel to their linear segment. In both cases, the hammer struck vertically on a 20-cm2 aluminum plate. The data logger is triggered by the closure of a circuit between the hammer and the plate. Then, the plate was moved to the next place of impact. Rotating each horizontal geophone by 90 degrees (to be oriented perpendicular to the trend of that survey line) and striking horizontally allows the acquisition of Love wave data (THF). The horizontal blows were produced by striking the end of a wooden beam, aligned perpendicular to each array, anchored by the weight of a field vehicle, or using the U.S Geological Survey’s shear-wave seismic source (Haines, 2007). For an acceptable signal-to-noise ratio, I used the so-called “vertical stacking” approach, which is a synchronized summation of multiple blows of the hammer (Table 1). I recorded the data for 0.5–2 s according to the site conditions. At each site, 3–6 offsets were selected for each data set. This gave the option of picking the dispersion curve over the clearest dispersion image.
Passive surveys
ReMi data sets were obtained passively, using the same linear array as was used for the vertical component of the Rayleigh wave (ZVF) (part of the active MASW method). ESAC data sets were simultaneously recorded on both arms of the “L”-shaped array, which included 24–48 vertical geophones. For both passive methods, field data were collected for much longer times (25–30 min) (Table 1).
Processing
The surface-wave analysis focuses on the interpretation of the experimental dispersion curve, which is used for the inversion process. The processing technique is based on the transformation of the experimental data from the time-space domain to the velocity-frequency domain. The different data sets were processed in the WinMASW® Academy software by Eliosoft by following these steps: (1) inputting the geometry of each data set (array origin coordinates, source offset, and receiver spacing). (2) For the active data sets, removing the non-relevant sections (first arrivals, coda—signals arriving after the surface waves, etc.). For the passive data sets, choosing a window length (typically 10 s). (3) For both active and passive data sets, calculating the dispersion image (Figure 5). (4) Inspecting the dispersion images from the various offsets, and picking the fundamental mode (and higher modes if they exist within the MASW data sets) at the maximum energy. In the ReMi technique (linear array), the dispersion curve is picked along the lower velocity bound (which is subjective), rather than at the maximum energy (Figure 5). (5) For each type of data set, constructed in shear wave velocity/layer thickness space, an initial model and search space is defined according to the dispersion curve. (6) Finding by the inversion process (described below), the mean

Example of data processing (results from Be’er Tuvia (BET)). In the upper section: the first column is the raw signal. The second column is the dispersion image after transformation from the time-distance domain. The third column is the marginal posterior probability density (MPPD) (Dal Moro et al., 2007), which is defined as the mean shear-wave profile. The green areas represent the model boundary. In the lower section is the “noise” analysis: A—amplitude spectrum, B—spectrograms for each component, and C—HVSR directivity.
“Noise” analysis was also processed by the WinMASW® Academy software. Amplitude spectra and spectrograms were calculated for each component (east-west, north-south, and up-down) (Figures 5A, 5B and in Supplemental material). The field experiment and the parameters of calculation of the azimuthal HVSR followed the recommendations of the site effects assessment using ambient excitations (SESAME) consortium (Bard et al., 2004) (Figure 5C and in Supplemental material).
Inversion
In the WinMASW® Academy software, the automatic inversion is conducted via genetic algorithms (Dal Moro et al., 2007). The parameterization of the model is based on the assumption of a 1D layered model consisting of a stack of homogeneous linear elastic layers over a half-space. Each layer is characterized by its thickness (except for the half-space), density, and two of the three following elastic properties: the shear-wave velocity, P-wave velocity, and Poisson’s ratio. The unknown model parameters of primary interest for our purpose are the thickness and
Results
Figures 6 to 11 display the results for each site, where each figure is divided into six sections: (a) the family of dispersion curves picked from each of the different methods, the average dispersion curve, and the frequency boundaries for inter-method comparison (between methods). (b) The inter-method variability of the theoretical dispersion curves. The variability of the results is evaluated through a coefficient of variation (CoV) which defined as:
where

The BET site results (see text for more information): (a) dispersion curves, (b) inter-method variability (dispersion curves), (c)

The BSH site results (see text for more information): (a) dispersion curves, (b) inter-method variability (dispersion curves), (c)

The CAS site results (see text for more information): (a) dispersion curves, (b) inter-method variability (dispersion curves), (c)

The HOD site results (see text for more information): (a) dispersion curves, (b) inter-method variability (dispersion curves), (c)

The MGH site results (see text for more information): (a) dispersion curves, (b) inter-method variability (dispersion curves), (c)

The RAM site results (see text for more information): (a) dispersion curves, (b) inter-method variability (dispersion curves), (c)
The numerator is the time-averaged shear wave velocity for the investigation depth range. The denominator
where
where
Table 2 summarizes the results of the six sites. The table lists the number of tests, range of depths for which all measurements are available, maximum CoV,
Summary of the results
CoV: coefficient of variation.
Figure 12 displays the noise analysis of each site: the HVSR overtime in which each curve is calculated over a time window sample of 20–30 s. HVSR directivity shows the azimuth of the signals, and the amplitude spectra illustrate the spectrum of the site.

Background “Noise” analysis. Each row represents a site. The first column is the HVSR over time, in which each curve is a sample representing a window of 20–30 s. The second column is the HVSR directivity, and the third column is the amplitude spectra.
It is worth mentioning that RVF refers to the vertical component (SV) of the shear wave and THF refers to the horizontal component (SH). Hence, anisotropy may explain some of the small variations between the two components which were integrated together in the analysis, that is, no separation between SH and SV.
Discussion
Our results indicate well-defined uncertainties in the measurements at the six different industry sites. The discussion below focuses on what is expected and measured in industrial sites in terms of the induced bias from the noise wavefield. I will also focus on the inter-method and intra-method comparisons and the depth of penetration for the various methods.
Misinterpretation in the identification of the dispersion curve modes (fundamental mode vs higher modes) can bias the results and is not covered in this study. In addition, the inversion process attempts to converge the solution to the local minimum. Thus, the introduction of known information (e.g. stratigraphic logs and water table depth) can constrain the inversion results and would certainly reduce the variability of the
The influence of the industrial seismic noise
Seismic noise in cities, and particularly in industrial areas, is driven by numerous processes such as transportation, machinery, humans, and so on. It is expected that these areas will be characterized by an unnatural wavefield, which may include sharp spikes within the amplitude spectra, amplitude variation within the spectrograms, and also a preferential direction of the HVSR measurement. Unsurprisingly, these phenomena are viewed in our observations (Figure 12).
The HVSR curves over time (Figure 12) indicate a significant variation between the samples at each site. For example, in the BSH site at ∼12 Hz, the HVSR varies for different sampling windows from 0.5 to 3.5, while from 1 to 4 Hz, it also varies significantly from 0.5 to 9. Also, in most of the sites, the HVSR directivity represents a preferential direction of the signal. For example, the HVSR directivity at the MGH site indicates seismic signals from different directions and frequencies; 2 Hz in azimuth between 30 and 60 degrees and ∼27 Hz in an azimuth mainly between 60 and 90 degrees. The amplitude spectra of all six sites are characterized by several sharp spikes that are probably due to noisy local machines.
This observation of HVSR over time, HVSR directivity, and amplitude spectra (also spectrograms within the Supplementary) indicate that some of the frequencies are biased due to industrial sources. However, a closer look at the biased frequencies at the HVSR over time and HVSR directivity reveal that most of the sites have biased frequencies lower than 4 Hz or higher than 20 Hz. The meaning of this is that they are beyond the uncertainty comparisons. Besides, despite the sharp spikes which influence most of the spectra at most of the sites, the data quality and especially the dispersion images are not compromised, even if the industrial noised is not eliminated (Figure 13).

The ReMi dispersion images of all six sites.
Inter-method variability
The agreement between the experimental dispersion curves at each site is very satisfactory, with an average CoV that varies from 2.5% to 12.6% (panel b in Figures 6 to 11 and Table 2). This confirms the robustness of the estimate of the experimental dispersion curve by the different methods, an observation consistent with previous studies (Cornou et al., 2006; Cox et al., 2014). CoVs for the dispersion relations are lower than CoVs for the
Compared to the maximum CoV values, those of the average shear-wave profiles

Inter-method variability—the normalized
It is expected that the CoV values would increase with depth due to lower resolution in deeper layers. However, since in all our sites, the
Intra-method variability
The calculated

Intra-method variability—the
The low CoV values within the intra-method variability reflect only the inversion processing stage. This is because I compare data from the same deployments for each site. I also use the same sources for active methods, and I did not process the data differently. Hence, the inversion processing stage is the only stage presenting variation. The surface wave inverse problem is mathematically ill-posed, and it is affected by solution non-uniqueness (Luke et al., 2003). Therefore, several different shear-wave profiles can produce the same theoretical dispersion curve.
Depth of investigation
Besides the instrumental limits (sensor and source), two main factors will define the depth of investigation: the

Maximum depth of investigation following the different methods.
Certainly, to provide a better constraint to the model, both in the shallow and the deep subsurface, a combination of the information from both active and passive data is essential (Garofalo et al., 2016b). However, this discussion is beyond the scope of this study.
Conclusion
In this study, I explored six different sites located in noisy environments. I acquired data by five different types of surface-wave methods and compared the results to estimate the inter-method and intra-method variabilities at each site. The purpose of the research was to apply common practice site characterization approaches while also evaluating approaches not commonly applied (RVF, THF). Our results indicate that, while the average CoV of the different dispersion curves is up to ∼13% at specific frequencies, the frequency-averaged CoV is about 8%. The same concept holds for the CoV of the
On the contrary, the average index
In our tests, the data available for deeper depths came from the passive tests, particularly from the ESAC method. These confirm the common wisdom that passive tests provide lower frequencies than active methods and therefore are more suitable for deeper seismic characterization.
It is possible that more robust survey design using best practices may result in much lower uncertainty. These may include longer arrays, larger energy sources for active source data acquisition, and use of passive arrays with better azimuthal averaging (T shape, fan shape, triangular, circular, etc.).
Although our study focuses on surface-wave measurements in noisy environments, it aligns with comparable shear-wave studies and falls within the acceptable uncertainties in the results. This suggests that such tests are also applicable in industrial zones where the signals processed may be affected by noise. This conclusion has practical implications for use in estimating seismic hazards in industrial zones.
Supplemental Material
sj-docx-1-eqs-10.1177_8755293020988029 – Supplemental material for Shear-wave velocity measurements and their uncertainties at six industrial sites
Supplemental material, sj-docx-1-eqs-10.1177_8755293020988029 for Shear-wave velocity measurements and their uncertainties at six industrial sites by Yaniv Darvasi in Earthquake Spectra
Footnotes
Acknowledgements
I thank HUJI’s Neev Center for Geoinfomatics for its equipment and computational facilities and the great contribution of labor and assistance from its students. I am especially grateful to Dr John K. Hall and Prof. Amotz Agnon, co-founders of the Center. I thank Amos Shiran and Hanan Alexander for advice and support. Finally, I especially thank Eli Ram for his field assistance.
Declaration of conflicting interests
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.
Funding
The author(s) received no financial support for the research, authorship, and/or publication of this article.
Supplemental material
Supplemental material for this article is available online.
References
Supplementary Material
Please find the following supplemental material available below.
For Open Access articles published under a Creative Commons License, all supplemental material carries the same license as the article it is associated with.
For non-Open Access articles published, all supplemental material carries a non-exclusive license, and permission requests for re-use of supplemental material or any part of supplemental material shall be sent directly to the copyright owner as specified in the copyright notice associated with the article.
