A review of near-surface QS estimation methods using active and passive sources

Seismic attenuation and the associated quality factor (Q) have long been studied in various sub-disciplines of seismology, ranging from observational and engineering seismology to near-surface geophysics and soil/rock dynamics with particular emphasis on geotechnical earthquake engineering and engineering seismology. Within the broader framework of seismic site characterization, various experimental techniques have been adopted over the years to measure the near-surface shear-wave quality factor (QS). Common methods include active- and passive-source recording techniques performed at the free surface of soil deposits and within boreholes, as well as laboratory tests. This paper intends to provide an in-depth review of what Q is and, in particular, how QS is estimated in the current practice. After motivating the importance of this parameter in seismology, we proceed by recalling various theoretical definitions of Q and its measurement through laboratory tests, considering various deformation modes, most notably QP and QS. We next provide a review of the literature on QS estimation methods that use data from surface and borehole sensor recordings. We distinguish between active- and passive-source approaches, along with their pros and cons, as well as the state-of-the-practice and state-of-the-art. Finally, we summarize the phenomena associated with the high-frequency shear-wave attenuation factor (kappa) and its relation to Q, as well as other lesser-known attenuation parameters.


3
Vol:. (1234567890) conditions were conducted by George Heinrich Otto Volger (1822Volger ( -1897 and Robert Mallet (1810Mallet ( -1881. Volger (1856Volger ( , 1858 described qualitatively the 1855 Visp (Switzerland) earthquake and produced a map showing areas characterized by equal seismic intensity (isoseismal map) with the collaboration of the cartographer August Petermann . Not long after, the Great Neapolitan earthquake of 1857 was described by Mallet (1862), who also developed a similar map. These studies were followed by Lawson (1908) on the great 1906 California earthquake, as well as Borcherdt (1970) and Borcherdt and Gibbs (1976), which related the effects of local geological and geotechnical characteristics, and in particular V S and Q S , on observed or recorded ground motions. Anderson et al. (1996) notably stated that site effects have an enormous influence on the character of ground motions, despite site conditions (in terms of their characteristic dimensions) rarely represent more than 1% of the source-to-site distance. By extension, a comprehensive understanding of the regional nearsurface conditions is also necessary to accurately assess the characteristics and spatial variability of the earthquake-induced ground motions.
In the last decades, most attention has been given to the development of methods that allowed reliable and robust estimates of the shear-wave velocity (V S ) beneath a site, by applying active (e.g., linear seismic reflection, surface waves, and borehole investigations) and/or passive sources (e.g., surface arraybased or single-station). The former methods (e.g., Xia et al. 2002a, b;Williams et al. 2003;Onnis et al. 2019) take advantage of the knowledge of the source location and of the receivers defined ad hoc for the experiment. The latter methods are either based on analysis of earthquake recordings at permanent or temporary seismological stations (e.g., Lermo and Chavez-Garcia 1993;Field and Jacob 1995;Lachet et al. 1996;), or on the use of seismic ambient noise (e.g., Aki 1957;Fäh et al. 2003;Okada 2003;Scherbaum et al. 2003;Tokimatsu 2004, 2005;Parolai et al. 2005;Köhler et al. 2007;Boxberger et al. 2011;Foti et al. 2011). In urban areas, active source methods, known mainly to yield reliable data in the high-frequency range, can be adversely affected when high-amplitude anthropogenic noise contamination occurs. These approaches are also associated with relatively higher operational costs. Passive-based approaches-in particular, those based on the seismic noise acquisition and analyseshave advantages, mainly that they are lower in cost than active source options, as well as for their noninvasive nature and short data acquisition times. More recently, site investigators, motivated by the need for accurate modeling of V S profiles and derivation of the average V S values in the first 30 m (V S30 ), have combined the aforementioned active and passive source methods (Richwalski et al. 2007;Odum et al. 2013;Yong et al. 2013Yong et al. , 2019Martin et al. 2021).
A complete assessment of site effects, however, can be obtained only if the effect of near-surface attenuation is also included in the site response analysis (e.g., Lai and Rix 1998;Parolai 2012). Site attenuation, as is understood today in a soil dynamics and ground motion framework, counteracts amplification effects that are primarily due to impedance contrasts within the soil column, unless in cases of extreme non-linear behavior and liquefaction (e.g., Fiegel and Kutter 1994;Bonilla et al. 2005). The absolute amplitude of ground shaking in the presence of local site resonance-and especially for the higher modeseven for well-defined velocity profiles, depends also rather strongly on the degree of attenuation in the materials of the soil profile. At relatively low frequencies (e.g., 1-10 Hz), attenuation effects have been studied in detail with reference to soil stiffness and damping ratio curves in fairly soft materials (e.g., Régnier et al. 2018).
Although extensive seminal work exists on soil/ rock material properties with respect to their attenuation and anisotropy (Barton 2006), it is only recently that interests have grown about their effects on ground motion prediction and scaling at higher frequencies (> 10 Hz) (Ktenidou and Abrahamson 2016). The quantities used to describe material attenuation may differ across different disciplines and frequency ranges, from t* (Solomon 1973;Singh et al. 1982;Cormier 1982) and kappa, or κ (Anderson and Hough 1984) in seismological terms to the quality factor (Q) (Futterman 1962;Knopoff 1964a, b) in geophysical notation. However, there is a general consensus that these quantities describe an overall decay of the amplitude of ground motion due to two physical mechanisms: (1) scattering by heterogeneities encountered by the waves along the seismic path, which is considered mostly frequency-dependent, and (2) intrinsic damping or anelastic attenuation due to internal friction within the material, which is often approximated (in the frequency range of engineering seismology interests, i.e., mainly from 0.1 to 10-20 Hz) as frequency-independent (i.e., hysteretic). However, there is still an ongoing debate regarding the frequency dependence of Q, in particular, on whether it is part of the intrinsic attenuation or if it is only related to the propagation (scattering) in the medium (e.g., Singh et al. 1982;Morozov 2009). This work focuses on Q or its geotechnical engineering counterpart D or ξ (damping ratio), while a short discussion of its relation to κ, t* is provided in the last Section of this article.
Despite the importance of attenuation in determining the modification of shaking caused by the wave propagation in the shallowest geological layers, the assessment of site attenuation, specifically that of the shear-wave quality factor (Q S ), has attracted less attention. This is probably due to difficulties in constraining damping effects with accuracy from seismic data. Several studies have relied on borehole earthquake recordings (e.g., Redpath and Lee 1986;Hauksson et al. 1987;Seale and Archuleta 1989;Assimaki et al. 2006Assimaki et al. , 2008Kinoshita 2008;Parolai et al. 2010), although other attempts were also successful when using active seismic sources generating body and surface waves (e.g., Wang et al. 1994;Rix and Lai 1998;Xia et al. 2002a, b;Foti 2004;Haase and Stewart 2005;Badsar et al. 2010;Xia 2014;Gao et al. 2018). Prieto et al. (2009) reported that it was possible to estimate the one-dimensional (1D) Q S structure at a regional scale from seismic noise. Albarello and Baliva (2009) deduced the average damping in the soil from seismic noise at a local scale. Other studies followed in this direction both at regional and local scales (e.g., Harmon et al. 2010;Tsai 2011;Weemstra et al. 2013;Parolai 2014;Magrini and Boschi 2021). Recently, Dikmen et al. (2016) and Haendel et al. (2019), based on previous approaches, analyzed seismic noise data in boreholes. Spatial resolution, intended as the characteristic surface length of soil from which Q S is estimated, and the range of frequencies investigated by most of the methods analyzed in this work are shown in Fig. 1.
In this paper, we provide an overview of the approaches proposed in the literature for estimating Q S in the near-surface layers (here, the term near-surface refers in general to tens to hundreds of meters depth from the surface, although this may be extended to a few km in the case of deep sedimentary basins) using surface and borehole arrays that record active and passive sources. Thus, this manuscript does not attempt to provide a full review for alternative attenuation models such as κ, t*, or other types of quality factors (e.g., coda Q). We also do not address related studies applied in laboratory settings. However, we briefly mention laboratory tests because they represent a key benchmark in the measurement of Q. A thorough review of publications on the attenuation is an enormous undertaking. We, nevertheless, present an earnest attempt at describing the state-ofknowledge. The paper was thought to be consulted if only for single Sections, in case the reader wants to focus just on one (sub)Chapter and not on the full manuscript.
We begin by defining Q. Then, we describe Q S approaches based on the analysis of signals generated by active sources and collected by sensors installed both at the Earth's surface and in boreholes. Next, we provide an overview of methods for estimating Q S based on the analyses of earthquakes in boreholes and seismic noise, the latter recorded by arrays both at the surface and in boreholes; we also present advantages/disadvantages therein. Finally, we discuss prospects in current and ongoing developments for improving methods for determining Q S .

Measures of energy dissipation in linear attenuating media
The mechanical response of geomaterials to lowstrain dynamic excitations is studied in different yet interacting disciplines such as geophysics, seismology, and geotechnical engineering. It follows that each of these disciplines has independently developed different terminologies and technical words, although referring to exactly the same physical phenomenon. This is the case for the energy loss suffered by a mechanical disturbance propagating in a geological material, for which various parameters have been introduced (e.g., O'Connell and Budiansky 1978;Ishihara 1996;Aki and Richards 2002). Some of them were inspired by quantities adopted in different scientific fields such as circuit theory in 1 3 Vol:. (1234567890) electrical engineering (Cole and Cole 1941). Most of these parameters are dimensionless and proportional to the ratio between the energy dissipated by a unit volume of geomaterial in one cycle of harmonic excitation and a measure of the corresponding stored energy. In essence, the proposed parameters attempt to provide a normalized estimate of the internal entropy density production, which is an indicator of the amount of energy dissipated by a unit volume of soil or rock mass undergoing a cyclic deformation. Although the energy-related interpretation of Q is independent from any formulation of constitutive modeling of material behavior, a link to the latter is greatly desired especially in computational geophysics. The linear theory of viscoelasticity is the simplest constitutive model to satisfactorily capture the most salient aspects observed in the mechanical response of geomaterials undergoing low-amplitude dynamic oscillations, which is the capability to store and simultaneously to dissipate strain energy over a finite period of time. Indeed, experimental evidence shows that geomaterials tend to exhibit a linear, yet inelastic response when subjected to dynamic excitations at strain levels below the linear cyclic threshold strain (Vucetic 1994;Kramer 1996). It has also been shown that the shape of the stress-strain loop predicted by the theory for a general viscoelastic material undergoing harmonic oscillations is elliptical (Pipkin 1986), a feature that compares well with the experimental data in geomaterials at low strain (Dobry 1970). A widely Fig. 1 Different spatial resolution and frequency ranges covered by the approaches described in this paper for measuring attenuation and the related shear-wave quality factor (Q S ). The Multi-channel Analysis of Surface Waves (MASW) and Spectral Analysis of Surface Waves (SASW) methods are noninvasive active techniques, based on the measurement of the dispersion properties of Rayleigh surface waves for estimating V S and Q S in the subsurface materials. SASW uses one active source and two sensors, whereas MASW adopts one active source and a linear array with a variable number of sensors (geophones) 1 3 Vol.: (0123456789) adopted definition of Q, which is hereinafter denoted by Q A , is the following: where ω is the angular frequency, E max is the maximum stored energy per unit volume of geomaterial during a cycle of harmonic excitation and ∆E diss is the amount of energy (per unit volume) dissipated by the material during the same cycle of loading (see Fig. 2). The latter is also equal to the amount of entropy produced due to unrecoverable mechanical work. A similar definition is used in geotechnical engineering for material damping ratio D A , which is related to Q A by the relation: Within the time domain, the constitutive relationships of an isotropic, linear, viscoelastic material are represented by a pair of integro-differential equations known as the Boltzmann's equations, while in the frequency domain they assume a simple algebraic form that resembles Hooke's law of linear elasticity. The only difference is that the two elastic constants, G and M, representing the shear and the constrained moduli, respectively, are replaced by the corresponding complex-valued shear and constrained (or oedometric) moduli Ĝ S ( ) = G S1 ( ) + iG S2 ( ) and Ĝ P ( ) = G P1 ( ) + iG P2 ( ) , respectively. If the linear theory of viscoelasticity is applied to Eq. (2.1), the quantity ∆E diss (ω) can be simply expressed for each mode of deformation in terms of the imaginary part (named the loss modulus) of the complex modulus, either Ĝ S ( ) or Ĝ P ( ) . However, the numerator of Eq. (2.1), namely the maximum stored energy E max (ω), depends not only on the real and the imaginary parts of the corresponding complex modulus, but also on their derivatives with respect to frequency (Tschoegl 1989). This is due to the phase lag existing among various energy storing mechanisms governing the response of viscoelastic materials during a harmonic excitation. As a result, when Q defined by Eq. (2.1) is expressed in terms of the material functions Ĝ S ( ) or Ĝ P ( ) , the ensuing result is an awkward expression that is inconvenient to use. The difficulties associated with the application of viscoelasticity theory to Eq. (2.1) are overcome if Q is defined as follows (Dain 1962;O'Connell and Budiansky 1978): where E ave is the average stored energy per unit volume of geomaterial during a cycle of harmonic excitation. This definition of Q, hereinafter denoted as Q B , has the advantage that E ave can be expressed, for each mode of deformation, in terms of the real part only (named the storage modulus) of either the complex-valued shear modulus Ĝ S (ω) or the complexvalued constrained (or oedometric) modulus Ĝ P (ω). Therefore, with the definition (2.3), Q B is linked to the material functions of linear viscoelasticity by the following simple relation: where G S1 and G P1 are the real parts of the shear and constrained complex moduli, respectively, while G S2 and G P2 are the corresponding imaginary parts, respectively. Equation (2.4) is sometimes used as a formal definition of Q, since the energetic . 2 Hysteretic loop of a cyclic stress-strain curve with definition of the quality factor Q A . The x and y axes are the engineering shear strain (γ) and the shear stress (τ), respectively. The hysteretic loop with the elliptical shape colored in light brown represents the energy dissipated in the geomaterial during one cycle of harmonic excitation and has an area equal to ∆E diss . The colored triangle is the maximum strain energy stored during that cycle and is defined by A ∆ = E max . Q is the quality factor and is obtained from the ratio of E max and ∆E diss multiplied by the constant 2π 1 3 Vol:. (1234567890) interpretation provided by Eq. (2.3) has been shown to hold only for a viscoelastic material that can be represented by a network of linear elastic springs and viscous dashpots (O'Connell and Budiansky 1978).
The quantity Q −1 , as defined by Eq. (2.4), is referred to in the literature as the loss tangent, since it is the ratio between the imaginary and the real parts of the complex modulus. Thus, it is the tangent of the loss angle which is the argument of either Ĝ S (ω) or Ĝ P (ω). The loss angle represents the phase angle by which the strain lags the stress in a steady-state harmonic motion due to internal friction.
It is worth remarking that when the energy losses in the geomaterial are small, the definitions (2.1), (2.3), and (2.4) of Q provide nearly the same results and the loss angle can be approximated by the loss tangent. It has been shown (e.g.,  that the energy losses may indeed be considered small when Q −1 is smaller or equal than 0.10 (Q ⩾ 10). In these cases, the distinction among the various definitions of Q become irrelevant, and the material is named weakly dissipative or low-loss or loss-less medium. Experimental evidence shows that nearsurface soil deposits (e.g., Vucetic and Dobry 1991;Ishibashi and Zhang 1993;Ishihara 1996) strained below the linear cyclic threshold strain (for normally consolidated clays with a plasticity index of 50, this parameter is of the order of 10 −4 ) exhibit values of Q higher than 10. For deep subsurface geomaterials, the threshold strain is larger (Vucetic 1994;Kramer 1996) and Q is significantly higher (e.g., Menq 2003).
Over time, other definitions of Q have been proposed in the literature, with some attaining moderate success. With reference to waves traveling in an attenuating medium, Futterman (1962) for instance proposed a definition of Q similar to Eq. (2.1), but with E max and E diss replaced by the peak kinetic energy density T max observed at a point and the drop of T max over a spatial wavelength, respectively. Aki and Richards (2002) distinguished between temporal and spatial Q, depending on the method used to measure this energy loss parameter. Temporal Q is determined from monitoring the decay in time of the amplitude of a standing wave of prescribed wavenumber, whereas spatial Q is obtained from the spatial decay of the amplitude of a wave of prescribed frequency. In the spatial Q definition, it is assumed that the direction of maximum attenuation coincides with the direction of propagation. In this case, for plane waves propagating in homogeneous low-loss media, there is no difference between temporal and spatial Q. Note that, when the direction of propagation is coincident with the direction of maximum attenuation, "simple" or "homogeneous" waves (Lockett 1962) are generated. "Non-simple" waves may arise as a result of boundary effects (e.g., reflection and refraction of monochromatic waves at a plane interface) combined with particular types of viscoelastic models (Christensen 2010). From a mathematical viewpoint, temporal Q is represented by a complex-value frequency, whereas spatial Q requires the wavenumber to be a complexvalue. These results originate from the application of either the Laplace or the Fourier transform, respectively, to the Boltzmann constitutive equations of linear viscoelasticity.
It is, however, conceptually relevant to distinguish the definition of a material parameter (Q) from its experimental determination. Clearly, Q defined by either Eq. (2.3) or Eq. (2.4) is a constitutive parameter (or better, a material function) of a specific material model. Thus, determining Q becomes a parameteridentification problem and different experimental techniques may be conceived to solve it in geophysics, seismology, and geotechnical engineering.

The relation between Q and material dispersion imposed by physical causality
In the theory of linear viscoelasticity, for each mode of deformation (say, e.g., shear), only one material function is required in the time domain to completely specify the mechanical response of the material. This may be, for instance, the shear relaxation G S (t), which is a real-valued function that represents the stress response (in shear) of a material subjected to a strain history specified as a Heaviside function. Since in the frequency domain the material response is described by the complex shear modulus, the real and imaginary parts G S1 (ω) and G S2 (ω) of Ĝ S (ω) cannot be independent. Their dependence can be established by applying the Fourier integral theorem to Boltzmann's equations after assuming the relaxation function G S (t) to obey at the principle of non-retroactivity, which states that in a viscoelastic material the stress at the time t is caused by a strain history that acted only in the past and not in the future (Fung 1965). This implies G S (t) = 0 for t < 0 and, if this condition is satisfied, the function G S (t) is said to be causal. The

3
Vol.: (0123456789) same concept applies also to G P (t), where P denotes P-waves. It can be shown that the Fourier transform of a causal, real-valued function G S (t) is a complex-valued function Ĝ S (ω), the complex shear modulus, where the real and the imaginary parts are Hilbert transforms pairs (Tschoegl 1989). Formally, this functional dependence between G S1 (ω) and G S2 (ω) or between G P1 (ω) and G P2 (ω) is represented by the so-called Kramers-Kronig integral equation (Christensen 2010) in honor of Ralph de Laer Kronig (1926) and Hans A. Kramers (1927), who studied the relation between dispersion and absorption of electromagnetic waves.
In a linear viscoelastic medium, Ĝ P (ω) and Ĝ S (ω) may be replaced by the complex-value velocity of propagation V P ( ) and V S ( ) , to which they are linked by the following expressions: where ρ is the mass density of the material that has been assumed homogeneous. For a harmonic viscoelastic plane S-wave exp( t −k S ⋅ x) propagating along the direction x , the complex-value wavenumber k S ( ) = V S ( ) can be decomposed as follows: where V S (ω) and α S (ω) are the real-valued physical velocity of propagation and the attenuation coefficient of the S-wave, respectively. A relation identical to Eq. (2.6) can be written for P-waves, as well. Equation (2.6) suggests two important remarks. The first is the frequency-dependence of V S (ω) and α S (ω) inherited by the constitutive parameter Ĝ S (ω). Thus, a viscoelastic material is inherently dispersive, i.e., the speed of propagation of viscoelastic waves is necessarily frequency-dependent (Futterman 1962). The second observation is that V S (ω) and α S (ω) are a pair of material functions alternative to the real and imaginary parts of the complex modulus Ĝ S (ω). The same applies also to V P (ω) and α P (ω) with respect to Ĝ P (ω). Equations (2.5) and (2.6), combined with the functional relationship between the real and the imaginary parts G S1 (ω) and G S2 (ω) of the complex modulus Ĝ S (ω), suggest that the Kramers-Kronig relation admits an alternative representation where G S1 (ω) and G S2 (ω) are replaced by V S (ω) and α S (ω) (Aki and Richards 2002): where V S (0) = lim →0 V S ( ) . A similar expression could be derived for P-waves after replacing Q S (ω) and V S (ω) with Q P (ω) and V P (ω). Equation (2.9) provides an explicit expression of dispersion-attenuation pairs for arbitrary dissipative, linear viscoelastic materials. By measuring the frequency-dependence of V S (ω), the first equation can be used to calculate the Q S (ω) spectrum. Alternatively, from the second equation dispersion V S (ω) can be calculated if measuring the Q S (ω) spectrum. Figure 3 shows the dispersion curves corresponding to a Gaussian and Rayleigh Q S (ω) spectra. It should be noted that, although both Q S (ω) spectra have a symmetrical shape, the corresponding dispersion curves are anti-symmetrical (as expected from the theory).
In seismology, a well-known, particular solution of the Kramers-Kronig relation is obtained assuming that Q S (ω) is hysteretic, namely rate-independent over the seismic bandwidth (i.e., 0.001-100 Hz). For shear waves, a dispersion relation consistent with the hysteretic assumption and based on the Hilbert transform pair proposed by Azimi et al. (1968) can be written as follows (Liu et al. 1976;Kjartansson 1979;Kennett 1983;Keilis-Borok 1989;Aki and Richards 2002): where ω ref denotes a reference angular frequency usually assumed equal to 2π. Equation (2.10) is better known as the Kolsky attenuation-dispersion model and has also been proposed for P-waves. Dispersion relation (2.10) predicts a V S that decreases monotonically with Q S for a fixed frequency. On the other hand, for a particular value of Q S , it predicts a log-linear increase of V S (ω) with frequency (high-frequency waves travel faster) as shown in Fig. 4. Equation (2.10) is often postulated a priori not only in seismology (e.g., Lomnitz 1957;Knopoff 1964a, b;Liu et al. 1976;Sipkin and Jordan 1979), but also in soil dynamics, owing to the rate independence of Q S (ω) exhibited by soils (and rocks) at low strains in the seismic bandwidth (Shibuya et al. 1995;Lo Presti and Pallara 1997;Ling et al. 2007; (2.10) Santamarina 2007;Tatsuoka et al. 2008;Tatsuoka 2009;Assimaki et al. 2012). However, other studies yielded opposite results and seem to demonstrate that in geomaterials Q S (ω) is sensitive to the loading frequency even in the seismic band of seismology (e.g., Murphy 1982;Spencer 1981;Berckhemer et al. 1982;Jackson et al. 1992Jackson et al. , 2002Satoh 2006) and geotechnical engineering (e.g., Kim et al. 1991;Leroueil and Marques 1996;Lin et al. 1996;d'Onofrio et al. 1999;Darendeli 2001;Matesic and Vucetic 2003;Rix and Meng 2005;Zambelli et al. 2007;Araei et al. 2012). Thus, the experimental results obtained so far are controversial and they do not allow to draw a definitive conclusion upon the presumed hysteretic/nonhysteretic nature of Q S (ω) in geomaterials at low strain. Fundamentally, this way of posing the problem is conceptually incorrect. The frequency dependence or independence of Q S (ω) should not be a matter to be postulated a priori. Indeed, the Q S (ω) spectrum is a material function that should be experimentally measured to fully characterize the response of geomaterials. Once Q S (ω) is known, the frequency-dependent V S (ω) can be computed by using Eq. (2.9). Dispersion curves calculated with Eq. (2.9). Influence of concavity of Q S (ω) spectrum on the computed dispersion curves (modified from Lai and Özcebe 2016b) Dispersion-attenuation relation (2.10) is a simple algebraic equation, and it is compatible with the physical principle of causality. However, a constant Q S over the entire frequency range ∈ [0, +∞[ would imply a frequency-independent shear (or compression) wave velocity and this would violate causality, since no Hilbert transform pair may satisfy Kramers-Kronig Eq. (2.7) with a constant Q S (Aki and Richards 2002). Therefore, frequency-dependence should be considered outside the seismic band even for a hysteretic Q S .
As an alternative to Eq. (2.10), sometimes Q S (ω), based on the good match of simulated versus observed ground motion recordings (e.g., Michaels 2006; Kawase 2019), is postulated to obey at specific frequencies the results predicted by such models as the Kelvin-Voigt model or the standard linear solid. However, the a priori enforcement of specific and oversimplified viscoelastic rheologies in the interpretation of attenuation measurements not only corresponds to assuming a specific distribution of relaxation times, which may not be supported by experimental data, but it may even lead to violation of physical causality (e.g., with the Kelvin-Voigt model). Assumptions like hysteretic Q S or other types of rheologies represent unnecessary constraints (Moczo and Kristek 2005) or even speculations regarding the interpretation of the mechanical response of geomaterials during wave propagation without a valid experimental validation. The good match of the numerical simulations with the ground motion recordings by itself is not a guarantee of virtuous modeling, since a parameter-estimation problem is known to be illposed because it suffers from the non-uniqueness of the solution; namely, the same set of ground motion recordings may in principle be obtained with different dispersion-attenuation models.

Geotechnical approaches in measuring Q
In geotechnical engineering, it is standard practice for near-surface site characterization to measure Q S (actually D S , namely damping ratio in the shear mode of deformation, Eq. (2.2) shows the relation between D S and Q S ) in a laboratory test conducted on a small soil specimen using a particular device known as the resonant column apparatus (ASTM D4015-15e1). A solid, cylindrical soil specimen is subjected to harmonic excitation by an electromagnetic driving system applied at the top of the sample (Drnevich 1985). Sometimes, hollow specimens are used to obtain a more uniform distribution of the shear strain amplitude along the cross-section of the sample.
The soil specimen can be excited in either the torsional or the longitudinal mode of vibration. Thus, the equipment allows to measure Q S and Q P at low strain (< 10 −5 ). Soils are barotropic materials (i.e., their mechanical response is sensitive to the intensity of the confining stress), thus prior to the Q measurement with the resonant column apparatus, the undisturbed soil sample is consolidated to reproduce the lithostatic pressure acting at the sampling depth. This is a standard practice in experimental soil mechanics when conducting laboratory tests. A typical configuration of the torsional resonant column apparatus is shown in Fig. 5, with the soil specimen fixed at the base and free to rotate at the top where the driving torque is applied.
Q is then obtained via the free-vibration decay method, where the amplitude reduction of successive oscillations exhibited by the soil specimen is monitored after setting the sample into free vibration. Alternatively, Q may be determined through the half-power bandwidth technique, which exploits the dependency of the shape of the frequency response curve on the magnitude of the energy loss. The same equipment is also used to determine the V S of the soil sample, although applying a different experimental procedure which is inconsistent with that adopted for determining Q. In fact, V S is measured assuming soil behaves in a linear elastic manner. Furthermore, the frequency dependence of V S (i.e., material dispersion) and Q S is disregarded, and this leads to a violation of physical causality (Lai and Özcebe 2016a, b). Note that when Q is less than 30 (Futterman 1962), the dispersion for V S starts to become not negligible. Nonconventional use of the resonant column apparatus to determine the frequency dependence of V S and Q S in undisturbed soil specimens has, however, been proposed through a direct measurement of the complexvalue shear modulus Ĝ S ( ) = G S1 ( ) + iG S2 ( ) Rix and Meng 2005;Khan et al. 2008).
One drawback that occurs when using the resonant column apparatus is the equipment-generated damping due to the driving system. Typically, in a resonant column device or in a combined resonant column and torsional shear apparatus, the driving system is constituted of a coil and a magnet. The magnet moving in the coil generates an electromotive force (back electromotive force) that is opposed to the driving movement. If not considered, this effect may influence damping measurement in soils and thus the comparison of laboratory versus in situ test results. For this reason, research has been conducted to determine the magnitude of the equipment-generated damping (e.g., Wang et al. 2003;Meng and Rix 2003;Sasanakul and Bay 2010). Based on the outcomes of these studies, nowadays, in advanced geotechnical laboratories, the problem has been eliminated or at least mitigated. Yet, it is a matter of concern when referring to damping values measured with old-fashion resonant column equipment.

Influence of the mode of deformation: Q S and Q P
The energy dissipated by geomaterials in harmonic excitations, and thus Q, vary not only with the amplitude of the induced strain, but also with the mode of deformation: for instance, shear mode (Q S ), com- Particularly relevant in geophysics and seismology are the differences, if any, between Q S (ω) and Q P (ω). Ultimately, these differences are connected to distinctive types of dissipation mechanisms that are activated when a S-rather than a P-wave propagates through a porous medium.
Using linear viscoelasticity, Toksöz and Johnston (1981) developed theoretical models for Q S (ω) and Q P (ω) based on assuming specific attenuation mechanisms associated to wave propagation in porous geomaterials, including friction, fluid-flow, viscous relaxation, and wave scattering. The deformation process in coarse-grained geomaterials involves internal rearrangements of the particles and contact re-orientation. In fluid-filled porous media and cracked rock, the phenomenon is further complicated by solid-fluid interaction and by the magnitude of the confining pressure. The energy loss occurring during P-and S-wave propagation is a function of the degree of saturation. Winkler and Nur (1982) made an interesting experimental study aimed at investigating the attenuation of P-and S-waves in rocks due to frictional sliding and fluid-flow mechanisms. They thoroughly analyzed the effects on attenuation of confining pressure, pore fluid, degree of saturation, strain amplitude, and frequency. The shear and constrained complex moduli of linear viscoelasticity were used as material functions, and Eq. (2.4) was adopted as a definition of Q. Partial water saturation significantly increases the energy loss of both P-and S-waves relative to that in dry rock, with Q P resulting lower than Q S . When complete saturation conditions prevail, Q P is greater than Q S . The ratio Q P /Q S is found to be a more sensitive indicator of partial gas saturation than the corresponding V P /V S ratio. In saturated rocks at low strain, "squirt-flow" (Mavko and Nur 1975) is the predominant loss mechanism of solid-fluid interaction (Mavko and Nur 1978;Palmer and Traviolia 1980;Dvorkin et al. 1995;Santamarina et al. 2001;Pride et al. 2004;Gurevic et al. 2010). The frequency of excitation plays an important role, as different relaxation times are associated to distinct yet interacting dissipation mechanisms and this occurs differently in shear and longitudinal modes of deformation. The geophysical literature is particularly rich of studies in this regard, especially in connection to material modeling, starting from the pioneering work of Biot who developed a theory of linear poroelasticity and visco-poroelasticity to describe the propagation of low-frequency and high-frequency waves in fluid-saturated porous continua (Biot 1956a, b). In the original articles, Biot did not provide details on how to determine the constitutive parameters of his model which includes Q P and Q S ; however, he did so in later work (Biot and Willis 1957). Since the Biot's theory accounts only for a single mechanism of energy dissipation, when compared with experimental data, this model tends to underestimate dispersion and attenuation, particularly at higher frequencies. Dvorkin and Nur (1993) developed a linear visco-poroelastic model that improved upon Biot's theory by accounting for the squirt-flow mechanism, which is connected to the heterogeneity of the porous medium at the microscopic scale of the individual pores. This combined model accounting for the Biot mechanism and squirt-flow is known in the literature as the BISQ theory (e.g., Nie and Yang 2008;Nie et al. 2012).
In Table 1, we report some measured values of Q P and Q S in near-surface deposits found in literature. Q P and Q S values in sandstones by Winkler and Nur (1982) are included in order to allow a comparison with the values obtained in unconsolidated soils.   3 Active source methods

Surface wave analysis
The spatial attenuation of surface waves at the free surface of a soil deposit, propagating away from its source, is caused in part to geometrical spreading, which depends on the source geometry, and in part to the inelasticity of the medium and, thus, to Q. Consequently, Love and Rayleigh waves can in principle be used to estimate Q in geomaterials, although this application of surface wave methods is not yet fully included in the state-of-the-practice of nearsurface site characterization (Foti et al. 2018). This is despite surface wave techniques offering several advantages for Q measurements over the borehole methods, including that they are economical because they are (by nature) non-invasive, and the adopted frequencies are better correlated to those in earthquake site response analyses when compared to those used in borehole tests. Also, the effects of poor coupling between the soil and the receiver, which may adversely affect the measurement of particle motion amplitudes, are easier to identify at the ground surface than inside a borehole. In seismology, it is customary to determine the velocity and Q structure of the Earth's crust and mantle from the inversion of long period dispersion and attenuation curves (e.g., Anderson and Archambeau 1964;Solomon 1975, 1979;Cheng and Mitchell 1981;Herrmann 2013). Surface waves generated by earthquakes have periods ranging from several seconds to 100 s and more, which allow Q and the velocity structure to be calculated at depths of tens if not hundreds of kilometers. Seismologists and geophysicists have adjusted their methods used in studying the Earth's structure to shallower depths, from near-surface up to few kilometers for hydrocarbons exploration (e.g., Mokhtar et al. 1988;Jongmans 1990;Malagnini et al. 1995). These studies used explosive sources to generate transient Rayleigh waves that have been recorded at large distances (up to 100 km). Also, it is worth mentioning the contribution by Li et al. (1995) who, under the assumptions of hysteretic Q and weak attenuations, used the dispersive and absorption properties of Love waves to estimate Q of a country rock and a coal seam.
Regarding the near-surface geophysical-geotechnical characterization, Spang (1995) set up a methodology for estimating hysteretic Q S of a soil deposit from the inversion of an attenuation curve of Rayleigh waves determined from Spectral Analysis of Surface Waves (SASW) measurements. The analyses were conducted assuming weak dissipation. Furthermore, the theoretical attenuation curve was associated with the fundamental mode of Rayleigh wave propagation, and when accounting for geometric spreading the effects of geometric dispersion were neglected. Continuing this pioneering work, Rix and Lai (1998) and Lai and Rix (1998) proposed a method for the simultaneous inversion of experimental Rayleigh phase velocity and attenuation curves to compute V S and the quality factor Q S of a soil deposit using an electro-mechanical shaker as a source, which allowed frequency control. The simultaneous (i.e., coupled) inversion of surface wave data is deemed superior over the corresponding uncoupled analysis because it considers the inherent coupling existing between Q and the velocity of propagation of seismic waves due to material dispersion (Section 2). The theory of linear viscoelasticity coupled with the elegant methods of complex analysis were used by the authors as the theoretical framework of their methodology. The inversion of the experimental dispersion and attenuation curves was conducted by accounting for multi-mode wave propagation through the concept of apparent (or effective) dispersion and attenuation curves, which is compatible with the use of harmonic sources (Lai et al. 2014). To improve the efficiency of the algorithm, elements of the Jacobian matrix involving the partial derivatives of the apparent Rayleigh phase velocity with respect to the medium parameters required for the solution of the Rayleigh inverse problem were computed analytically.
The uncoupled approach, for which the measurement and inversion of surface-wave attenuation data to obtain Q S are separate from the measurement and inversion of Rayleigh-wave velocity data to compute V S , was applied by Rix et al. (2000) at the US Geotechnical Experimentation Site of Treasure Island where independent laboratory measurements of Q S were available for comparison. As a further advancement of this methodology, Rix et al. (2001) and  developed a test procedure to measure and invert surface wave dispersion and attenuation data simultaneously from a single set of acquisitions. Experimentally, the method is based on measuring the experimental transfer function from the traces recorded at the receivers and the signature of the source represented by an electro-mechanical shaker. The proposed technique introduced consistency in surface wave methods between phase velocity and attenuation measurements by using the same sourcereceiver configuration array, a step followed by the implementation of a joint inversion. Foti (2003) and Foti (2004) further enhanced the method by defining a transfer function from the deconvolution of the seismic traces recorded at the receivers, thereby eliminating the need of the source's signature which potentially may introduce uncertainties for the poor coupling between the source and the underlying subsurface.
The successful inversion of Rayleigh attenuation data to estimate Q was also undertaken by Xia et al. (2002a, b) and Xia et al. (2012) using an uncoupled, constrained inversion algorithm (i.e., damped leastsquares) based on a trade-off between data misfit and roughness of the vector of model parameters represented by the dissipation factors of soil layers. Besides Q S , when the ratio V S /V P is greater than 0.45-a situation which is common in the oil industry and occasionally encountered in near-surface soil deposits-the authors also attempted to estimate Q P . In determining Q P , it was found that most contributions to Rayleigh-wave attenuation coefficients are in a relatively high frequency range, while contributions from Q S are in a lower frequency range. In computing the Rayleigh-wave attenuation coefficients, the correction of the measured amplitudes for cylindrical divergence (i.e., geometric attenuation) was made without considering geometric dispersion.
A joint estimation of Q S and V S from the inversion of dispersion and attenuation surface wave data was also proposed by Misbah and Strobbia (2014) using an array of possibly unevenly spaced receivers, coupled with a set of multiple shots at different locations. The method is based on measuring the complex-value Rayleigh wavenumber of multiple normal modes. The interference among the various modes was considered without the need of a priori specifying a multimode Rayleigh spreading function, so that multiple modal attenuation curves could be extracted. The possibility of considering multiple shot points into a receiver array also allowed the increase of the modal resolution without increasing the receiver aperture. The feasibility of estimating the near-surface quality factor Q S in a layered ground model from the inversion of the Love waves' attenuation coefficients was studied by Xia et al. (2013) and Xia (2014). Attenuation coefficients of Love waves are independent of Q P , which makes the inversion of attenuation coefficients of Love waves to estimate Q S simpler than that of Rayleigh waves, since fewer parameters are required.
Despite the encouraging results of the aforementioned studies, an accurate measurement of both Rayleigh and Love attenuation coefficients remains a challenge, given there are no straightforward procedures to separate material and geometric attenuation taking also into account the role played by the different modes of propagation. Also, it should be remarked that since attenuation measurements rely on accurate estimates of the amplitude of particle motion, it is essential that accurate vertical emplacements and physical coupling of each receiver are checked carefully in advance because the uncertainty of attenuation measurements contributes to the uncertainty of the Q S structure (Rix et al. 2001;Foti et al. 2018). Badsar et al. (2010) and Badsar et al. (2011) proposed an alternative approach for determining Q S in shallow soil layers from surface wave measurements that was not based on the spatial decay of Rayleighwave amplitudes, but on transforming the surface wave field into a frequency-wavenumber spectrum. Then, the modal attenuation curves are derived from the width of the f-k spectrum (f and k denote the cyclic frequency and circular wavenumber, respectively) using the half-power bandwidth method. The latter is a technique widely used in the fields of engineering structural dynamics (Clough and Penzien 1993) to determine the modal damping ratio of structures idealized as multi-degree of freedom systems (with widely spaced resonance frequencies) from the shape of properly defined transfer functions. The method does not require the calculation of geometric attenuation, and despite the fact that Q S is determined by inverting the attenuation curve associated to the fundamental Rayleigh mode of propagation, the authors claim that in irregular soil profiles their technique yields more accurate results than methods based on computing the attenuation coefficients by separating the effects of geometric from intrinsic attenuation. This is because the occurrence of higher Rayleigh modes does not affect the fundamental mode attenuation curve. Yet, they admit that their method relied on the use of a few critical parameters whose values significantly affect the results. Thus, the authors suggest that their most appropriate values should be determined from a parametric study.
Finally, among the most recent studies is that by Gao et al. (2018), who, based on numerical examples and a real case study, proposed to estimate the attenuation coefficients based on separating the attenuation curves associated with different modes of propagation of both Love and Rayleigh waves. The Q S structure is then calculated by inverting the attenuation curve associated to the fundamental mode. The authors claim that in irregular soil deposits where multimode wave propagation is important, prior separation of the various attenuation coefficients provides more accurate estimation of the measured attenuation coefficients and, correspondingly, of the Q S structure. However, in this study, the attenuation coefficients are computed without considering geometric dispersion.
Borehole investigations provide the opportunity to assess the low-strain quality factor Q free from the undesirable effects of sample disturbance affecting laboratory measurements. At the same time, they involve a larger, more representative volume of soil when compared with the small-sized specimens used in lab experiments. However, the measurement of Q in a soil deposit via borehole methods is not a standard practice when compared to measuring the velocity structure, and most of the work conducted so far was conceived within research projects. The obtained results have only been moderately successful, mainly due to the difficulties of separating geometric and intrinsic attenuation, but also because of the influence of wave scattering. Furthermore, most of the adopted techniques assume a frequency-independent Q in a certain frequency band and in their measurements the authors do not always distinguish between Q P and Q S .
In the above references, Q was computed from borehole measurements using a variety of techniques associated to either cross-hole or downhole data (i.e., vertical seismic profiling) including amplitude decay with distance, spectral slope, spectral ratio, wavelet and phase modelling, matching technique, inverse Q −1 filtering, and pulse rise-time. None of these methods were clearly demonstrated to be superior. Some of them appeared more suitable than others depending on soil stratigraphy, depth, seismic source, level of ambient noise, and other factors.
More recent contributions include the work by Hall and Bodare (2000), who used cross-hole seismic tests with a configuration of four aligned boreholes to compute the dispersion function V S (ω) from the phase of cross-power spectra using the classical two-station method. The latter was also used to estimate Q after assuming rate-independent (i.e., hysteretic) dissipation and the geometric spreading law proportional to r −1 , with r being the distance from the source. The testing site was located in Sweden and the investigated soil consisted of a thick (75 m) soft clayey layer. The methods of determining V S (ω) and Q were verified through numeral simulations conducted with the finite element technique. On average, the value of Q in the clayey layer resulted to be about 20 (Q S −1 about 0.05). Pujol et al. (2002) studied the attenuation of S-waves at three sites in Arkansas and Tennessee (USA) using data recorded in boreholes up to a depth of 60 m and a source represented by a compressedair-driven hammer. The lithology encountered was typical of fluvial, floodplain deposits with interbedded layers of medium to coarse sands with lenses of sandy clay, sandy silt, clayey silt, and silty clay of variable thickness. The attenuation of wave amplitude was estimated using the spectral ratio technique with a fixed depth and variable frequency (higher than 10-20 Hz). The values of Q S resulting from the study ranged from 18 to 44 (Q S −1 from 0.02-0.03 to 0.06) depending on the depth range analyzed. These values are smaller than typical values reported in the literature for these geomaterials, which are of an order of around 100 (Q S −1 = 0.01). An important conclusion from the study is also that low velocity does not necessarily imply low Q S , as is generally assumed. This result has been confirmed recently by Boore et al. (2020) and is further discussed below.

3
Vol.: (0123456789) Michaels (1998) and (2006) proposed a methodology to determine low-strain Q S (ω) and Q P (ω) based on simultaneously inverting dispersion and attenuation curves of P-and S-waves in downhole seismic testing, assuming soils behaving according to the Kelvin-Voigt rheological model. The amplitude decay with distance was corrected with a geometric spreading law assumed proportional to r −1 . The fieldtesting site was located in Colorado (USA) and the soil deposit comprised a partially saturated silty sand layer of about 4.5-m thickness overlying a gravelly sand which extended to a depth of about 10 m. The author of the study found Q P (ω) lower than Q S (ω), suggesting the explanation that while in shear (i.e., Q S ) inertial coupling is the key dissipation mechanism, in compression (i.e., Q P ) this mechanism is diffusion. Thus, a low-density pore fluid like air can result in high levels of dissipation for P-waves, but not for SH-waves. On the other hand, a dense pore fluid is required to increase dissipation in shear. The frequency dependence of the dissipation factors determined in this experiment is biased by the assumed Kelvin-Voigt rheology.
Using Seismic Cone Penetration Test (SCPT) data, Karl et al. (2006) estimated the variation of Q S with depth by applying the spectral ratio method to the Fourier transforms of the horizontal acceleration time histories that have been recorded during the experiment at two testing sites in Belgium. The soil deposits investigated at the first site consisted of a 10-m-thick layer of clayey silt to silty clay, while the second site was characterized by a layer of about 10 m thick of fine sand to silty sand. The obtained results are characterized by a significant scatter with a coefficient of variation exceeding 60%. At one of the two sites, Q S at low strain resulted unexpectedly low with some values in the range between 5 and 10 (Q S −1 between 0.10 and 0.20). At the other site, the average value was about 25 (Q S −1 = 0.04) and was confirmed by independent laboratory measurements. A major assumption made by the authors in their analyses was frequency-independent Q S in the frequency range of interest.
More recently, Crow et al. (2011a) performed downhole tests in Ontario (Canada) to measure dispersion functions of S-waves and low-strain Q S in soft clayey silts (V S < 250 m/s) using a frequency-controlled vibratory source in the range 10-100 Hz. The dissipation factor Q S was obtained using the spectral ratio method with measured values ranging from 125 to 250 (Q S −1 from 0.004 to 0.008), thus indicating that the investigated soil is a very low loss material at low strains. Results also showed negligible frequency dependence of Q S between 10 and 100 Hz, thereby supporting the assumption of hysteretic Q S at least in the investigated frequency band. Average values of Q S measured in the laboratory by Crow et al. (2011b) with the resonant column apparatus were lower (Q S = 66.7; Q S −1 = 0.015), perhaps due to sample disturbance. Also, some frequency dependence was noticed both in Q S and V S in the range between 65 and 85 Hz.
The interdependency of Q S (ω) and V S (ω) stated by Eqs. (2.9), representing the solution of Kramers-Kronig relation, is the framework for setting up a procedure of determining Q S (ω) and Q P (ω) from the measurement of the dispersion functions V S (ω) and V P (ω), respectively. A preliminary application to demonstrate the feasibility of the method was performed by Lai and Özcebe (2012) and Lai and Özcebe (2016a, b) who attempted to estimate Q S (ω) from borehole cross-hole seismic data. Although the preliminary results obtained at two sites in Italy were encouraging, there are technical difficulties that need to be overcome, the most important of which concerns the limited frequency bandwidth at which the dispersion measurements were made. This prevented a broadband characterization of V S (ω) and, thus, of Q S (ω). Yet, the method offers several advantages over the conventional cross-hole test based on picking the first arrivals of P and S phases. Firstly, Q S (ω) is obtained from the inversion of the dispersion function V S (ω) calculated using the full waveforms of the recorded seismograms. Consequently, physical causality of V S (ω) and Q S (ω) is automatically fulfilled, as Q S (ω) is determined from the exact solution of the Kramers-Kronig relation. Secondly, in interpreting the borehole measurements, no a priori assumption is made concerning specific rheological behaviors like the controversial hystereticity of Q S or the presumed single relaxation time enforced by the Kelvin-Voigt model. With this method, the estimated Q S and the shape of the corresponding function Q S (ω) is solely determined by the shape of the dispersion curve V S (ω), as predicted by the Kramers-Kronig relation. Thirdly, the method replaces the measurement of the amplitude characteristics of seismic signals, which is 1 3 Vol:. (1234567890) critical due to the influence of geometrical spreading with the measurement of velocity dispersion.
Using a frequency-controlled S-wave source vibrator, Koedel and Karl (2020) were able to estimate Q S from multi-channel spectral analysis of seismic downhole data. The authors estimated the shear dissipation factor by fitting a hysteretic (i.e., rate-independent) dispersion model to an experimentally measured velocity dispersion function V S (ω). The latter is extracted from a phase velocity-frequency spectrum calculated using the same f-k methodology adopted in Multi-channel Analysis of active Surface Wave (MASW) testing. The test site is located in the city of Hannover (Germany). The downhole experiment was conducted in a three-layer soil deposit consisting of fine and coarse sands down to 6-m depth overlying an intermediate gravelly layer between depths of 6 and 19 m and a cretaceous chalky claystone down to the final examined depth of 100 m. Although no independent laboratory data are available at the testing site, the obtained results (e.g., average Q S of 20 at depths of 5-35 m and of 7.7 at depths of 35-90 m) are comparable to the experimental data reported in the literature for similar soils.
Boore et al. (2020) applied two approaches discussed by Gibbs et al. (1994) to extract Q S from downhole receiver measurements using a repeatable shear-wave surface source. One approach adopted the ratio of Fourier spectra at two different depths, whereas the other used the change in Fourier spectral amplitudes as a function of depth for individual frequencies. The procedure was applied at 22 boreholes in the San Francisco Bay area in central California and the San Fernando Valley of southern California to determine the average Q S over depth intervals ranging from about 10 to 245 m (at one site), with most maximum depths being between 35 and 90 m. The average Q S values ranged from 6.25 to 50 (Q S −1 from less than 0.02 to almost 0.16) with little dependence on soil type and grain size for sites in sediments. Interestingly, the average Q S values for sites with average V S greater than about 450 m/s, including, but not limited to, rock sites, turned out to be generally lower than for sites with lower average value of V S . This was already observed by others, including Pujol et al. (2002) and Boxberger et al. (2017). Although the obtained values of Q S are generally compatible with low-strain dissipation factors from laboratory measurements which are in the range of 25-50 (Q S −1 in the range 0.02-0.04), no specific laboratory tests were performed at the sites for confirmation. The estimate of Q S versus depth was made at frequencies greater than 10 Hz. Considering the experimental procedure, it is uncertain how the measured values of Q S will actually reflect the combined effect of scattering and intrinsic attenuation due to the inelasticity of the medium. This, however, is true for most indirect surface or borehole methods, since only laboratory techniques can eliminate the effects of wave scattering.

Earthquake recordings in boreholes
Some early studies investigated regional earthquakes by means of a single sensor at depth in order to analyze wave attenuation. For instance, Rautian et al. (1978) examined the spectral content of P-and S-waves at three stations placed in bedrock within tunnels from 15-to 30-m depth. They retrieved Q P and Q S by leveraging the high seismicity of the area considered and different rocks to the North and South of which the region of research was placed. From their results, they stated that attenuation was not determinant for the modeling of the ratio of the S-to-P-wave spectra, as their differences were inconsequential.
Stork and Ito (2004) determined a frequency-independent Q with the spectral analysis by considering the best-fitting Brune's source model (Brune 1970) in one 800-m-deep borehole sensor. To check the validity of assuming a frequency-independent Q, they estimated frequency-dependent Q P and Q S and attenuation model for small-magnitude-range events by means of the extended coda normalization method of Yoshimoto et al. (1993). This last methodology considers that the coda spectral amplitude is proportional to the spectral amplitude of P-and S-waves. The Q model obtained following Yoshimoto et al. (1993) showed that it was not always appropriate to assume a frequency-independent Q on a wide frequency range.
Vertical arrays of seismometers or strong motion sensors provide recordings of earthquake signals from different depths and at the surface, allowing (in principle) an in situ estimation of the medium's characteristics over the frequency range of engineering 1 3 Vol.: (0123456789) interest. Joyner et al. (1976) compared observed and predicted spectral ratios between different depth levels in a four-level 186-m-deep borehole, by considering a viscoelastic and linear Haskell plane-layer model. They assumed a frequency-independent Q and found good agreement when considering surface and bedrock spectra. From spectral ratios between surface and bedrock recordings, Archuleta (1986) noted that the near-surface sediment amplification was almost homogeneous in the frequency range 2-50 Hz. Malin et al. (1988) considered velocity spectra and particle motion of ultra-microearthquakes as functions of depth and found attenuation of S-waves at high frequencies near the surface. Attempts to retrieve mean Q P and Q S values by Blakeslee and Malin (1991) from spectral ratios in the frequency range 1-80 Hz highlighted that Q could also be fully masked by both amplification in the near-surface and resonance effects. They observed a strong oscillatory behavior when Q estimations were conducted considering only low-frequency (< 24 Hz) recordings. This determined an uncertain assessment of the spectral ratio decay rate and resulted in Q values that could not consider, for instance, the lack of a high-frequency free-surface reflection in the downhole signal. Therefore, they concluded that it was important not to estimate Q from low-frequency recordings only, due to the significant role of low-Q materials at higher frequencies.
Alternatively, fitting the high-frequency part of the spectral ratio (f > 20 Hz), which was assumed to be less affected by down-going reflected phases to a first order, by simple exponential decay was proposed by Aster and Shearer (1991). Jongmans and Malin (1995) examined a 938-m-deep borehole. They estimated both a section of the 1D Q S profile between 298-and 938-m depth and the mean Q S for the entirety of the borehole depth by considering the average slope of the spectral ratio above 20 Hz. They observed that in the depth interval from 0 to 298 m, attenuation could not be calculated from a simple constant Q S model of the spectral ratio, since scattering effects of internal low-velocity zones, as well as interference, increased the shallow high-frequency signal. Thus, as the sensor depth reduced, a marked increase in amplitude was generated. Abercrombie and Leary (1993) found Q P and Q S from the comparison between uphole and 2.5-kmdeep downhole recordings. They observed that at the surface seismic amplitudes were attenuated at higher frequencies compared to downhole, whereas at lower frequencies they were amplified. They carried out a fit for the downhole earthquake spectra to a standard Brune's source model (Brune 1970) and observed that downhole spectra did not need to be completely corrected for intrinsic attenuation. Abercrombie (1997) estimated a 1D Q S profile by calculating the spectral ratio between surface and deep recordings, avoiding considering the interference between the up-going wave and the surface reflection at the borehole site. Abercrombie (1998) reviewed studies on attenuation from earthquake analyses in boreholes (specifically, Abercrombie and Leary 1993;Malin et al. 1988;Aster and Shearer 1991;Blakeslee and Malin 1991;Jongmans and Malin 1995;Abercrombie 1997). These studies relied on a frequency-independent near-surface Q. Abercrombie (1998) suggested that the attenuation between two borehole sensors and the frequency interval to calculate the gradient of the spectral ratio determine the accuracy of Q values obtained with the spectral ratio method. The author specified that while Q can be efficiently estimated by means of the spectral ratio technique for frequencies greater than 10 Hz, for frequencies below 10 Hz, Q estimation is difficult due to the dominance of nearsurface attenuation. This is because resonance effects, as well as amplifications, become pervasive. She also indicated that in spectral ratios the problem of surface and shallow layers' reflections recorded from borehole sensors can be either corrected by means of simple reflectivity codes (e.g., Steidl 1996), or avoided by carefully selecting window lengths (Abercrombie 1997). Safak (1997) highlighted that down-going waves reflected at the surface might affect the downhole recordings, especially for shallow boreholes. Therefore, the simple spectral ratio method, applied generally to the surface and the downhole recordings, cannot lead to a robust estimation of Q S because of the transmission coefficients at the boundaries of the layers. To overcome this drawback, when possible (i.e., for a deep enough borehole sensor), the spectral ratio is taken between the up-going and down-going pulses in the downhole seismogram (e.g., Kinoshita 1983Kinoshita , 2008Hauksson et al. 1987;Fukushima et al. 1992).
Attenuation studies by Archuleta et al. (1992Archuleta et al. ( , 1993 in a downhole array were followed by computations of synthetic seismograms in the frequency band 0-10 Hz by Bonilla et al. (2002), and their comparisons with those recorded from a small nearby earthquake. P-, S-wave velocities, and Q P and Q S profiles were calibrated through a trial-and-error procedure to gain the best fit in both the time and amplitude domains between synthetics and observed data, by referring to geotechnical and seismic surveys. Similarly, Assimaki et al. (2006) and Assimaki et al. (2008) proposed an inversion procedure that estimates the best borehole model in terms of shear-wave velocity, attenuation, and mass density, by optimizing the correlation between observed and synthetic seismograms.
Furthermore, under the condition that the orientation of the sensor is correctly known, Q S might be estimated by an inversion procedure that optimizes the fit either between the observed and the calculated, for a certain model, amplitude spectral ratios (Seale and Archuleta 1989) or between the observed and theoretical temporal propagator for a layered medium (Trampert et al. 1993). Parolai et al. (2010) proposed to estimate the attenuation in the shallow geological layers considering the deconvolution of the wavefield recorded in a borehole with that recorded at the surface. The first method required the Fourier transform of the deconvolved wavefield to be fitted with a theoretical transfer function valid for the vertical or nearly vertical propagation of S-waves. The second method was based on the spectral fitting of the Fourier transform of only the acausal part of the deconvolved wavefield with a theoretical transfer function. Both methods work more accurately when the medium between the two sensors does not present large vertical velocity contrasts. To overcome this limitation, Parolai et al. (2012) proposed a linear inversion of the spectra of a deconvolved wavefield for estimating the Q S structure below a site that also considers the V S velocity profile. However, they indicated that while the V S profile can be well constrained, the Q S is less reliably constrained. Parolai et al. (2013) applied an inversion on the deconvolved wavefield between the sensors inside the borehole based on the approaches of Parolai et al. (2010Parolai et al. ( , 2012, by using the non-stationary ray decomposition of Kinoshita (2009). It should be remarked that all methods based on the wavefield deconvolution are assuming a linear transfer function and are, therefore, suitable only for linear behavior of soils. Both Parolai et al. (2012) and Parolai et al. (2013) argued that the number of layers of the subsoil model should be consistent with the number of the borehole sensors used, regardless of the local stratigraphic features. Raub et al. (2016) forward modeled the deconvolved seismograms in time domain for obtaining the effective Q P and Q S values, instead of deconvolving the wavefield in frequency domain as Parolai et al. (2010Parolai et al. ( , 2012. Standard spectral ratio techniques were not applicable, due to the strong interference effects between up-going and down-going waves they reported. This complexity, together with the assumption of a single homogeneous layer above each sensor in their model, led them to apparent Q values. Although Q S included both intrinsic and scattered attenuation, as well as impedance effects, it was similar to the results of Parolai (2010). Fukushima et al. (2016) retrieved Q S in sediments by combining the methods of Fukushima et al. (1992) and Trampert et al. (1993). Up-going (incident) and down-going (surface-reflected) waves were separated by the deconvolution of the seismogram recorded at the bottom of the borehole with the seismogram from the ground surface. Q S values were, then, obtained from the transfer function between up-going and down-going waves. They applied the combined method to boreholes deeper than 300 m, with a theoretical S-wave two-way travel time larger than 0.5 s. Riga et al. (2019) applied the method of Fukushima et al. (2016), both in its original form and with some modifications, to synthetic tests. They concluded that the technique of Fukushima et al. (2016) could be enhanced to improve the accuracy in estimating the attenuation, if up-going and down-going wave windows were selected based on travel times rather than on signals' coherence measures, as originally proposed. They also showed that the method could be potentially extended to shallower borehole arrays. Their results were valid for downhole sensors placed both above the bedrock interface and, for a certain frequency range, at the soil-bedrock interface (or close to it). For real cases, Riga et al. (2019) suggested carrying out simulations first to determine the best analyzable frequency band.
Recently, Seylabi et al. (2020) used a sequential data assimilation method (Evensen 2009) based on the Ensemble Kalman Inversion (EKI) (Iglesias et al. 2013) for obtaining soil V S and damping ratio from the joint inversion of dispersion data as input. They first refined the method by means of synthetic data analyses, then applied it to real data from a downhole array. The algorithm allowed them to set a high number of layers for increasing the precision. They showed that inverting dispersion data together with acceleration time series can improve the results.

Seismic noise arrays
In the last decade, several studies showed that it is possible to retrieve the attenuation of Rayleigh waves, and thus, that of the shear waves, starting from the analysis of seismic noise data collected from arrays of seismometers. Most studies have proposed either a modification of the SPatial Auto-Correlation analysis (SPAC) introduced by Aki (1957), or the analysis of the results of the interferometric approach (Claerbout 1968). Note that in both cases the effect of the phase difference between two signals is estimated either looking at their frequency and space dependent correlation, when the wavelengths are much larger than the SPAC array dimension, or at the estimated band-limited Green's function, when the analyzed wavelengths are much shorter than the interstation distances (Fig. 6). Prieto et al. (2009) showed that it was possible, at a regional scale, to estimate the attenuation of surface waves using seismic noise recordings. They estimated attenuation by demonstrating that the spatial coherency of the ambient seismic field is related to the Green's function. They also inverted data from the network seismic stations for the 1D Q S structure, based on the assumption that Q P /Q S is equal to 2. Later, Lawrence and Prieto (2011) and Prieto et al. (2011) adapted and extended the method of Prieto et al. (2009) to calculate lateral variations in the attenuation structure by means of attenuation tomography. Prieto et al. (2011) highlighted that the attenuation measures gained using the approach of Prieto et al. (2009) may be contaminated by some factors. Among them, there is the focusing and defocusing of Rayleigh waves, which can change amplitudes when travelling through 3D heterogeneous media (e.g., Dalton and Ekström 2006), and the non-uniformity of the source of seismic noise that can affect both phase velocity and attenuation with a certain percentage (Harmon et al. 2010). Prieto et al. (2011) tried to reduce the effects of focusing during the measures. Regarding the second contamination, they explained that this was minimum for a yearlong coherency stack and the results of Prieto et al. (2009) could be well overlapped with those obtained by Yang and Forsyth (2008) analyzing surface waves with longer periods from earthquake recordings. Lin et al. (2011) demonstrated that the spatially averaged attenuation observed with ambient noise and regional seismic event measurements observed from data recorded by the US Transportable Array were consistent. The cross-correlation with spectral whitening used by them was similar to the application of coherency supported by Prieto et al. (2009). Tsai (2011), in a theoretical framework study based on Tsai (2009Tsai ( , 2010, demonstrated that the derivation of the attenuation proposed by Prieto et al. (2009) was appropriate only for uniformly distributed noise sources when homogeneous attenuation was assumed. In fact, this did not Fig. 6 3D model representing the distance between a pair of seismic stations in relation to the wavelength (λ) of a seismic wave (blue color) for both SPAC (orange color) and seismic interferometry (green color) methodologies. Usually, in the SPAC analysis the stations' separation is much shorter than λ (< < λ), and vice versa in seismic interferometry (> > λ) 1 3 Vol:. (1234567890) work for far-field isotropic noise sources. It follows that without a perfectly diffusive seismic noise field, it is not easy to constrain the attenuation (Harmon et al. 2010;Tsai 2011). Nakahara (2012) showed that the expression proposed by Prieto et al. (2009) was a good approximation already for the ratio between the attenuation factor κ and the wavenumber k 0 . Nakahara (2012) also observed a restricted validity of the exponential decay model, and the substantial influence of the direction of noise sources on attenuation values. Similar observations were acquired from numerical simulations and experiments by Cupillard and Capdeville (2010), Cupillard et al. (2011), Weaver (2011), and Liu and Ben-Zion (2013. Other numerical simulations were conducted by Weaver (2013). Weaver (2013) estimated attenuation coefficients and site amplification factors by employing the radiative transfer equation to both causal and anti-causal side, and the stationary phase approximation (Snieder 2004) from a linear array of no less than five stations. Lawrence et al. (2013) showed that, for a wide range of noise source distributions, the coherency of the noise correlation functions matches a Bessel function decaying exponentially with a specific attenuation coefficient. Zhang and Yang (2013), starting from the correlation of the coda of correlation (C 3 ) method developed by Stehly et al. (2008) and the temporal flattening method of Weaver (2011), estimated attenuation using noise and compared the results with earthquake data. They demonstrated that the C 3 method reduced bias and improved the attenuation evaluation for noise around the primary microseisms peak (18 s). On the contrary, less reliable results were gained by them for noise at the secondary microseisms peak (8 s). Weemstra et al. (2013) estimated the attenuation and Q of surface waves using recordings from an array with an aperture of several kilometers but did not attempt any 1D Q S inversions. They also examined pros and cons of attenuation estimations by means of seismic noise, similarly to Prieto et al. (2011). Weemstra et al. (2014) highlighted an inconsistency in the required proportionality constant of Weemstra et al. (2013) and reinterpreted the results of those authors. Weemstra et al. (2014) adopted the model proposed by Tsai (2011) to get the difference between the spectral whitening after cross-terms removal and the whitening before it. They found that for inter-receiver distances greater than two wavelengths, the resultant coherency was proportional to the complex coherency, if cross-terms were not deleted before spectral whitening. For minor distances, conversely, the prompt decay of coherency due to cross-terms in the normalization process could have led to erroneous shallow attenuation interpretations. Menon et al. (2014) showed that due to the slowness inhomogeneity, azimuthally averaging the coherence is not equivalent to homogenizing the medium. In fact, the former introduces apparent attenuation in the coherence because of interference. Bowden et al. (2015) showed that a source term added to the 2D wavefront tracking approach of Lin et al. (2012) could resolve spatially both intrinsic attenuation and scattering of a dense array if the seismic noise field is reasonably omni-directional. Liu et al. (2015) relied on theoretical constructs to apply conditions on amplitude decay due to attenuation to form a linear least-square inversion for interstation Rayleigh-wave quality factor values in southern California. Allmark et al. (2018) estimated both phase velocity and Q tomography at the Ekofisk Field by implementing the methods of Bloch and Hales (1968) and Liu et al. (2015). Q values gained by Allmark et al. (2018) were not dominated by the phase velocities estimates and showed that Q and compression were non-linearly correlated. They found that some Q structures could be related with local known geological characteristics. Boschi et al. (2019) developed a mathematical expression for the multiplicative factor that correlates the normalized cross correlations to the Rayleigh-wave Green's function, and the obtained value, equal to 1, confirmed the speculation of Prieto et al. (2009). A numerical validation of Boschi et al. (2019Boschi et al. ( , 2020)'s method for quantifying the attenuation of Rayleigh waves from the cross-correlation of real seismic noise data was conducted by Magrini and Boschi (2021). They first validated the algorithm through synthetic tests, then applied it successfully to broad-band recordings. Their results were similar to those of previous studies in the same area and were also reliable when data did not support a perfectly diffuse wavefield. They pointed out that source spectrum and attenuation can be greatly underestimated when no noise sources from the near field of receivers are detected.
In terms of more local scale, fewer attempts have been done to retrieve the attenuation by using data 1 3 Vol.: (0123456789) collected by seismic arrays. Albarello and Baliva (2009) estimated the attenuation in soil by seismic noise measurements but did not obtain Q S profiles. Parolai (2014), taking advantage of the previous studies, showed that it is possible to reliably estimate Q S in the shallow-most geological layers by using seismic noise recordings from microarrays with an extension of a few tens of meters. Boxberger et al. (2017) confirmed these results considering data sets of seismic noise recorded by microarrays of seismic stations in different geological environments of Europe and Central Asia. Furthermore, they highlighted a slight inverse correlation between V S and Q S , in agreement with Campbell (2009). They proposed that in unconsolidated materials the role of the grain matrix and grain contacts play a major role with respect to V S and Q S . Recently, similar results were obtained by Boore et al. (2020) from analyses of active shear-wave surface source recordings in borehole data. Based on the idea of Boxberger et al. (2017), Nardone et al. (2020) retrieved the mean regional 1D Q S model of Ischia volcanic island (Italy) up to a depth of 2 km. Their results were realistic and comparable with those found at both Stromboli volcano and the volcanic area of Campi Flegrei by other authors. Haendel et al. (2016), adopting the method of Zhang and Yang (2013), estimated the Love wave quality factor (Q L ) at frequencies above 1 Hz from a two-circle array with the largest diameter of 1.8 km. Q L curves were modeled using Q S values calculated from previous laboratory and geophysical measures. They found some differences between observations and theoretical Q L curves that can be ascribed either to local scattering or poor representation of the subsoil below the stations from the considered profiles.
Lastly, Meng et al. (2021) investigated both controlled and uncontrolled vehicle-induced ground motion using two short linear arrays. Their results were consistent with those of Liu et al. (2015) and Q S values from shallow logging data of Fletcher et al. (1990). They noted that in their experiment involving controlled vehicle-induced ground motion, Q values depended in part on the phase velocities adopted. They also highlighted that phase velocities could vary across damage areas (Zigone et al. 2019) giving rise to overestimated Q values at the sensors for a certain frequency, and that their non-homogeneity could perturb the vehicle-produced radiation from the isotropic pattern taken.

Seismic noise recordings in boreholes
Recently, recordings of ambient vibrations in vertical arrays placed inside boreholes were investigated for retrieving near-surface attenuation. Haendel et al. (2019) conducted research to estimate Q from seismic noise by means of DeConvolution Interferometry (DCI), based on previous applications of DCI to earthquake data in boreholes (e.g., Trampert et al. 1993;Parolai et al. 2010Parolai et al. , 2012Fukushima et al. 2016;Raub et al. 2016). Haendel et al. (2019) analyzed seismic noise and small local magnitude (M L ) earthquake recordings in a three-component 87-m-deep borehole and a surface station in the West Bohemia/Vogtland area, at the Czech-German border. They used the interferometric technique by deconvolving seismic motion recorded at depth with the one recorded at the surface for both earthquakes and seismic noise data. For both types of data, they separated the wavefield into up-and down-going waves and obtained Q S from the transfer function of up-and down-going waves in a deconvolved trace, based on the method proposed by Fukushima et al. (2016). Haendel et al. (2019) noted that DCI applied to 1 h of seismic noise recordings led to stable deconvolution results and, consequently, Q S values were determined in the frequency range from 5 to 15 Hz. It was, thus, simpler than considering several earthquakes for obtaining an equal stable deconvolution. Furthermore, with seismic noise, the non-zero incidence of earthquake waves could be excluded. DCI of seismic noise, however, did not allow them to derive Q S above 15 Hz, differently from DCI of earthquakes.
A few years before, Dikmen et al. (2016) focused their study on man-made noise (also called anthropic noise), that is vibrations generated by human activities and structures. This kind of noise, which is recorded by seismic sensors together with that produced by natural sources (e.g., wind and ocean waves), was discussed in previous publications (e.g., O'Connell 2007;Latorre et al. 2014). On this basis, Dikmen et al. (2016) considered traffic-induced seismic noise at the 140-m-deep Ataköy (Turkey) vertical array, to calculate principally κ and, secondly, Q. Seismic noise induced by known sources (i.e., a major highway, a subway line with relative station, a shopping mall, and a sports arena) at a close distance was recorded at the engineering bedrock level. Results of noise data collected in one day were compared by them with the ones obtained by the analysis of a strong motion event, highlighting a difference attributed to near-surface attenuation. They concluded that traffic-induced seismic noise could be useful for estimating subsurface attenuation. However, strong traffic-induced seismic noise level was required, and quiet hours could not provide a sufficient level of noise.
Based on the aforementioned studies, seismic noise analyses in boreholes appear promising for deriving Q S in the near-surface and/or taking advantage of lesser attenuated recordings at depth. Further investigations should be conducted in this respect to enhance what has been achieved so far.

High-frequency S-wave attenuation
In the previous Sections, methodologies to estimate the quality factor for shear waves are presented in detail. In this Section, we briefly describe alternative attenuation characteristics of the Earth crust, with a focus on the high-frequency decay parameter kappa (κ) and its relation to Q. In particular, we highlight how some of the first seminal papers influenced more recent advances. We note here the distinction between this high-frequency decay parameter, for which we will use the notation of the Greek letter κ-with the appropriate subscripts where necessary-and the wavenumber mentioned in previous Sections, for which we have used the notation k 0 or k S . Anderson and Hough (1984) was the seminal study that first defined κ as the slope of the Fourier amplitude spectrum of acceleration of the S-wave window when plotted in a log-linear scale and when the regression is performed in an appropriate frequency band, above the source corner frequency and below the frequencies where noise becomes dominant (e.g., see Fig. 7). Prior to that, the observation had been made of the spectrum deteriorating rapidly above a certain frequency f max , though its interpretation was not straightforward: Hanks (1979) first considered the decay as a path effect, until Papageorgiou (1981) suggested a source effect after showing that the decay was still present after path correction; this led Hanks (1982) to reinterpret f max as a site effect, while later researchers also suggested effects with source origins (Papageorgiou 2003). A comprehensive review on κ as site attenuation and its relation to Q and damping can be found in Campbell (2009). Anderson and Hough (1984) measured κ from S-wave windows recorded at different sites (with varying site conditions), from recordings of several events with varying distances. They found two κ components, a zero-distance one that remained roughly constant per site and insensitive to the different events, and the positive derivative of κ with distance. The former was related to the attenuation of S-waves in their short, vertical propagation beneath the site within the shallower geological layers, and the latter was related to their incremental attenuation along their longer, horizontal propagation from source to site within the crust. Because the quantity originally defined by Anderson and Hough (1984) as κ includes the effect of distance, and in order to disambiguate Fig. 7 (a) Generic shape of the acceleration Fourier amplitude spectrum in log-log space, with plateau between source corner frequency and f max . (b) Definition of κ as the slope of the spectrum in the frequency range from f 1 to f 2 , in loglinear space (modified from Ktenidou et al. 2014) the various notations used in κ studies, the taxonomy of Ktenidou et al. (2014) recommended the use of the notation κ r for the total effective slope of the spectrum, and κ 0 for its site-specific component at zero distance. The derivative of κ r with distance R will be referred to in this work as dκ/dR, so that: where R e is the epicentral distance. (dκ/dR)R e was initially approximated by a straight line.  later confirmed that κ 0 varied with local site stiffness, while dκ/dR was a deeper, regional effect. In the (much later) asymptotic model of Ktenidou et al. (2015), κ 0 for stiff formations is considered to be also bound by regional upper crust and source characteristics, while in the nationwide maps of New Zealand by Van Houtte et al. (2018), hard-rock κ 0 was correlated to crustal attenuation.
The underlying assumption in the κ r model as a straight line trending with frequency is that the path attenuation can be modeled as frequency-independent Q S within the range κ r is computed, i.e., from f 1 to f 2 . In reality, a Q S (f) = Q o f n model, as generally derived from the inversion of seismological data, would incur a curvature upon the linear trend of κ r with frequency, and hence would cause concavity, i.e., a decreasing slope with increasing frequency. This was observed, e.g., by Edwards et al. (2015) and later modeled by Haendel et al. (2020).
The underlying assumption in the dκ/dR model as a straight line with distance is that Q S is constant and frequency-independent across the path, within the crust below the sediments. Below the site, the assumption is that Q depends strongly on depth, increasing from the near-surface sediments to the crustal value. There is a rough expression to compute the overall regional crustal Q S from the slope of κ r with epicentral distance (Eq. 8 of Ktenidou et al. 2015): where β is the crustal V S . Some studies have used this equation to compute an averaged Q across the frequency range from f 1 to f 2 and compare it to regional studies of Q S (f). However, this usually leads to an extrapolation to higher frequencies with respect to the maximum frequency usually considered in Q S (f) models.
(5.1) Anderson and Hough (1984) suggested the linear dκ/dR model as a convenient simplification rather than a necessity. Indeed,  later found non-linear distance dependencies corresponding to multi-layer (rather than constant with depth) crustal Q models. Anderson (1991) extended the model by simply imposing a smooth variation with no other prior assumption. Alternative concave or convex models were later used (e.g., Fernández et al. 2010;Gentili and Franceschina 2011), while more recent studies sometimes found a lack of distance dependence in the first 40-90 km, depending on the region (Kilb et al. 2012;Ktenidou et al. 2017Ktenidou et al. , 2021, an observation that has been referred to as the "hockey stick" effect. One of its advantages is that the long-distance recordings cannot bias the site κ 0 value, which is determined only by short-distance data (such bias on κ 0 had been quantified, e.g., by Edwards et al. 2015). In lieu of this zero-distance extrapolation based on the κ r distribution, an alternative approach has long existed, which avoids assumptions as to dκ/dR and allows for the consideration of the frequency dependence of Q. The alternative is the correction of the spectrum with a prior Q(f) = Q o f n model from literature. This was done, e.g., by Silva and Darragh (1995), yielding equivalent κ 0 values with no need to extrapolate to zero distance. The downside is that Q(f) power laws are not often derived out to very high frequencies and thus may need to be extrapolated themselves. Another alternative, just if very abundant data exist, is to only use recordings from very short distances (e.g., Van Houtte et al. 2014).
Unlike Q and damping ratio, κ r has units of time and is thus comparable to the site attenuation factor t* of Solomon (1973). Fletcher and Boatwright (2020) noted the advantage of t* being how it additive from one layer to the next, and of being related to both the ray path and the site conditions, and we believe the same advantages hold for κ r . However, although κ r is usually related to epicentral distance, t* is typically related to hypocentral distance in order to trace the ray path. Fletcher and Boatwright (2020) give the relation: where, as in Eq. (5.2), β is the crustal V S and R h is the hypocentral distance. Comparing Eqs. (5.2) and (5.3) and ignoring the difference in the distance metric, t* is again seen to be roughly equivalent to κ r .  defined the conditions under which these two parameters are equivalent, namely assuming a ω −2 source spectrum, that a frequency-independent contribution to Q is separable, that the decay is mostly due to the direct wave, and that site resonances are small. The issue of shallow site resonances possibly biasing κ r was first noted by Anderson and Hough (1984). Parolai and Bindi (2004) quantified such bias theoretically in the presence of soil layers, while Ktenidou and Abrahamson (2016) demonstrated that there may be inherent effects even in rock κ 0 values, due to near-surface site-specific broadband amplifications, or due to generic crustal amplification down to several kilometers.
Resonance effects due to shallow or deeper impedance contrasts can bias κ 0 up or down, causing unpredictable trade-offs between site amplification and site attenuation depending on the chosen frequency band f 2 -f 1 . On the other hand, the correlation between κ 0 and the overall site stiffness is more predictable. Anderson et al. (1996) first implied monotonic negative correlations between t* and V S , i.e., lower site attenuation for stiffer sites, after which several empirical correlations between κ 0 and V S30 were proposed, starting with Silva et al. (1998) and including Chandler et al. (2006), Van Houtte et al. (2011), and others. Although most such correlations suggest an indefinite decrease of κ 0 for very hard rock, the conceptual model of Ktenidou et al. (2015) suggested instead a regionally bound asymptotic stabilization for high V S . All correlations are statistically unreliable with large uncertainties, and some of the reasons, as outlined in the review paper of Ktenidou et al. (2014), include these being related to regional differences, methodological issues, and inherent uncertainties in the V S estimates themselves. Among the various reasons behind the large uncertainties in the κ 0 values, one is the inter-analyst biases, i.e., the differences in technique between different analysts. After quantifying the possible effect of such biases, and in order to minimize them, Ktenidou et al. (2013) proposed a step-by-step flowchart procedure and some criteria to homogenize future calculations.
Similar to Q, κ 0 as energy loss can also be related to two physical mechanisms: frequency-dependent scattering and frequency-independent material absorption (Dainty 1981). From the point of view of scattering, and particularly for high frequencies, even before the effects of κ 0 , stratigraphic filtering was well known (i.e., that wave propagation through fine layering can filter out high frequencies and may increase the apparent attenuation through short-period multiples; O'Doherty and Anstey 1971; Richards and Menke 1983). The notion of scattering κ 0 attenuation was more recently examined theoretically by Parolai et al. (2015) and Pilz and Fäh (2017), and also empirically by Ktenidou et al. (2015) and Parolai (2018). From the point of view of absorption,  proposed a seismological framework in which the frequency-independent part of κ related to intrinsic energy loss depended on temperature and pressure conditions in the upper crust. From a geotechnical engineering point of view, Silva (1997) related Q in the surface layers with the soil damping ratio ξ as: e.g., a soil with 5% damping corresponds to a Q of 10 (this likely holds for a frequency range between around 0.1-10 Hz). Fernandez-Heredia et al. (2012) first related κ 0 to the soil damping for surface recordings. After Van Houtte et al. (2011) opened the way to computing κ 0 for downhole arrays, Ktenidou et al. (2015) computed the difference in κ 0 between the top and bottom sensor of a vertical array (κ 0_surface − κ 0_bedrock ) as an indicator of the total effective attenuation-due to both intrinsic absorption and frequency-dependent scattering-in the soil column between the stations. Recently, Xu et al. (2020) did so for an extended downhole dataset. For strong ground shaking, an unavoidable aspect of damping is its increase due to nonlinearity. Dimitriu et al. (2001) first related κ 0 to non-linearity, finding a positive correlation between PGA and κ 0 , as stronger shaking increased the absorption κ 0 component. Van Houtte et al. (2014) (perhaps surprisingly) found the opposite, i.e., κ 0 decreasing with peak ground acceleration; this was attributed either to the stiffness degradation and non-linearity smoothing out the small-scale heterogeneities contributing to scattering κ 0 component, or to the fast-shearing rate decreasing the material damping.
It is clear from all of the above that distinguishing the various components of κ 0 and resolving its trade-offs with distance and site amplification can be challenging (note that trade-off issues related to the seismic source were deliberately left out of this discussion as peripheral to its topic, even though the uncertainty in stress drop and thus corner frequency is a key issue for correct κ estimates, e.g., Ktenidou et al. 2017). One of the key factors that limits data availability and complicates resolution of all such trade-offs is the level of noise in the recording. This was recognized ever since Anderson and Hough (1984) examined digitization noise levels for old analog data. In modern data, site noise can impede higher-frequency calculations for both κ and Q. One way forward was recently suggested by Pikoulis et al. (2020), utilizing stochastic noise modeling to correct the signal spectra and render them usable out to higher frequencies than was possible up to now.

Summary and discussion
In this manuscript, we reviewed the state-of-the-art methods for estimating the shear-wave quality factor (Q S ) in the near-surface geological layers. After an overview of the various definitions of the quality factor Q, we compiled diverse methodologies used by different disciplines (applied geophysics, geotechnical engineering, and seismology) under a single unified framework, using a structure that allows readers interested only in a particular approach to find all the necessary information in a single Chapter/paragraph. This review paper illustrates the wide range of efforts that have been conducted, providing reliable Q S estimates. This is indeed an extremely challenging problem because of the variety of interacting phenomena that affect the spatial and temporal attenuation of a mechanical wave propagating in multicomponent systems, such as geomaterials, which are inherently inhomogeneous and anisotropic at multiscales. If Q is considered a constitutive parameter only when it is associated with the energy dissipated by a unit volume of soil or rock undergoing cyclic deformations and thus, only related to the inelasticity of the material, then in situ experimental measurements are difficult because wave attenuation may also be influenced by scattering and geometric spreading. Laboratory tests, on the other hand, do not suffer from these phenomena; however, they are performed on small soil samples which are not representative of the large volume of soil at a site of interest and, in addition, may also be affected by undesirable effects such as sample disturbance and equipment-generated damping. For these reasons, there is currently no single method that can be considered superior to the others, as each of them has its own advantages and shortcomings.
Our work highlights the need for further harmonization and consistency in the terminology used to achieve better communication across different disciplines. We hope this article will contribute towards facilitating this process and help shedding light onto the world of attenuation, where the numerous different definitions of what is often the same phenomenon are sometimes confusing, especially for practitioners. In Table 2, we summarize the methodologies discussed throughout the manuscript, together with their main pros and cons.
This review paper also indicates that while the effect of the apparent attenuation can be reliably estimated, the distinction between intrinsic and transmission (i.e., reflection, refraction, and scattering) effects still requires dedicated efforts. This becomes an issue of major importance when empirically estimated Q values are used in numerical simulations of ground motion or site response. In fact, substituting the intrinsic quality factor with an apparent one (which includes all components, hence leading to double counting) might strongly affect the amplitude and duration of the obtained time series, with obvious impact on downstream calculations. This overview also highlights the wide range of applications of the Q S factor, as estimated in different disciplines. We recommend caution when using Q values estimated with techniques intended for specific purposes that were not intended for estimating Q. Furthermore, the frequency range used for the Q estimation and the strain level exerted on the material-both in the laboratory and in the field-should always be considered.
Although this review may inadvertently overlook studies contributing to near-surface Q S estimation, we hope it may foster further developments which may in turn fill the existing research gaps that we suggest herein.    Li and Liu (2015) and Xue et al. (2020) Seismic refraction In situ Active Free-surface Method not yet well-established for estimating Q Bibliographic references are, e.g., Socco et al. (2005) and Kurtuluş and Sertçelik (2010) Earthquake Acknowledgements The authors wish to thank Frank Rackwitz (Technische Universität Berlin) for his help and suggestions, and Carla Barnaba (National Institute of Oceanography and Applied Geophysics -OGS) for her useful comments. Kevin Fleming kindly revised the English. Thanks to the Editor and the two anonymous reviewers whose comments and valuable suggestions improved the quality of the manuscript. Figure 1 was produced using the Generic Mapping Tools (GMT; Wessel and Smith 1998) and the open-source vector graphics editor Inkscape. Figure 2 was created using PowerPoint. Figures 3 and 4 were realized using MATLAB with the assistance of Ali Guney Özcebe (Eucentre, Pavia). Figure 5 was produced with the open-source vector graphics editor Inkscape, whereas Figure 6 was created with the Python's libraries Matplotlib (Hunter 2007) and NumPy (Harris et al. 2020). Figure 7 was produced using the Generic Mapping Tools (GMT; Wessel and Smith 1998).

Funding
The Consortium of Organizations for Strong Motion Observation Systems (COSMOS) and a group of North American power utility companies, consisting of Southern California Edison and Pacific Gas and Electric, identified the need for these guidelines for best-practices and provided funding and encouragement to facilitate the project. This material is also based upon work supported by the US Geological Survey under Cooperative Agreement No. G17AC00058. The views and conclusions contained in this document are those of the authors and should not be interpreted as representing the opinions or policies of the US Geological Survey. Mention of trade names or commercial products does not constitute their endorsement by the US Geological Survey. Open access funding were provided by the National Institute of Oceanography and Applied Geophysics-OGS within the CRUI-CARE agreement.

Data availability Not applicable.
Code availability Not applicable.

Conflict of interest
The authors declare no competing interests.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.