
==== Front
Heliyon
Heliyon
Heliyon
2405-8440
Elsevier

S2405-8440(24)12208-6
10.1016/j.heliyon.2024.e36177
e36177
Research Article
Inversion method of NHV based on novel model parameter space acquisition strategy and its application in the Tonghai basin site in Yunnan, China
Wang Jixin a
Rong Mianshui rongmianshui@bjut.edu.cn
ab⁎⁎
Li Xiaojun lixiaojun@bjut.edu.cn
ab⁎
a Key Laboratory of Urban Security and Disaster Engineering of China Ministry of Education, Beijing University of Technology, Beijing, 100124, China
b Institute of Disaster Prevention, Sanhe, Hebei, 065201, China
⁎ Corresponding author. Key Laboratory of Urban Security and Disaster Engineering of China Ministry of Education, Beijing University of Technology, Beijing 100124, China. lixiaojun@bjut.edu.cn
⁎⁎ Corresponding author. Key Laboratory of Urban Security and Disaster Engineering of China Ministry of Education, Beijing University of Technology, Beijing, 100124, China. rongmianshui@bjut.edu.cn
19 8 2024
15 9 2024
19 8 2024
10 17 e361777 11 2023
11 8 2024
12 8 2024
© 2024 The Authors
2024
https://creativecommons.org/licenses/by-nc-nd/4.0/ This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).
The imaging of subsurface soil velocity structures from ambient noise inversion is a difficult problem. Few recording points and a simplified 1-D layered profile lead to important non-uniqueness. From our point of view, improving the reliability of processing methods of the observed data to obtain noise horizontal-to-vertical spectral ratio (NHV) curves and setting a complete model parameter space are important tasks to reduce the non-uniqueness of inversion. In this study, using a local site near the border of the Tonghai Basin, China, as a case study, we first demonstrate how to identify and mitigate the influence of industrial sources using surface observations to obtain more reliable NHV curves. Then, a new strategy to determine model parameter space is proposed, that is, stratifying soil layers based on the number of NHV peaks and determining the shear wave velocities, thicknesses, and their ranges based on the empirical relationship between sedimentary thickness and resonant frequency (h-fr). Subsequently, combining the model parameter space acquisition strategy with the NHV inversion, a novel NHV inversion approach is developed and applied to obtain the 2-D VS profile of the investigated Tonghai site. The inverted 2-D VS profile aligns favorably with the frequency-depth conversion results of the measured NHV curves (NHV-profiling) and the measured borehole profiles, affirming the reliability of the proposed NHV inversion method. Finally, by comparing the empirical transfer functions from the strong-motion recordings, we validated the applicability of the inverted models for characterizing site effects. The model parameter space acquisition strategy proposed in this paper and the analysis procedure of the observed data are also applicable to other study areas, which can provide a referable approach to quickly and effectively acquire the soil layer velocity structure of the site.

Keywords

Noise horizontal-to-vertical spectral ratio (NHV)
Initial model
Resonance frequency
Shear wave velocity
Empirical relationship
==== Body
pmc1 Introduction

Numerous seismic damage observations and theoretical research results have demonstrated that subsurface soil velocity structures in a region significantly influence seismic wave propagation, especially reflections and refractions of seismic waves in subsurface layers markedly amplifying ground motions, which is one of the essential causes of severe damage to many buildings and facilities [[1], [2], [3], [4]]. Accurately determining subsurface soil velocity structures is crucial for reasonably evaluating site seismic effects and effectively implementing engineering seismic prevention.

Methods for determining subsurface sedimentary soil velocity structures include drilling, shallow seismic prospecting, electromagnetics, and ambient noise inversion. Among these, the ambient noise inversion method has the advantages of being simple, fast, economical, and less site-constrained than other methods. In the past two decades, this method has developed rapidly. For example, since 2001, European scientists have carried out the SESAME (Site Effects Assessment using Ambient Excitations) project [5] to develop noise-based methods for imaging subsurface velocity structures and estimating site effects. Arai and Tokimatsu [6] proposed a formulation for the Green's function based on the surface wave assumption, representing it as a summation of Rayleigh and Love wave modes. They used the observed NHV at six sites to invert and estimate shallow soil layer profiles down to 100 m depth. Based on the body waves assumption, Herak [7] developed a ModelHVSR tool to invert shallow velocity structures using ambient noise. Bignardi et al. [8] developed the OpenHVSR program for simulating and inverting large datasets of NHV, aimed at building 2-D/3-D subsurface profiles, assuming either body or surface waves. García-Jerez et al. [9] proposed the HV-inv program based on the diffuse field assumption (DFA), which considers the NHV as the result of the combination of body and surface wavefields. Guo et al. [10] analyzed the polarization characteristics of Rayleigh waves using one year of continuous ambient noise records from 30 movable array stations along the east coast of the United States, based on the fundamental mode Rayleigh waves, which improved the reliability of NHV methods in site effect analysis. Rong et al. [11] combined the inversion algorithm of genetic and simulated annealing with the forward NHV calculation based on diffuse field assumption and proposed a global inversion method for soil layer velocities. From the above discussions, it can be seen that NHV inversion has continuously developed with the understanding of noise wavefields, from initially considering only the contribution of body or surface waves alone to being able to take into account the combined effect of different wave types. However, existing NHV inversion methods still have significant non-uniqueness problems.

Determining an appropriate initial model parameter space is crucial for reducing non-uniqueness in inversion results. The process of inverting soil velocity structures can be simply described as follows: 1) calculating theoretical NHV curves based on the initial model (predicted model parameters); 2) comparing the theoretical NHV curves with observed ones and calculating their errors; 3) randomly varying model parameters within a preset model parameter space to generate predicted models - when the errors of NHVs between predicted model and observation converges to an acceptable preset value, the predicted model is taken as the inverted result; otherwise, iterative repetition is performed until convergence. Two main categories of inversion methods are primarily used in the inversion process: local linearized and nonlinear global optimization methods. Relevant analysis indicates that local linearized methods are easily trapped in local minima, and the reliability of their results largely depends on the selection of the initial model. In theory, the convergence of nonlinear global optimization methods is unaffected by the initial model, but practical results show that excessively crude initial models significantly reduce the convergence speed. Also, if computation time is limited, the obtained results may still only be local optimal solutions [11,12].

Additionally, according to research by Scherbaum et al. [13] and Piña-Flores et al. [14], using the NHV alone cannot accurately explain the trade-off relationship between layer thickness and shear wave velocity, resulting in non-unique inversion results. Piña-Flores et al. [14] pointed out that using a priori information (such as borehole, geological, and geophysical data) or analyzing the shape of the observed NHV itself facilitates determining the upper and lower bounds of the target values, thereby improving the solving efficiency. Research by Perton et al. [15] showed that inverting NHV curves at multiple locations and constraining thickness and velocity variations in the profile direction helps to obtain a coherent and consistent final subsurface structure. To improve the uniqueness and convergence of the solution, Bora et al. [16] utilized existing standard penetration test (SPT) data and earlier soil measurement data to determine the model parameter space (including the number of layers, thicknesses, and shear wave velocity ranges). They finally combined the shear wave velocity profiles obtained from NHV inversion with existing soil information to acquire a better subsurface 3D shear wave velocity model. Thomas et al. [17] manually determined the combinations of layer thicknesses and velocities required to generate results that fit the observed NHV curves well at a particular site, considering independent constraints from in-situ thickness measurements. They used the resulting velocity structure as an initial model for inversion at other sites. In the inversion strategy, no constraints were imposed on layer thicknesses, but each layer's shear wave velocity (VS) was allowed to vary only within 15 % of the initial model. When studying shallow VS profiles in Taiwan, Chen et al. [18] alleviate the non-uniqueness of inversion by utilizing ambient noise arrays and velocity structural data from engineering geological repositories. In order to obtain VS of the Bengal Basin through NHV inversion, Farazi et al. [19] utilized existing boreholes and VP obtained from active seismic surveys to constrain the initial model parameter space. Overall, the initial model and its model parameter space significantly impact the inversion results. Existing research relies heavily on prior geological survey data when setting up initial models, and there is a problem in setting the range of model parameter variations. Therefore, reducing the reliance on prior knowledge and rapidly determining the initial model and its parameter variation range is the key step in NHV inversion.

The accuracy of noise data processing has an essential impact on the non-uniqueness of inversion. If the target NHV curves for inversion intrinsically contain significant discrepancies, it will lead to considerable deviations in the inversion results. The observed data recorded on site may encompass various types of noise; the vibrations induced by industrial activities could be misconstrued as the resonant frequencies of the site, resulting in erroneous interpretations in subsequent NHV analysis. Thomas et al. [17] found that there were significant differences in the peak values of NHV between daytime and nighttime when using a dense seismic array to study the regional site effects. They proposed that the impact of industrial sources should be taken into account and their influence mitigated before calculating NHV curves. Field et al. [20] reported some early studies describing the effects of certain human activities on NHV, while more recently, Bard and SESAME-team [5], Dal Moro [21,22], and Mohamed et al. [23] delineated the presence of industrial peaks. Before performing NHV analysis, Setiawan et al. [24] applied the random decrement technique to the noise data to detect industrial noise sources, validating the efficacy of ambient noise data utilized in subsequent NHV analysis. Dal Moro [25] conducted a systematic investigation on identifying industrial sources in observed noises, and his results show that industrial sources can propagate for tens of kilometers and contaminate the NHV. Then, he proposed two simple approaches to mitigate the effects of industrial noise. The research mentioned above has demonstrated the difficulty in avoiding the impact of industrial sources during NHV observations. Detecting the presence of industrial sources and mitigating their impact is necessary before calculating NHV curves.

In this study, we first demonstrate how to identify and mitigate the influence of industrial sources from surface observations and obtain more reliable NHV curves. Then, we suggest a novel model parameter space acquisition strategy based on surface observation data. This strategy stratified the layers according to the number of peaks in the observed NHV curves and determined the model's initial shear wave velocity and thickness ranges based on the observed peak frequencies combined with the universal h-fr empirical relationship, thereby significantly reducing dependence on a priori information. Based on the above two aspects, an NHV inversion method for shallow subsurface velocity structure based on a novel model parameter space acquisition strategy is proposed, which can effectively mitigate inversion non-uniqueness. This method successfully obtains a 2D velocity profile for the strong motion observation array (SMOA) site at the margin of the Tonghai Basin, Yunnan, China, through ambient noise array observation. Finally, the applicability and reliability of the proposed method are validated by comparing the inverted velocity profile with borehole profiles and results from other methods.

2 Data acquisition and processing

The study focuses on the Tonghai SMOA site, depicted in Fig. 1, situated in Sijie Town, Tonghai County, Yunnan Province, in southwest China (24.13°N, 102.73°E). The primary objective of this investigation is to acquire a 2-D subsurface shear-wave velocity (VS) profile of this area using the NHV inversion method. To this end, we established a continuous ambient noise observatory array over 24 h from November 29 to November 30, 2019. The array comprises 12 stations, with inter-site distances of 28m, 59m, 37m, 52m, 32m, 80m, 46m, 81m, 105m, 55m, and 23m respectively. Each station recorded data at a sampling rate of 500 Hz. The observations were conducted using three-component short-period velocity-type seismometers (frequency range: 5s - 240Hz, model: IES-3C). Given the site-specific conditions, we positioned the sensors in a linear array, as depicted in Fig. 1. In order to minimize the impact of wind-induced noise and human interference, the sensors were securely placed in trenches and independently deployed at each observation point in the array to collect data. Time synchronization was automatically achieved through the reception of GPS satellite signals.Fig. 1 Distribution of the temporary ambient noise array observation points (black circles) along the Tonghai strong motion observation stations (red solid circles) in Yuxi City. The black triangles are the drilling boreholes, and the red rectangle of the inset (lower left) shows the location of the research area in China. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

Fig. 1

2.1 Sources of recorded ambient noise data

The noise signals generated by natural phenomena such as air pressure, wind, and tide have a dominant frequency below 1 Hz, while high-frequency noise sources (>1 Hz) are primarily from human activities, followed by multiple scattering of ocean noise sources in the subsurface media [26]. In regional studies, researchers generally use low-frequency ambient noise below 1 Hz, focusing primarily on the intensity, components, frequency range, and directionality of the noise rather than the cause of the noise. For smaller-scale subsurface studies, since the high-frequency ambient noise used is primarily determined by the temporal and spatial characteristics of human activities near the array, both surface wave and body wave components in the ambient noise are important research targets. Therefore, analyzing the causes of noise sources is crucial for utilizing different noise components [27].

When collecting noise data, it is imperative to ensure that no specific vibration sources (such as vehicular transportation, construction activities, etc.) affect the observation site within a certain area. Optimal conditions for collecting high-quality noise signals include selecting calm weather conditions with weak winds for on-site data collection. Considering these factors, for the 24-h data recorded by the noise array, this study extracted the three-component recordings from 2 a.m. to 3 a.m. on November 29, 2019, during clear weather conditions with calm winds and exceptionally quiet surroundings. The data is resampled to 100Hz and subjected to demeaning and detrending processing. Fig. 2 shows the preprocessed noise waveform records of 12 stations and the FFT amplitude spectra. The amplitude spectra in Fig. 2(b) indicate that the main energy is concentrated above 1Hz, indicating that the source of the recorded environmental noise data is mainly related to human activities. As shown in Fig. 3, we calculated the power spectral density of the horizontal components at 12 stations (Fig. 3(a)), provided the corresponding mean power spectral density curves and their standard deviations (Fig. 3(b)), and presented the distribution of the horizontal component power spectral density along the measurement line (Fig. 3(c)). Additionally, we calculated the power spectral density of the vertical components at the 12 stations (Fig. 3(d)), provided the corresponding mean power spectral density curves and their standard deviations (Fig. 3(e)), and presented the distribution of the vertical component power spectral density along the measurement line (Fig. 3(f)). The power spectral density curves for the horizontal and vertical components in the 0.1–0.4 Hz frequency band are similar. Haubrich et al. [28] believed that the microseismic signal in this frequency band results from coupling between seismic waves and solid earth vibration, and there is a main peak between 0.1 and 0.4 Hz. Tanimoto [29] believed that the microseismic component in this frequency band is composed of Rayleigh waves excited by nearshore water pressure perturbations. In the frequency band above 0.5 Hz, the power spectral density curves exhibit significant changes, indicating the influence of local site effects at the experimental site. From the power spectral density distribution in the frequency spectrum, it can be observed that the ambient noise power spectral density at sites between 200 and 600 m is greater than that at sites between 0 and 100 m for frequencies above 0.5 Hz. Based on the correlation between site effects and layer thickness, this suggests that the sediment layer thickness within the 200–600 m range may be greater than that within the 0–100 m range.Fig. 2 Noise observation recording of the 12 stations. (a) Waveforms of noise recordings; (b) Amplitude spectra of noise recordings.

Fig. 2

Fig. 3 Power spectral density analysis of ambient noise data at the experimental site. (a) Power spectral density curves of the horizontal components at 12 stations; (b) Average power spectral density curve of the horizontal components (solid line) and corresponding standard deviation (dashed lines); (c) Map of the power spectral density of horizontal components along the survey line; (d) Power spectral density curves of the vertical components at 12 stations; (e) Average power spectral density curve of the vertical components (solid line) and corresponding standard deviation (dashed lines); (f) Map of the power spectral density of vertical components along the survey line.

Fig. 3

Before calculating the NHV curve, it is necessary to mitigate the influence of industrial sources, as vibrations from structures such as machinery, buildings, and trees can contaminate the signal and lead to erroneous NHV interpretations. The guidelines for implementing the H/V spectral ratio technique on ambient vibration measurements, processing, and interpretation [5] offer various ways to verify the origin of peaks, especially those caused by human activity. The random decrement technique (RDT) is a special averaging procedure used to determine step and impulse responses based on random responses under stable random vibration conditions. An effective method for checking is to apply the random decrement technique to the ambient vibration records to obtain "impulse response" near the target frequency: if the corresponding damping is very low (e.g., below 1 %), it can be almost certain to be of artificial origin, and frequency should not be considered in the interpretation [24,30,31]. Based on this, we performed the random decrement technique analysis on the extracted 1-h observed ambient noise data. As depicted in Fig. 4(b–c), the random decrement technique analysis results of the vertical observation records from stations TH06 and TH07 demonstrate that their damping ratios are far below 1 % while the frequencies remain unchanged, indicating the existence of industrial sources in the vertical observation records of stations TH06 and TH07. More details regarding the random decrement technique analysis of all the observation records can be found in Appendix A of the supporting information. In order to further evaluate the validity of the recorded data, we performed the time-frequency analysis of the extracted 1-h observed noise data. The result shows (Fig. 4(d-f)) that the vertical component recorded by the TH05, TH06, and TH07 stations continuously has noise with a vibration frequency of about 9Hz. Detailed time-frequency analysis results of all stations can be found in Appendix B of the supporting information. In addition, have a look at the raw spectra from each individual window: if they all exhibit a sharp peak (often on the three components together) at this particular frequency, there is a 95 % chance that this is anthropic "forced" ambient vibration and it should not be considered in the interpretation [5]. Fig. 4(g–i) illustrates that this narrow peak appears at the same 9Hz frequency in the three-component spectra at stations TH05, TH06, and TH07. Through a comprehensive analysis that included a random decrement technique, time-frequency analysis, and spectrum analysis, we concluded that the industrial source exists in the noise of TH05, TH06, and TH07 stations, with a vibration frequency of approximately 9 Hz.Fig. 4 Example results of identification of data from an industrial source. (a–c) The damping effect on the noise records of the TH05, TH06, and TH07 Z components. f represents frequency, z represents the damping factor, and N win represents the number of windows analyzed. The red line represents the fitted damping, while the black line represents the average damping. (d–f) Results of time-frequency analysis for Z-component of TH05, TH06, and TH07. (g–i) Power spectral density curves for the three components at TH05, TH06, and TH07 stations. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

Fig. 4

2.2 Calculation of NHV and reliability analysis

The open-source GEOPSY program [32] was utilized to process noise data in this study. The amplitude spectra of three components within selected time windows were computed, and the wavefield statistics in the frequency domain were analyzed. The width of these time windows is contingent on the target frequency band and record length, as the window length determines the minimum resolvable frequency [5]. Following the SESAME guidelines recommending at least ten cycles per window, the minimum frequency (in Hz) was defined as 10/lw, where lw is the window length in seconds. In this study, to ensure the reliability of the minimum frequency of 0.1 Hz, we set lw to 100 s. For each time window, a taper with a width of 5 % is applied to avoid potential truncation effects. Albarello and Lunedei [33] proposed several alternative methods for averaging horizontal components to calculate NHV spectral ratios and concluded that vector summing is the most physically logical averaging method. In light of this, this study employs vector summing to combine the amplitude spectra of the horizontal components (H=HN2+HE2). An anti-trigger algorithm [34] is applied to the noise recording sequence to include only signals with quasi-stationary amplitudes in the calculation. The selection criteria are based on comparing the short-term average (STA) and long-term average (LTA) amplitudes, with the typical values of STA, LTA, Min STA/LTA, and Max STA/LTA set to 1, 30, 0.2, and 2.5 s, respectively. The upper bound of the STA/LTA ratio aims to exclude spikes produced, for example, by operators walking close to the sensors. We overlap the time windows by 50 % for all the considered cases. In the frequency analysis, we apply a smoothing function with a b-value of 35 (using the Konno-Ohmachi approach) to ensure consistent data points for both low and high frequencies [35,36]. For the TH05, TH06, and TH07 stations affected by industrial sources, we apply band-stop filtering to mitigate specific industrial source frequencies before calculating the NHV curves. Taking the TH09 observation point as an example, Fig. 5 demonstrates the data processing steps for obtaining the NHV curve. First, the original noise data (Fig. 5(a)) is subjected to random decrement technique analysis (Fig. 5(b)) and time-frequency analysis (Fig. 5(c)). Subsequently, the NHV curve (Fig. 5(d)) and its corresponding directional NHV spectrum (Fig. 5(e)) are calculated. Due to space limitations, NHV curves obtained for the remaining stations can be found in Appendix C of the supporting information.Fig. 5 Different steps of data processing at the TH09 observation point. (a) Raw data; (b) Influence of the damping on ambient noise recordings of the three components (EW, NS, and UD). f is frequency, z is the damping factor, and N win is the number of windows analyzed. The red line is the fit damping, and the black line is the average damping; (c) Represents the time-frequency spectrum of the three components; (d) the NHV curve; and (e) The corresponding directional NHV spectrum. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

Fig. 5

Since noise observation results are susceptible to various factors such as the environment, site conditions, and instrumentation, verifying the reliability of the NHV curves obtained from the observations is necessary. The guidelines for implementing the ambient vibration H/V spectral ratio technique [5] provide criteria for assessing the reliability and peak clarity of the NHV curves (refer to SESAME H/V User Guidelines "2.1 Criteria for reliability of results"). Based on these criteria, this study examined the observed NHV curves of the Tonghai SMOA site (as shown in Table 1), and the results confirm that the NHV curves obtained in this study meet the reliability and peak clarity criteria.Table 1 Fundamental resonance frequency, maximum horizontal-to-vertical-spectral ratio, SESAME reliability criteria checks, and dominant direction.

Table 1No.	Latitude (°)	Longitude (°)	f0	A0	Check for reliability	Check for clarity	Dominant direction	
i	ii	iii	i	ii	iii	iv	v	vi		
1	24.16062	102.69442	3.847	7.05	PASS	PASS	PASS	PASS	PASS	PASS	PASS	FALL	PASS	None	
2	24.16037	102.69446	3.445	4.55	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	0°–110°	
3	24.15985	102.69435	1.488	3.89	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	None	
4	24.15952	102.69436	1.240	3.71	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	None	
5	24.15906	102.69428	1.034	3.89	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	None	
6	24.15877	102.69428	1.051	4.02	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	150°–260°	
7	24.15805	102.69424	0.926	5.59	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	None	
8	24.15763	102.69424	1.057	5.13	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	150°–200°	
9	24.15692	102.69402	0.750	4.18	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	None	
10	24.15596	102.69411	0.706	4.86	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	160°–210°	
11	24.15547	102.69408	0.695	5.25	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	120°–190°	
12	24.15526	102.69410	0.646	3.10	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	PASS	50°–130°	

In order to assess site conditions, including subsurface layer anisotropy and topographic effects, this study employed NHV analysis to investigate the directional characteristics of seismic effects at the experimental site. The noise recordings obtained from the experimental array contain the effects of the actual three-dimensional site conditions. The NHV analysis, besides accurately estimating the fundamental resonance frequency of the site, also allows for the study of the dominant directions of seismic effects on the site. Del Gaudio et al. [37] were the first to study the directional characteristics of environmental noise using the NHV method. By comparing the results with actual strong motion observation data, they concluded that NHV analysis based on site environmental noise could effectively identify the dominant directions of seismic effects on the site. Since then, the research results of numerous scholars have also confirmed this conclusion [[38], [39], [40], [41], [42], [43], [44], [45]]. Given this, we utilized directional NHV spectra to investigate the dominant direction of seismic effects at the experimental site, following the judgment criteria proposed by Del Gaudio et al. [37,38]. Additionally, we conducted a preliminary examination of 2-D effects at the site by analyzing the correlation between resonance frequency and azimuth angle, as illustrated in Fig. 5(e). At the TH09 measurement point, the directional features were not prominent due to amplification within the 0°–180° range, which made it challenging to determine the dominant direction of the seismic effects. The determination results of the dominant directions for the 12 measurement points in this study are listed in Table 1. On the other hand, the directional NHV spectra indicate almost no correlation between the resonance frequency and azimuth angle at the TH09 measurement point, suggesting the absence of 2-D effects. Due to space limitations, further details regarding such analyses are included in Appendix D of the supporting information.

In this section, we explored the causes of noise sources and mitigated the influence of industrial sources on NHV. We investigated the dominant directions of seismic effects at the observation sites and examined the impact of 2-D effects at these sites. These investigations ensure the reliability of the subsequent NHV inversion of the site velocity structure and provide a basis for evaluating the uncertainty of the inversion results.

3 NHV inversion method based on diffuse field assumption

The Diffuse Field Assumption (DFA) model posits that the noise wavefield behaves as a diffuse field, adhering to the principle of energy equipartition, thereby comprehensively considering the combined effects of body waves and surface waves. When the seismic source and receiver are located at the same position, a proportional relationship exists between the power spectrum and the imaginary part of the Green's function [46]. The DFA model is consistent with the assumption of noise-based interferometric imaging [47], thus gaining wider acceptance [9,14,15,18,48]. This study employs the DFA model to interpret NHV curves, with the specific formula as follows [46]:(1) NHV(ω)≡P1(ω)+P2(ω)P3(ω)=|F1(ω)|2+|F2(ω)|2|F3(ω)|=Im[G111D(0,0;ω)]+Im[G221D(0,0;ω)]Im[G331D(0,0;ω)].

In this equation, Pi(ω) represents the power spectral density of the ground surface records, Fi(ω) denotes the corresponding Fourier spectrum, and Im[Gii1D(0,0;ω)] signifies the imaginary part of the displacement Green's function of the source and receiver of the surface observation point at frequency ω and component i. In the computation, G11 denotes the north-south component of the ambient noise records, G22 represents the east-west component, G33 indicates the vertical component, and Pi(ω)∝Im[Gii1D(0,0,ω)] for i = 1, 2, 3.

In the above equation, the imaginary part of Green's function is related to each soil layer's VS、VP、h, and ρ. Due to the complexity of its function form, which involves solving a matrix composed of VS、VP、h, and ρ of each layer, as well as a wave number integral of the same soil layer, a detailed derivation process and calculation expression for the imaginary part of the layered site Green's function are provided in the appendix of Sánchez-Sesma et al. [46].

The theoretical NHV can be computed using equation (1) based on the DFA theory. In this study, we utilized the HV-inv program (https://w3.ual.es/GruposInv/hv-inv) developed by García-Jerez et al. [9] to perform forward computation and inversion analysis of NHV. The inversion is achieved by finding the minimum value of the objective function. In this study, the objective function ΦNHV is defined as:(2) ΦNHV=∑i=1n(NHVobs(ωi)−NHVth(ωi))2σNHV2(ωi),

where NHVobs is the observed NHV at frequency ω, NHVth is the NHV calculated based on the DFA theory at frequency ω, and σNHV is the standard deviation of NHVobs. The model with the minimum objective function value is adopted as the inverted best model.

4 New initial model and model space acquisition strategy

This study introduces a novel strategy to obtain initial models for inversion. At first, we employ 1-D shear wave propagation theory to establish empirical relationships between resonance frequency (fr) and sediment thickness (h). Meanwhile, the empirical relationships from previous studies for different regions have been analyzed to offer references for deriving empirical relationships suitable for our study area. Subsequently, we examine the characteristics of the measured NHV curve and ascertain the count of interfaces (N interfaces excluding the free surface) by identifying the distinct peaks present in the NHV curves (N peaks) [14,49,50]. Finally, we utilize the peak frequencies corresponding to each clear peak, in conjunction with the derived empirical relationships, to establish the upper and lower limits of thickness and shear wave velocity for each soil layer. This process allows us to construct the initial model essential for NHV inversion.

4.1 Relationship between the resonance frequency and sediment thickness

The thickness (h) of the sediments can be estimated from the empirical relationship between the site's fundamental resonance frequency (fr) and the average shear wave velocity (V‾s) of the sediments: h = V‾s/4 fr [51,52]. However, such estimates only provide a rough approximation of the sedimentary thickness. Ibs-von Seht and Wohlenberg [53] considered the increase in VS with depth to better estimate sedimentary thickness. It proposes an empirical relationship between sedimentary thickness (h) and resonance frequency (fr):(3) h=[VS0(1−x)4fr+1]1/(1−x)−1

When h≫1 and VS0(1−x)≫4fr, equation (3) can be expressed as:(4) h=AfrB

i.e.,A=[VS0(1−x)4]1/(1−x);B=−11−x

In the above equation, coefficient A represents the soil layer thickness at a frequency of 1 Hz, while coefficient B governs the slope of the regression line. When B equals −1, the average shear wave velocity (V‾S) remains constant with depth. Regression parameter A is linked to surface velocity and exhibits significant regional variations, rendering it a characteristic parameter for each area. Conversely, parameter B pertains to the exponential component of the S-wave velocity curve change with depth, showing consistency across studies and minimal regional variation. This implies that the absolute velocities of sedimentary layers exhibit substantial differences among different areas, while the curve of shear wave velocity with depth in different areas has a similar change pattern [49,54,55].

We collected and analyzed relevant research that applied the bivariate regression model to different regions [49,50,[53], [54], [55], [56], [57], [58], [59], [60], [61], [62], [63], [64], [65], [66], [67], [68], [69], [70], [71], [72], [73], [74], [75]]. Based on the regression parameters (A and B values) and their frequency band range provided by previous studies, the corresponding values of surface shear wave velocity VS0 and depth-dependent parameter x were determined, as shown in Appendix E of the supporting information. Looking back at the above derivation process, the validity of equation (4) relies on the conditions h≫1 and VS0(1−x)≫4fr. The condition of h≫1 is easy to satisfy, and the assumption of VS0(1−x)≫4fr is more stringent. In this study, we statistically analyzed the corresponding VS0 and x values in equations (3), (4), as shown in Fig. 6. Based on equations (3), (4), the root-mean-square error of the surface shear wave velocity (VS0) obtained from both is 20 m/s, and the root-mean-square error of the depth-dependency parameter (x) is 0.029. The average value of VS0 obtained based on equation (3) is 214 m/s with a standard deviation of 108 m/s, and the average value of x is 0.258 with a standard deviation of 0.134. The average value of VS0 obtained based on equation (4) is 231 m/s with a standard deviation of 109 m/s, and the average value of x is 0.236 with a standard deviation of 0.127. Comparing the empirical relationships given by different researchers indicates that a regional regression model is a function of local surface conditions and regional geological features. Surface conditions vary little at small spatial scales, while different soil/bedrock types are reflected in different resonance frequencies and bedrock depths [50]. Since VS0 is a characteristic parameter in each region (obtainable through simple means), it exhibits significant variations across different regions, while x shows minimal variability across different regions. For this purpose, this study collected VS0 from eight boreholes in the Tonghai Basin, with an average VS0 of 160 m/s and a standard deviation of 12 m/s, as shown in Fig. 1. By utilizing h = V‾s/4 fr and equation (4), the interrelationship between the average shear wave velocity of the overlying soil layer (V‾s), the resonant frequency (fr), and the thickness (h) of the overlying soil layer can be expressed as:(5) V‾s=4Afr(B+1),

(6) V‾s=4A(−1/B)h(B+1)/B.

Fig. 6 VS0 and x derived from published resonance frequency-sediment thickness relationships.

Fig. 6

Given that the condition for the VS0(1−x)≫4fr assumption to hold is quite stringent, this study provides the statistical parameter x based on formula (3) using the grid search method. According to the above statistical analysis, the surface shear wave velocity VS0, which represents the Tonghai basin, is 160 ± 12 (m/s), and the depth-dependent parameter x is 0.258 ± 0.134, while the regression parameters are A = 96.48 and B = −1.348.

Thus, the S-wave velocity-depth function for the Tonghai Basin is as follows:(7) VS(z)=(160±12)(1+Z)0.258±0.134

The empirical relationship between the resonance frequency and the thickness of the covering soil layer in the Tonghai Basin is as follows:(8) h=96.48fr−1.348

The relationship between the resonant frequency and the average shear wave velocity of the covering soil layer in the Tonghai Basin is as follows:(9) V‾s=385.92fr−0.348

The relationship between the thickness of the soil layer and the average shear wave velocity of the covering soil layer in the Tonghai Basin is as follows:(10) V‾s=118.63h0.258

Researchers have typically derived the relationship between resonant frequency and sediment thickness for regional basins within the frequency range of 0.1–10 Hz, focusing on alluvial sites characterized by significant impedance contrast between soft sediments and bedrock [76]. Fig. 7(a) compares the empirical h-fr relationships given by previous researchers. Similarities and parallelism exist among most empirical relationship curves in the log-log coordinate system, illustrating parameter B's small variability among different regions. Fig. 7(b) shows the h-fr empirical relationship obtained in this study, which is very similar to the one given by Ibs-von Seht and Wohlenberg [53], and the reason for this is that the geological conditions of both have certain similarities. To further assess the validity of the obtained h-fr relationship equation, this study utilized borehole data from eight locations uniformly distributed across the Tonghai Basin. The peak frequencies of the H/V spectral ratio curves were calculated for the experimental sites using the approximate formula proposed by Tuan et al. [77] in a multilayer model, as shown in Table 2. The peak frequencies and corresponding borehole depths obtained from calculations based on the borehole model, along with the sediment layer thicknesses obtained from inversion and the observed peak frequencies, are plotted in Fig. 7(b), as shown by the blue and green triangles. It can be seen that the h-fr relationship proposed in this study demonstrates good reliability.Fig. 7 (a) Comparison of proposed relationship (Eq. (8)) between estimated bedrock depth (sediment thicknesses) and resonance frequency with published relationships following appendix E of the supporting information. The resonance frequency range between 0.1 and 10 Hz is plotted as published relationships are developed mainly up to 10 Hz. (b) Comparison between the present study regression relationship for the Tonghai basin and Ibs-von Seht and Wohlenberg [53] derived a regression relationship for the Lower Rhine embayment. The triangular markers in blue represent the results from drill data, while those in green depict the results derived from inversion data. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

Fig. 7

Table 2 Information on borehole stations in the Southwest China Tonghai basin was collected from a previous study [78].

Table 2Sta	Z (m)	VS0 (m/s)	VSZ (m/s)	fr(Hz)	
B01	65	181	745	1.69	
B02	49	159	321	1.13	
B03	80	159	321	0.69	
B04	93	159	760	0.63	
B05	80	145	398	0.90	
B06	82	152	370	0.83	
B07	83	153	357	0.64	
B08	30	173	456	2.39	
fr based on Tuan et al. [77].

4.2 Sensitivity analysis of NHV curves to stratified interfaces

In the context of 1-D modeling, we investigated the sensitivity of NHV curves to stratified interfaces within the frequency range of engineering interest (0.1–20Hz) using DFA and 1-D numerical simulations. Taking a three-layer model (with the last layer representing the bedrock half-space) as an example, NHV curves exhibit two prominent peaks. As illustrated in Fig. 8(a), varying the thickness of the first layer affects only the resonance frequency and amplitude of the high-frequency peak of the NHV curve; as depicted in Fig. 8(b), changing the thickness of the second layer affects only the resonance frequency and amplitude of the low-frequency peak of the NHV curve; as shown in Fig. 8(c), altering the shear wave velocity (VS) of the first layer affects only the resonance frequency and amplitude of the high-frequency peak of the NHV curve; as demonstrated in Fig. 8(d), modifying the VS of the third layer affects only the resonance frequency and amplitude of the low-frequency peak of the NHV curve. In summary, we believe that the number of significant peaks in NHV is related to the number of strong wave impedance interfaces underground under the assumption of 1-D site conditions.Fig. 8 Sensitivity analysis of NHV curves to stratified interfaces. (a) Sensitivity of NHV to changes in the thickness of the first layer; (b) Sensitivity of NHV to changes in the thickness of the second layer; (c) Sensitivity of NHV to changes in the VS of the first layer; (d) Sensitivity of NHV to changes in the VS of the third layer. The color bar represents the 1-D models and their respective NHV curves. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

Fig. 8

4.3 Acquisition of initial model and model space

To present the complete framework of the proposed novel strategy for initial model acquisition in this study, Fig. 9 depicts the flowchart of the strategy. The figure illustrates the main steps involved, which are as follows.1) Processing observed noise data (identifying and mitigating industrial sources; calculating and testing the NHV curve);

2) Identify the peak value and frequency of the NHV curve;

3) Performing a statistical analysis of the h-fr relationship in different research areas and providing the depth dependence parameter x and its standard deviation SD based on the grid search method.

4) Using a simple method, it is acquiring the shear wave velocity VS0 and its standard deviation SD within the shallow surface 1-m depth of the study area.

5) Using the values of x±SD and VS0±SD obtained from 3) and 4), we determine the parameters of h-fr, hmin-fr, hmax-fr, VS(z)-z, VS(z)min-z, and VS(z)max-z.

6) Based on 2) and 5), we obtained the values of hmin-hmax and VS(z)min-VS(z)max for each soil layer.

Fig. 9 Flowchart diagram of the novel initial model acquisition strategy.

Fig. 9

Taking the TH09 observation point as an example, we can see from Fig. 5(d) that there are three clear peaks with corresponding resonant frequencies of 4.632 Hz, 2.127 Hz, and 0.75 Hz, respectively, indicating a four-layer model. We take the empirical relationship of shear wave velocity with a depth corresponding to VS0+SD and x + SD (VS(z)max=172(1+Z)0.392) as the upper limit of shear wave velocity variation with depth and the empirical relationship of layer thickness with resonance frequency (hmax=214.39fr−1.645) as the upper limit of layer thickness variation with resonance frequency. Similarly, we use the empirical relationship of shear wave velocity with a depth corresponding to VS0-SD and x-SD (VS(z)min=148(1+Z)0.124) as the lower limit of shear wave velocity variation with depth and the empirical relationship of layer thickness with resonance frequency (hmin=53.03fr−1.142) as the lower limit of layer thickness variation with resonance frequency. Therefore, we obtain the upper and lower limits of each soil layer's thickness and shear wave velocity sequentially, based on the three observed resonance frequencies. Rong et al. [11] indicate that the VP, VS, and h parameters of soil layers and the VS parameter of bedrock are the main factors affecting NHV, while the influence of the density (ρ) of soil layers and bedrock on NHV can be ignored. In the initial model construction of this study, the Poisson's ratio in the soil layer was fixed at 0.33, and the bedrock layer was fixed at 0.25. The range of shear wave velocity in the bedrock layer was fixed at 400–1500 m/s. The upper and lower limits of compression wave velocity for each layer were given by the relationship between VP, VS, and ν, i.e., VP=2(1−ν)/(1−2ν)VS. The density value of each layer was given by the empirical relationship, which can be expressed as follows [79,80]:(11) ρ=1.227+1.53VS−0.837VS2+0.207VS3−0.0166VS4,

Table 3 presents the initial model for the TH09 observation site. Leveraging the obtained resonant frequencies (fr) and combining the relationships of h, hmax, and hmin with fr and the relationships of VS(z), VS(z)max, and VS(z)min with depth (z) derived above, Fig. 10 presents the parameter space for the test site, including the cover layer thickness (Fig. 10(a)), surface shear wave velocity (Fig. 10(b)), average shear wave velocity of the cover layer (Fig. 10(c)), and shear wave velocity at the bottom of the cover layer (Fig. 10(d)).Table 3 The initial model with TH09 taken as an instance.

Table 3Layer	hmin(m)	hmax(m)	VPmin(m/s)	VPmax(m/s)	VSmin(m/s)	VSmax(m/s)	ρmin(kg/m3)	ρmax(kg/m3)	Poisson's ratio ν	
1	9	17	391	1060	197	534	1497	1836	0.33	
2	13	45	433	1733	218	873	1523	2053	0.33	
3	52	282	502	3374	253	1700	1564	2287	0.33	
4	∞	∞	693	2598	400	1500	1718	2253	0.25	

Fig. 10 Parameter space analysis for the experimental site: (a) thickness of the covering layer, (b) shear wave velocities at the surface, (c) average shear wave velocities of the covering layer, and (d) shear wave velocities at the base of the covering layer.

Fig. 10

5 Inverted 2-D subsurface profiles and discussion

In engineering and seismological applications, two main features of NHV are principally considered: the peaks of NHV (corresponding frequencies and amplitudes) and the shape of the NHV curve. Numerous in situ investigations show that NHV maxima amplitudes are commonly associated with significant impedance contrasts in the soil layers [81,82]. More importantly, the shape of the NHV curve serves as the target spectrum for inverting the VS profiles of sites [8,9,32,83].

5.1 Experimental observation of 1-D/2-D resonances

Owing to the 2-D effects of the site having a different impact on the horizontal components, this study analyzed not only the directional NHV curves but also focused on the peak frequency and corresponding amplitude patterns of H/V, N/V, and E/V along the measured profile. As shown in Fig. 11, within the range of 0–175 m along the profile, the peak amplitude of E/V for a single station is always greater than that of N/V, while within the range of 175–570 m, the peak amplitude of N/V is always greater than that of E/V for a single station, which may be related to the predominant direction of the site effects. To quantify the disparities between the values of f0 and A0 provided by H/V, N/V, and E/V, we utilized the goodness-of-fit index S proposed by Anderson [84], expressed as:(12) S(p1,p2)=exp[−(p1−p2min(p1,p2))2],

where p1 and p2 are the values to be compared: S reaches its maximum value of 1 when they perfectly agree and rapidly decreases to 0 as the difference between them increases. Table 4 presents the goodness-of-fit parameters between H/V, N/V and E/V. Regarding the resonant frequencies, the values given by H/V, N/V, and E/V are generally consistent for each observation site, with goodness of fit parameters all above 0.97. For the peak amplitudes, except for the TH12 site, the goodness-of-fit parameters for the remaining sites are all above 0.93. The low fitting quality at the TH12 observation site is analyzed to be caused by the reduced amplitude of the original EW component record, possibly due to instrument malfunction. Based on the results in Fig. 11 and Table 4, treating the local sites under each observation point as 1-D sites is reasonable.Fig. 11 H/V (where H is the quadratic mean of E and N), N/V and E/V peak frequency and amplitude along the surveyed sections. The shaded areas on the surveyed profiles represent the confidence intervals.

Fig. 11

Table 4 The goodness-of-ﬁt parameter S by Anderson [84] for the fundamental resonance frequency f0 and the amplitude A0 from H/V, N/V, and E/V. Misﬁt values relative to f0 and A0 are reported in the matrix's upper and lower triangular parts.

Table 4f0/A0	H/V	N/V	E/V	H/V	N/V	E/V	H/V	N/V	E/V	H/V	N/V	E/V	H/V	N/V	E/V	H/V	N/V	E/V	
	TH01	TH02	TH03	TH04	TH05	TH06	
H/V		0.99	1		0.99	1		0.98	1		0.97	1		0.98	1		1	0.97	
N/V	1		0.98	1		0.99	1		0.98	1		0.93	1		0.99	1		0.95	
E/V	1	1		1	0.99		1	1		0.99	0.97		1	1		1	1		
	TH07	TH08	TH09	TH10	TH11	TH12	
H/V		1	0.97		1	0.98		1	0.99		1	0.96		1	0.98		0.08	0.95	
N/V	1		0.95	1		0.97	1		0.99	1		0.94	1		0.98	1		0.01	
E/V	1	1		0.98	0.97		1	1		1	1		1	1		1	1		

5.2 2-D S wave velocity model

After identifying the origins of ambient noise sources, discriminating and mitigating industrial source impacts, and evaluating the 1-D resonant behaviors at the observation site, this study conducted inversion analyses on the observed NHV curves based on DFA theory. Fig. 12 shows the NHV inversion results for 12 ambient noise observation points at the Tonghai SMOA site. The distribution of these ambient noise observation points is shown in Fig. 1. Notably, the inversion results all have relatively small objective function values, indicating that the theoretical NHV curves with minimum objective function values fit the observed NHV curves well, matching both the fundamental resonant frequencies and their peak values.Fig. 12 NHV inversion results at twelve sites (locations shown in Fig. 1). Panels (a–l) display the estimated shear wave velocity profiles at survey points TH01 through TH12. In the left panel of each figure, the red line indicates the observed NHV curve, the black line shows the calculated NHV curve, and the gray shading represents the standard deviation (SD) of the observed curve. The right panel presents the best-fit inverted velocity profile, depicted as a solid black line. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

Fig. 12

Considering the small size of the experimental site and the limited number of observation points, using the local polynomial spatial interpolation method is appropriate. This method directly employs the available data points to construct polynomial functions for interpolation, eliminating the need to determine model parameters beforehand. It can predict surface values, assess prediction accuracy using standard error measures, and generate predictive surfaces. Additionally, local polynomial interpolation offers excellent smoothness and approximation properties while being less sensitive to noisy data. In this study, the local polynomial spatial interpolation method was applied to interpolate spatially the 12 1-D shear wave velocity profiles obtained from inversion, resulting in a 2-D velocity profile, as depicted in Fig. 13(a).Fig. 13 (a) The VS map along the survey line from the noise array results. Black triangles indicate the projection of the measurement points. Borehole logs of B02, B03, and B04 are shown in the top; (b) Row NHV curves with picked maxima; (c) Depth domain NHV, where the white line is the strong wave impedance interface, produces the maximum peak amplitude of NHV.

Fig. 13

In order to verify the reliability of the inverted 2-D shear wave velocity profile, we first utilized the measured NHV curves (see Fig. 13(b)), normalized their amplitudes, and combined them with the relationship between resonance frequency and sediment thickness obtained in Section 4.1 (h=96.48fr−1.348) to convert them into depth-domain NHV curves and grid them, as shown in Fig. 13(c). The white line indicates the strong wave impedance interface that generates the maximum peak amplitude of the NHV. We can see that the strong wave impedance interface has some similarities with the inverted bedrock surface. In addition, as shown in Fig. 13(a), the inverted 2-D shear wave velocity profile also matches well with the three borehole data, further indicating that the inverted 2-D shear wave velocity profile has certain reliability.

To further validate this study's inversion results, we first analyze the characteristics of the site's fundamental resonance frequency and amplification effects through weak motion data recorded at the third strong motion station (53CD3). Next, we conducted a linear site response analysis for the inversion model of the noise observation point (TH07) near the third station. Finally, we compare observational results with theoretical calculations to verify the reliability of the inversion model.

In this study, 15 seismic records observed at the third strong motion observation point (53CD3) were selected. Based on the strong motion observation records, linear site response analysis was conducted for the site at the third observation point using the surface-to-borehole spectral ratio (SBSR) method. As shown in Fig. 14, the theoretical SBSR obtained using the inversion strategy proposed in this study is compared with the observed SBSR, showing minor discrepancies in the characteristics of fundamental resonance frequency and amplification effects. At 3Hz, the observed SBSR exhibits an additional peak compared to the theoretical SBSR, possibly due to the 35m distance difference between the strong motion observation point at 53CD3 and the noise observation point at TH07, which may lead to local velocity structure deviations between the two. Further comparison with observed strong motion records indicates the reliability of the soil model parameters obtained through the NHV inversion method proposed in this study, based on the new model parameter space acquisition strategy.Fig. 14 Comparison between the theoretical SBSR of the inverted model based on ambient noise station TH07 and the observed SBSR from station 53CD3. The black solid line represents the mean of the observed SBSR, the gray shaded area shows plus/minus one standard deviation, and the red solid line corresponds to the theoretical SBSR for the inverted model. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

Fig. 14

6 Conclusions

To fully explore the site information contained in ambient noise records, this study utilized the noise data from the Tonghai SMOA site to elucidate the origin of ambient noise sources through power spectrum analysis. Subsequently, based on the characteristics of industrial sources (long propagation distance and monochromatic correlated signals), random decrement techniques, time-frequency analysis, and spectral analysis techniques were implemented to identify noise data and effectively detect industrial noise sources. Additionally, bandstop filtering was applied to mitigate the influence of industrial noise on the NHV, thereby providing a more reliable target spectrum for subsequent NHV inversion. Finally, directional NHV analysis was used to investigate the dominant direction of seismic effects in the experimental site, and the 1-D/2-D resonance behavior of the site was understood from two aspects: the peak values (frequency and amplitude) of NHV and the overall shape of NHV curves. The NHV curves developed using this procedure lay the essential groundwork for ensuring the reliability of the implicit 1-D shear wave velocity inversion under the DFA theory.

On the other hand, we analyzed and derived an empirical relationship between resonance frequency (fr) and sediment thickness (h) using the 1-D shear wave propagation theory. We then systematically gathered and organized empirical relationships from previous researchers in diverse regions, subjecting their parameters to a detailed analysis. Following this, we established an empirical relationship between shear wave velocity and depth (VS(z)=(160±12)(1+Z)0.258±0.134), tailored specifically for the Tonghai Basin. Concurrently, we established an empirical relationship between resonant frequency and sediment thickness (h=96.48fr−1.348) in the same context. By closely scrutinizing the characteristics of the measured NHV curve and integrating the acquired empirical relationships, we proposed a straightforward method for constructing the initial model essential for NHV inversion. Building upon this foundation, we conducted NHV curve inversion using DFA theory and successfully established a 2-D VS model for the Tonghai SMOA site. Furthermore, we verified the reliability of the inversion results through comparison with the depth-domain NHV profile and three drilling profiles. Lastly, by comparing the empirical transfer functions from the strong-motion recordings, we validated the applicability of the inverted models for characterizing site effects.

The research findings suggest that the NHV inversion method developed in this study offers a cost-effective and efficient means of determining parameters for subsurface soil velocity models. Nonetheless, persistent challenges in current research, such as limited understanding regarding the influence of underground rock-soil layer boundaries, terrain, and soil media characteristics on ambient noise, along with the non-uniqueness issues associated with sedimentary soil velocity structure inversion, have hindered the method's further advancement and application. In order to address these challenges, it is necessary to undertake additional research, such as expanding the single-station NHV inversion based on diffuse field assumption to multi-point joint inversion, establishing a numerical model based on the 2-D VS model obtained from inversion, and using finite element analysis to compare the numerical results with actual seismic records to verify the accuracy of the inversion model. Furthermore, exploring the impact of site conditions on surface ambient noise should be incorporated into the research agenda.

Data availability statement

The datasets used and analyzed during the current study are available from the corresponding author on reasonable request.

CRediT authorship contribution statement

Jixin Wang: Writing – original draft, Visualization, Validation, Methodology, Investigation, Formal analysis, Data curation, Conceptualization. Mianshui Rong: Writing – review & editing, Supervision, Resources, Methodology, Investigation, Conceptualization. Xiaojun Li: Writing – review & editing, Supervision, Resources, Funding acquisition, Conceptualization.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Appendix A Supplementary data

The following are the supplementary data to this article:Multimedia component 1

Multimedia component 1

Multimedia component 2

Multimedia component 2

Multimedia component 3

Multimedia component 3

Multimedia component 4

Multimedia component 4

Multimedia component 5

Multimedia component 5

Acknowledgments

This work was supported by the 10.13039/501100001809 Natural Science Foundation of China [grant numbers 52192675 ], the 10.13039/501100012166 National Key Research and Development Program of China [grant numbers 2023YFC3000183 ], and the 10.13039/501100013314 111 Project , China [grant numbers D21001 ]. The strong motion data for this study are provided by Yunnan Earthquake Agency.

Appendix A Supplementary data to this article can be found online at https://doi.org/10.1016/j.heliyon.2024.e36177.
==== Refs
References

1 Aki K. Impact of earthquake seismology on the geological community since the Benioff zone Geol. Soc. Am. Bull. 100 1988 625 629 10.1130/0016-7606(1988)100<0625:IOESOT>2.3.CO;2
2 Ansal A. İyisan R. Yıldırım H. The cyclic behaviour of soils and effects of geotechnical factors in microzonation Soil Dynam. Earthq. Eng. 21 2001 445 452 10.1016/S0267-7261(01)00026-4
3 Cárdenas-Soto M. Chávez-García F.J. Regional path effects on seismic wave propagation in Central Mexico Bull. Seismol. Soc. Am. 93 2003 973 985 10.1785/0120020083
4 Hadley P.K. Askar A. Cakmak A.S. Subsoil geology and soil amplification in Mexico valley Soil Dynam. Earthq. Eng. 10 1991 101 109 10.1016/0267-7261(91)90040-7
5 Bard P.-Y. SESAME-team, Guideline for the implementation of the H/V spectral ratio technique on ambient vibrations-measurements Processing and Interpretations 2004 SESAME European Research Project EVG1-CT-2000-00026, D23.12
6 Arai H. Tokimatsu K. S-wave velocity profiling by inversion of microtremor H/V spectrum Bull. Seismol. Soc. Am. 94 2004 53 63 10.1785/0120030028
7 Herak M. ModelHVSR—a Matlab® tool to model horizontal-to-vertical spectral ratio of ambient noise Comput. Geosci. 34 2008 1514 1526 10.1016/j.cageo.2007.07.009
8 Bignardi S. Mantovani A. Abu Zeid N. OpenHVSR: imaging the subsurface 2D/3D elastic properties through multiple HVSR modeling and inversion Comput. Geosci. 93 2016 103 113 10.1016/j.cageo.2016.05.009
9 García-Jerez A. Piña-Flores J. Sánchez-Sesma F.J. Luzón F. Perton M. A computer code for forward calculation and inversion of the H/V spectral ratio under the diffuse field assumption Comput. Geosci. 97 2016 67 78 10.1016/j.cageo.2016.06.016
10 Guo Z. Aydin A. Huang Y. Xue M. Polarization characteristics of Rayleigh waves to improve seismic site effects analysis by HVSR method Eng. Geol. 292 2021 106274 10.1016/j.enggeo.2021.106274
11 Rong M. Wang J. Li X. Liu A. Kong X. Li H. Inversion method for velocity structure of soil lavers by ambient noise and its application Chin. J. Geophys. 66 2023 530 545 10.6038/cjg2022Q0165 (in Chinese)
12 Rong M. Li X. Fu L. Improvement of the objective function in the velocity structure inversion based on horizontal-to-vertical spectral ratio of earthquake ground motions Geophys. J. Int. 224 2020 1 16 10.1093/gji/ggaa347
13 Scherbaum F. Hinzen K.-G. Ohrnberger M. Determination of shallow shear wave velocity profiles in the Cologne, Germany area using ambient vibrations Geophys. J. Int. 152 2003 597 612 10.1046/j.1365-246X.2003.01856.x
14 Piña-Flores J. Perton M. García-Jerez A. Carmona E. Luzón F. Molina-Villegas J.C. Sánchez-Sesma F.J. The inversion of spectral ratio H/V in a layered system using the diffuse field assumption (DFA) Geophys. J. Int. 208 2017 577 588 10.1093/gji/ggw416
15 Perton M. Spica Z. Caudron C. Inversion of the horizontal-to-vertical spectral ratio in presence of strong lateral heterogeneity Geophys. J. Int. 212 2018 930 941 10.1093/gji/ggx458
16 Bora N. Biswas R. Malischewsky P. Imaging subsurface structure of an urban area based on diffuse-field theory concept using seismic ambient noise Pure Appl. Geophys. 177 2020 4733 4753 10.1007/s00024-020-02547-4
17 Thomas A.M. Spica Z. Bodmer M. Schulz W.H. Roering J.J. Using a dense seismic array to determine structure and site effects of the two towers earthflow in Northern California Seismol Res. Lett. 91 2020 913 920 10.1785/0220190206
18 Chen C.-T. Kuo C.-H. Lin C.-M. Huang J.-Y. Wen K.-L. Investigation of shallow S-wave velocity structure and site response parameters in Taiwan by using high-density microtremor measurements Eng. Geol. 297 2022 106498 10.1016/j.enggeo.2021.106498
19 Farazi A.H. Hossain MdS. Ito Y. Piña-Flores J. Kamal A.S.M.M. Rahman MdZ. Shear wave velocity estimation in the Bengal Basin, Bangladesh by HVSR analysis: implications for engineering bedrock depth J. Appl. Geophys. 211 2023 104967 10.1016/j.jappgeo.2023.104967
20 Field E.H. Hough S.E. Jacob K.H. Using microtremors to assess potential earthquake site response: a case study in Flushing Meadows, New York City Bull. Seismol. Soc. Am. 80 1990 1456 1480 10.1785/BSSA08006A1456
21 Dal Moro G. Some aspects about surface wave and HVSR analyses: a short overview and a case study Boll. Geofis. Teor. Appl. 52 2011 241 259 10.4430/bgta0007
22 Dal Moro G. Surface Wave Analysis for Near Surface Applications 2015 Elsevier Amsterdam ; Boston
23 Mohamed A. Lindholm C. Girgis M. Site characterization and seismic site response study of the Sahary area, South Egypt Acta Geodyn. Geomater. 12 2015 427 436 10.13168/AGG.2015.0032
24 Setiawan B. Jaksa M. Griffith M. Love D. Seismic site classification based on constrained modeling of measured HVSR curve in regolith sites Soil Dynam. Earthq. Eng. 110 2018 244 261 10.1016/j.soildyn.2017.08.006
25 Dal Moro G. On the identification of industrial components in the horizontal-to-vertical spectral ratio (HVSR) from microtremors Pure Appl. Geophys. 177 2020 3831 3849 10.1007/s00024-020-02424-0
26 Hillers G. Campillo M. Lin Y.-Y. Ma K.-F. Roux P. Anatomy of the high-frequency ambient seismic wave field at the TCDP borehole J. Geophys. Res. Solid Earth 117 2012 B06301 10.1029/2011JB008999
27 Xu Y. Luo Y. Methods of ambient noise-based seismology and their applications Chin. J. Geophys. 58 2015 2618 2636 10.6038/cjg20150803 (in Chinese)
28 Haubrich R.A. Munk W.H. Snodgrass F.E. Comparative spectra of microseisms and swell Bull. Seismol. Soc. Am. 53 1963 27 37 10.1785/BSSA0530010027
29 Tanimoto T. Excitation of microseisms Geophys. Res. Lett. 34 2007 L05308 10.1029/2006GL029046
30 Dunand F. Bard P.Y. Chatelain J.L. Guéguen P. Vassail T. Farsi M.N. Damping and Frequency from Randomec Method Applied to In-Situ Measurements of Ambient Vibrations: Evidence for Effective Soil Structure Interaction 2002 Elsevier Science Ltd London
31 Badreldin H. Abu El-Ata A. El-Hadidy M. Cornou C. Abd El-Aal A.E.-A.K. Lala A.M. Active and passive seismic methods for site characterization in Nuweiba, Gulf of Aqaba, Egypt Soil Dynam. Earthq. Eng. 172 2023 108002 10.1016/j.soildyn.2023.108002
32 Wathelet M. Chatelain J.-L. Cornou C. Giulio G.D. Guillier B. Ohrnberger M. Savvaidis A. Geopsy: a user-friendly open-source tool set for ambient vibration processing Seismol Res. Lett. 91 2020 1878 1889 10.1785/0220190360
33 Albarello D. Lunedei E. Combining horizontal ambient vibration components for H/V spectral ratio estimates Geophys. J. Int. 194 2013 936 951 10.1093/gji/ggt130
34 Withers M. Aster R. Young C. Beiriger J. Harris M. Moore S. Trujillo J. A comparison of select trigger algorithms for automated global seismic phase and event detection Bull. Seismol. Soc. Am. 88 1998 95 106 10.1785/BSSA0880010095
35 Konno K. Ohmachi T. Ground-motion characteristics estimated from spectral ratio between horizontal and vertical components of microtremor Bull. Seismol. Soc. Am. 88 1998 228 241 10.1785/BSSA0880010228
36 Picotti S. Francese R. Giorgi M. Pettenati F. Carcione J.M. Estimation of glacier thicknesses and basal properties using the horizontal-to-vertical component spectral ratio (HVSR) technique from passive seismic data J. Glaciol. 63 2017 229 248 10.1017/jog.2016.135
37 Del Gaudio V. Coccia S. Wasowski J. Gallipoli M.R. Mucciarelli M. Detection of directivity in seismic site response from microtremor spectral analysis Nat. Hazard. Earth Sys. 8 2008 751 762 10.5194/nhess-8-751-2008
38 Del Gaudio V. Wasowski J. Muscillo S. New developments in ambient noise analysis to characterise the seismic response of landslide-prone slopes Nat. Hazard. Earth Sys. 13 2013 2075 2087 10.5194/nhess-13-2075-2013
39 Del Gaudio V. Luo Y. Wang Y. Wasowski J. Using ambient noise to characterise seismic slope response: the case of Qiaozhuang peri-urban hillslopes (Sichuan, China) Eng. Geol. 246 2018 374 390 10.1016/j.enggeo.2018.10.008
40 Del Gaudio V. Zhao B. Luo Y. Wang Y. Wasowski J. Seismic response of steep slopes inferred from ambient noise and accelerometer recordings: the case of Dadu River valley, China Eng. Geol. 259 2019 105197 10.1016/j.enggeo.2019.105197
41 Del Gaudio V. Wasowski J. Hu W. Capone P. Venisti N. Li Y. Ambient noise and ERT data provide insights into the structure of co-seismic rock avalanche deposits in Sichuan (China) Bull. Eng. Geol. Environ. 80 2021 7153 7170 10.1007/s10064-021-02346-8
42 Hartzell S. Leeds A.L. Jibson R.W. Seismic response of soft deposits due to landslide: the mission peak, California, landslide Bull. Seismol. Soc. Am. 107 2017 2008 2020 10.1785/0120170033
43 Stolte A.C. Cox B.R. Lee R.C. An experimental topographic amplification study at Los Alamos national laboratory using ambient vibrations Bull. Seismol. Soc. Am. 107 2017 1386 1401 10.1785/0120160269
44 Ma N. Wang G. Kamai T. Doi I. Chigira M. Amplification of seismic response of a large deep-seated landslide in Tokushima, Japan Eng. Geol. 249 2019 218 234 10.1016/j.enggeo.2019.01.002
45 Kazemnia Kakhki M. Del Gaudio V. Rezaei S. Mansur W.J. Directional variations of site response in a landslide area using ambient noise analysis via Nakamura's and polarization-based method Soil Dynam. Earthq. Eng. 141 2021 106492 10.1016/j.soildyn.2020.106492
46 Sánchez-Sesma F.J. Rodríguez M. Iturrarán-Viveros U. Luzón F. Campillo M. Margerin L. García-Jerez A. Suarez M. Santoyo M.A. Rodríguez-Castellanos A. A theory for microtremor H/V spectral ratio: application for a layered medium: theory for microtremor H/V spectral ratio Geophys. J. Int. 186 2011 221 225 10.1111/j.1365-246X.2011.05064.x
47 Wapenaar K. Retrieving the elastodynamic Green's function of an arbitrary inhomogeneous medium by cross correlation Phys. Rev. Lett. 93 2004 254301 10.1103/PhysRevLett.93.254301
48 Sánchez-Sesma F.J. Modeling and inversion of the microtremor H/V spectral ratio: physical basis behind the diffuse field approach Earth Planets Space 69 2017 92 10.1186/s40623-017-0667-6
49 Bignardi S. The uncertainty of estimating the thickness of soft sediments with the HVSR method: a computational point of view on weak lateral variations J. Appl. Geophys. 145 2017 28 38 10.1016/j.jappgeo.2017.07.017
50 Stanko D. Markušić S. An empirical relationship between resonance frequency, bedrock depth and VS30 for Croatia based on HVSR forward modelling Nat. Hazards 103 2020 3715 3743 10.1007/s11069-020-04152-z
51 Dobry R. Oweis I. Urzua A. Simplified procedures for estimating the fundamental period of a soil profile Bull. Seismol. Soc. Am. 66 1976 1293 1321 10.1785/BSSA0660041293
52 Kramer S.L. Geotechnical Earthquake Engineering 1996 Prentice Hall Upper Saddle River, New Jersey
53 Ibs-von Seht M. Wohlenberg J. Microtremor measurements used to map thickness of soft sediments Bull. Seismol. Soc. Am. 89 1999 250 259 10.1785/BSSA0890010250
54 Gosar A. Lenart A. Mapping the thickness of sediments in the Ljubljana Moor basin (Slovenia) using microtremors Bull. Earthq. Eng. 8 2010 501 518 10.1007/s10518-009-9115-8
55 Rupar L. Mapping the thickness of Quaternary sediments in the Iska alluvial fan (central Slovenia) using microtremor method Acta Geodyn. Geomater. 2020 177 190 10.13168/AGG.2020.0013
56 Delgado J. López Casado C. Giner J. Estévez A. Cuenca A. Molina S. Microtremors as a geophysical exploration tool: applications and limitations Pure Appl. Geophys. 157 2000 1445 1462 10.1007/PL00001128
57 Delgado J. López Casado C. Estévez A. Giner J. Cuenca A. Molina S. Mapping soft soils in the Segura river valley (SE Spain): a case study of microtremors as an exploration tool J. Appl. Geophys. 45 2000 19 32 10.1016/S0926-9851(00)00016-1
58 Parolai S. Bormann P. Milkereit C. New relationships between Vs, thickness of sediments, and resonance frequency calculated by the H/V ratio of seismic noise for the Cologne Area (Germany) Bull. Seismol. Soc. Am. 92 2002 2521 2527 10.1785/0120010248
59 Hinzen K.-G. Weber B. Scherbaum F. On the resolution of H/V measurements to determine sediment thickness, a case study across a normal fault in the Lower Rhine Embayment, Germany J. Earthq. Eng. 8 2004 909 926 10.1080/13632460409350514
60 Garcia-Jerez A. Characterization of the sedimentary cover of the zafarraya basin, southern Spain, by means of ambient noise Bull. Seismol. Soc. Am. 96 2006 957 967 10.1785/0120050061
61 Motamed R. Ghalandarzadeh A. Tawhata I. Tabatabaei S.H. Seismic microzonation and damage assessment of Bam city, Southeastern Iran J. Earthq. Eng. 11 2007 110 132 10.1080/13632460601123164
62 D'Amico V. Picozzi M. Baliva F. Albarello D. Ambient noise measurements for preliminary site-effects characterization in the urban area of Florence, Italy Bull. Seismol. Soc. Am. 98 2008 1373 1388 10.1785/0120070231
63 Birgören G. Özel O. Siyahi B. Bedrock depth mapping of the Coast South of İstanbul: comparison of analytical and experimental analyses Turk. J. Earth Sci. 2009 10.3906/yer-0712-3
64 Özalaybey S. Zor E. Ergintav S. Tapırdamaz M.C. Investigation of 3-D basin structures in the İzmit Bay area (Turkey) by single-station microtremor and gravimetric methods: 3-D basin structures of the İzmit Bay area Geophys. J. Int. 186 2011 883 894 10.1111/j.1365-246X.2011.05085.x
65 Sukumaran P. Parvez I.A. Sant D.A. Rangarajan G. Krishnan K. Profiling of late Tertiary–early Quaternary surface in the lower reaches of Narmada valley using microtremors J. Asian Earth Sci. 41 2011 325 334 10.1016/j.jseaes.2011.02.011
66 Poggi V. Fäh D. Burjanek J. Giardini D. The use of Rayleigh-wave ellipticity for site-specific hazard assessment and microzonation: application to the city of Lucerne, Switzerland: Rayleigh wave for micronization Geophys. J. Int. 188 2012 1154 1172 10.1111/j.1365-246X.2011.05305.x
67 Paudyal Y.R. Yatabe R. Bhandary N.P. Dahal R.K. Basement topography of the Kathmandu Basin using microtremor observation J. Asian Earth Sci. 62 2013 627 637 10.1016/j.jseaes.2012.11.011
68 Johnson C. Lane J. Statistical comparison of methods for estimating sediment thickness from Horizontal-to-Vertical Spectral Ratio (HVSR) seismic methods: an example from Tylerville, Connecticut, USA Symposium on the Application of Geophysics to Engineering and Environmental Problems 2016 2016 Society of Exploration Geophysicists and Environment and Engineering Geophysical Society Denver, Colorado, USA 317 323 10.4133/SAGEEP.29-057
69 Maresca R. Berrino G. Investigation of the buried structure of the Volturara Irpina Basin (southern Italy) by microtremor and gravimetric data J. Appl. Geophys. 128 2016 96 109 10.1016/j.jappgeo.2016.03.010
70 Tün M. Pekkan E. Özel O. Guney Y. An investigation into the bedrock depth in the Eskisehir Quaternary Basin (Turkey) using the microtremor method Geophys. J. Int. 207 2016 589 607 10.1093/gji/ggw294
71 Sant D.A. Parvez I.A. Rangarajan G. Patel S.J. Bhatt M.N. Sanoop Salam T.A. Subsurface profiling along Banni Plains and bounding faults, Kachchh, Western India using microtremors method J. Asian Earth Sci. 146 2017 326 336 10.1016/j.jseaes.2017.06.002
72 Molnar S. Cassidy J.F. Castellaro S. Cornou C. Crow H. Hunter J.A. Matsushima S. Sánchez-Sesma F.J. Yong A. Application of microtremor horizontal-to-vertical spectral ratio (MHVSR) analysis for site characterization: state of the art Surv. Geophys. 39 2018 613 631 10.1007/s10712-018-9464-4
73 Mascandola C. Massa M. Barani S. Albarello D. Lovati S. Martelli L. Poggi V. Mapping the seismic bedrock of the Po Plain (Italy) through ambient‐vibration monitoring Bull. Seismol. Soc. Am. 109 2019 164 177 10.1785/0120180193
74 Thabet M. Site-Specific Relationships between bedrock depth and HVSR fundamental resonance frequency using KiK-NET data from Japan Pure Appl. Geophys. 176 2019 4809 4831 10.1007/s00024-019-02256-7
75 Kumar P. Mahajan A.K. New empirical relationship between resonance frequency and thickness of sediment using ambient noise measurements and joint-fit-inversion of the Rayleigh wave dispersion curve for Kangra Valley (NW Himalaya), India Environ. Earth Sci. 79 2020 256 10.1007/s12665-020-09000-8
76 Nelson S. McBride J. Application of HVSR to estimating thickness of laterite weathering profiles in basalt Earth Surf. Process. Landf. 44 2019 1365 1376 10.1002/esp.4580
77 Tuan T.T. Vinh P.C. Malischewsky P. Aoudia A. Approximate formula of peak frequency of H/V ratio curve in multilayered model and its use in H/V ratio technique Pure Appl. Geophys. 173 2016 487 498 10.1007/s00024-015-1098-6
78 He Z. Hu G. Lu L. Zhang W. Ye T. Shen K. The shallow velocity structure for the Tonghai basin in Yunnan Chin. J. Geophys. 56 2013 3819 3827 10.6038/cig20131123 (in Chinese)
79 Brocher T.M. Empirical Relations between Elastic Wavespeeds and density in the earth's crust Bull. Seismol. Soc. Am. 95 2005 2081 2092 10.1785/0120050077
80 Shen W. Ritzwoller M.H. Crustal and uppermost mantle structure beneath the United States J. Geophys. Res. Solid Earth 121 2016 4306 4342 10.1002/2016JB012887
81 Caielli G. De Franco R. Di Fiore V. Albarello D. Catalano S. Pergalani F. Cavuoto G. Cercato M. Compagnoni M. Facciorusso J. Famiani D. Ferri F. Imposa S. Martini G. Paciello A. Paolucci E. Passeri F. Piscitelli S. Puzzilli L.M. Vassallo M. Extensive surface geophysical prospecting for seismic microzonation Bull. Earthq. Eng. 18 2020 5475 5502 10.1007/s10518-020-00866-4
82 Albarello D. Herak M. Lunedei E. Paolucci E. Tanzini A. Simulating H/V spectral ratios (HVSR) of ambient vibrations: a comparison among numerical models Geophys. J. Int. 234 2023 870 878 10.1093/gji/ggad109
83 Foti S. Parolai S. Albarello D. Picozzi M. Application of surface-wave methods for seismic site characterization Surv. Geophys. 32 2011 777 825 10.1007/s10712-011-9134-2
84 Anderson J.G. Quantitative measure of the goodness-of-fit of synthetic seismograms International Association for Earthquake Engineering 2004
