A Procedure for Developing Uncertainty-Consistent Vs Profiles from Inversion of Surface Wave Dispersion Data
AA Procedure for Developing Uncertainty-Consistent Vs Profilesfrom Inversion of Surface Wave Dispersion Data
Joseph P. Vantassel and Brady R. Cox The University of Texas at Austin
July 21, 2020
Abstract
Non-invasive surface wave methods have become a popular alternative to traditional invasive forms ofsite-characterization for inferring a site’s subsurface shear-wave velocity (Vs) structure. The advantageof surface wave methods over traditional forms of site characterization is that measurements made solelyat the ground surface can be used routinely and economically to infer the subsurface structure of asite to depths of engineering interest (20-50m), and much greater depths ( > a r X i v : . [ phy s i c s . g e o - ph ] J u l Introduction
Surface wave methods have become a popular alternative to traditional invasive forms of site characteriza-tion for measuring a site’s subsurface shear wave velocity (Vs) structure. Vs, which is directly related tothe site’s small-strain shear stiffness ( G max = ρV s , where ρ refers to mass density), is a critical parameterin many seismic hazard analyses, such as ground response analyses (Foti et al., 2009; Rathje et al., 2010;Teague and Cox, 2016; Passeri et al., 2019) and liquefaction-triggering (Andrus and Stokoe, 2000; Kayenet al., 2013; Wood et al., 2017). The performance of surface wave testing is traditionally broken into threestages: acquisition, processing, and inversion (Foti et al., 2015). The acquisition stage involves non-invasively(i.e., from the ground surface) measuring surface waves as they propagate through a site. This can be ac-complished using active-source methods, where surface waves are generated by the experimenters, and/orpassive-wavefield methods, where sensors are left undisturbed to record ambient surface waves. The process-ing stage transforms these time-domain recordings to measurements of the site’s dispersive properties, whichin this context describes how the site’s surface wave phase velocity changes as a function of frequency (orequivalently wavelength). Importantly, this measurement of the site’s dispersive properties should includesite-specific estimates of frequency-dependent uncertainty, discussed in detail later. For clarity of expression,we will refer to the measurement of the site’s dispersive properties from the processing stage as the experi-mental dispersion data. Importantly, when surface wave processing is performed by an experienced analyst,the experimental dispersion data is robustly determined, with uncertainties that typically range from 5%– 10% coefficient of variation (COV) between analysts (Cox et al., 2014; Garofalo et al., 2016). This con-sistency led Griffiths et al. (2016a; 2016b) and Teague et al. (2018a) to aptly refer to the experimentaldispersion data as part of a site’s “signature”. The final stage of surface wave testing is inversion. In thisstage, numerical search algorithms are used to identify one-dimensional (1D) layered earth models whosetheoretical dispersion curves, determined through an analytical forward problem, best fit the experimentaldispersion data. It is important for the reader to understand that, while throughout this paper we focusprimarily on Vs, as it has the greatest impact on the theoretical dispersion curve (Wathelet, 2005) and is ofprimary importance to subsequent engineering analyses, the computation of a theoretical dispersion curve(and therefore surface wave inversion) requires the definition of an entire 1D ground model. These groundmodels are composed of a stack of layers described by their thickness (H), Vs, compression-wave velocity(Vp), and ρ . The models whose theoretical dispersion curves best fit the experimental dispersion data, asdetermined by a misfit function, are considered to be the most likely representations of the site’s subsurfaceconditions.As the primary focus of this paper is on the consistent propagation of experimental dispersion datauncertainty into the Vs profiles resulting from inversion, it is important to first discuss the types and potentialsources of these uncertainties. Uncertainties can be broadly grouped into two categories commonly used inprobabilistic seismic hazard analysis (PSHA): epistemic and aleatory. Epistemic uncertainty describes thoseunknowns which stem from a lack of knowledge. In surface wave testing, sources of epistemic uncertaintyinclude, among others, the selected dispersion processing wavefield transformation and the number of 1Dsubsurface layers used during inversion. Aleatory uncertainty, or sometimes referred to as aleatory variability,is the result of inherit randomness within the quantity being measured. When discussing site characterizationin particular, aleatory uncertainty is classically used to define how the material properties in the subsurfacechange in three-dimensional space (i.e., randomness with space, or spatial variability). Of course, this couldalso be considered epistemic uncertainty because it stems from a lack of knowledge rather than from anyinherit randomness within the subsurface (i.e., randomness with time). Being able to consider subsurfacevariability as either epistemic or aleatory in nature hints that, for most real-world applications, the twouncertainty “categories” are not as distinct as their names may indicate. Regardless, it is to be expectedthat all experimental dispersion data contains some amount of uncertainty and that this uncertainty is in partepistemic (dependent largely on the quantity and quality of the data acquired) and part aleatory (dependenton the complexity of the subsurface and the spatial extents of the surface wave arrays). And, while it isdifficult-to-impossible to separate epistemic and aleatory uncertainties in the experimental dispersion data,it is of paramount concern that the surface wave analyst attempt to quantify the combined uncertainty on asite-by-site basis (discussed briefly next) and then propagate that uncertainty through the inversion process2n order to develop uncertainty-consistent Vs profiles.As this study begins under the assumption that site-specific estimates of dispersion uncertainty havebeen established, how this might be accomplished deserves a brief discussion. However, before doing so, it isimportant to acknowledge that while many previous researchers have developed procedures for accountingfor uncertainty in experimental data (Lai et al., 2005; Foti et al., 2009; Cox and Wood, 2011; Teague et al.,2018a) no single procedure has been accepted widely into practice, and therefore the example presentedhere is but one of many potential alternatives that may be used. To facilitate a more practical discussion,the example presented here is shown in Figure 1. Figure 1a illustrates a measurement of a site’s dispersiveproperties at some location termed “Location A”. This measurement is directly linked to a given location (i.e.,Location A), experimental setup (e.g., array configuration), wavefield recording (e.g., source offset or noisetime window), and wavefield processing method (e.g., frequency-wavenumber transformation). The light-colored portions of Figure 1a show high surface wave power and the dark-colored portions low surface wavepower. The maximum power at each frequency is selected as the representative experimental dispersion datafor that location, setup, wavefield recording, and processing method. To then develop meaningful statisticsthat incorporate the uncertainties previously discussed, the procedure shown in Figure 1a should be modifiedand/or repeated in a systematic manner to encompass reasonable combinations of these contributing sources.As shown schematically in Figure 1b, these combinations may include performing the test at various locationsacross a site (i.e., Location B), using various source locations for active-source experiments (i.e., Offset X),various time-windows for passive-wavefield experiments (i.e., Window Y), and multiple wavefield processingmethods (i.e., Method Z). Finally, once these sources of uncertainty have been accounted for in the formof experimental dispersion data (i.e., Figure 1b) they may be summarized into a statistical representationat each frequency (or wavelength), as shown in Figure 1c. A previous study by Lai et al. (2005) showedthat the uncertainty in experimental dispersion data for active-source experiments based on multiple sourceimpacts was normally distributed, therefore, it is common to represent the Rayleigh wave velocity (Vr) ateach frequency with a mean and standard deviation. However, it is important to note that at some sites withsignificant lateral variability, such as the Garner Valley site examined by Teague et al. (2018b), it may bemore appropriate to develop alternative dispersion data sets (one per location) with their own experimentaluncertainty and invert them separately rather than trying to represent all of the dispersion data with asingle statistical distribution. As a more in depth discussion of developing experimental dispersion data isbeyond the scope of this paper, the reader is encouraged to carefully consider the impact of how definingtheir experimental dispersion uncertainty will affect the resulting distribution of Vs profiles and whetherthose profiles will correctly communicate the measured uncertainty.While not the primary focus of this work, it is important to briefly mention some of the methods currentlybeing used to account for Vs uncertainty in seismic hazard studies and highlight why developing uncertainty-consistent Vs profiles from experimental dispersion data is of such critical importance. At present, mostimportant seismic design projects utilize some form of Vs randomization to account for Vs uncertainty inground response analyses, where a baseline Vs profile is randomized (i.e., manipulated) using a statisticalprocedure. The most commonly used of these randomization approaches is that proposed by Toro (1995),leading to its adoption in the guidelines for hazard-consistent one-dimensional ground response analyses(Stewart et al., 2014) and the design of critical facilities (EPRI, 2012). However, multiple studies (Griffithset al., 2016a,b; Teague et al., 2018b; Teague and Cox, 2016) have expressed serious concern that the blindapplication of this type of Vs randomization can result in highly unrealistic ground models that do not fitthe experimental site signature and an inability of such models to accurately predict ground response. Thereader will note that a new Vs randomization approach has recently been proposed by Passeri et al. (2020).However, the authors have not yet evaluated it rigorously to determine if it produces Vs profiles that aremore consistent with the site signature. Importantly, this paper, while not discussing Vs randomizationdirectly, offers an alternative to Vs randomization by propagating measured, site-specific uncertainty intosuites of Vs profiles obtained from surface wave testing. These suites of profiles can then be used eitherdirectly to address Vs uncertainty in subsequent engineering analyses (e.g., ground response analyses) or,if the use of Vs randomization is still desired, to better inform its many unknown input parameters withsite-specific and uncertainty-consistent values. While this study will focus primarily on the development of3 V r ( m / s ) (a) Location A 3 10 (b) Location ALocation BOffset XWindow YMethod Z 3 10 (c) Exp. Disp. Data
Figure 1: Schematic illustrating one possible procedure for developing experimental dispersion data withmeasures of uncertainty in terms of the site’s Rayleigh wave velocity (Vr). This involves: (a) the processingof recorded waveforms from a single location, experimental array setup, source offset/noise time window, andwavefield transformation method, (b) the modification and repetition of the procedure illustrated in (a) toproduce estimates of the site’s dispersion uncertainty, and (c) the synthesis of these estimates of uncertaintyinto a statistical representation called the site’s experimental dispersion data.the former, it is important that the reader be aware of the wider implications of the proposed method ongreater engineering practice.This study begins by examining common approaches found in the literature for developing uncertainty-consistent Vs profiles from the inversion of experimental surface wave dispersion data. Using a syntheticdataset, these approaches from the literature are shown to yield suites of subsurface models whose theoreticaldispersion curves severely underestimate the experimental dispersion data’s uncertainty, resulting in Vsprofile which underestimate the site’s uncertainty. This is then followed by the presentation of a newprocedure for rigorously propagating measured experimental dispersion uncertainty through the inversionprocess to obtain suites of Vs profiles that more accurately represent Vs uncertainty at the site. Thisnew procedure is applied to the same synthetic dataset as the literature-based approaches and is shownto be quantitatively superior, and able to precisely propagate the experimental dispersion uncertainty intothe resulting Vs profiles. The new procedure is then extended to incorporate epistemic uncertainty in theinversion’s layering parameterization using the same synthetic example, and again shows excellent results.The study concludes with the application of the method at a real site where it is shown to produce suites ofVs profiles which agree favorably with a borehole Vs profile while simultaneously capturing the experimentaldispersion data’s uncertainty.
To illustrate the problem this paper proposes to solve, we first examine several approaches from the literaturethat have been used to account for Vs uncertainty in surface wave inversion. All of the approaches mentionedbelow have a basic commonality, in that they first search through large numbers (often tens-of-thousands tomillions) of trial layered-earth models to find a significant number of models with an acceptable fit to theexperimental dispersion data. The acceptability of a model is typically judged using a misfit function thatquantifies the goodness-of-fit between the theoretical dispersion curve for a given model and the experimentaldispersion data. Due to the non-uniqueness of the inverse problem and the experimental dispersion data’suncertainty, it is often possible to find thousands- to tens-of-thousands of acceptable models that rangefrom subtly to significantly different. From these acceptable models, a subset are selected to account for4ncertainty in Vs. The discussion below will focus primarily on how different studies obtained suites ofacceptable models, and how they then selected from those to estimate Vs uncertainty.Before presenting the various approaches that have been used to account for Vs uncertainty in surfacewave inversion, it is important to briefly discuss the details surrounding surface wave inversion, as thesedetails will be important for understanding the results presented later. First, all of the approaches presentedin the literature utilize large numbers of trial models from a global-search inversion algorithm to find suitesof acceptable models. Global-search algorithms vary in their implementation and may search for accept-able models in several different ways, including randomly (i.e., pure Monte-Carlo), with the aid of someoptimization algorithm (i.e., pure optimization), or by using a combination of the two. By far, the mostpopular tool for performing global-search surface wave inversion is the Dinver module (Wathelet et al., 2004)of the open-source software Geopsy (Wathelet et al., 2020). As a testament to its popularity, Dinver wasthe inversion algorithm of choice in all but one of the works discussed below. Of particular relevance to thediscussion below, and worth discussing here, is the misfit function proposed by Wathelet et al. (2004) andimplemented in Dinver. The Wathelet et al. (2004) misfit function can be described as a root-mean-squareerror normalized by the experimental dispersion uncertainty (Yust et al., 2018). This gives the misfit functiona useful physical interpretation, as it represents on average, across all frequencies/wavelengths, how far (innumber of standard deviations) a theoretical dispersion curve strays from the experimental dispersion data.For example, a theoretical dispersion curve with a misfit of 1.0 can be understood as a curve that on averageis one standard deviation away from the mean. With these details in mind, it is now possible to discuss theapproaches presented in the literature to account for Vs uncertainty in surface wave inversion.To develop suites of acceptable layered-earth models, Wathelet et al. (2004) and Wathelet (2008) per-formed pure-optimization inversions considering many (over 100,000) trial models. For convenience of ex-pression, and in order to be consistent with common vernacular, we will refer to a pure optimization inversionas a minimum misfit of zero (M0) inversion. M0 inversions are named as such because the goal of the in-version algorithm is to find a model with a misfit equal to zero (i.e., one whose theoretical dispersion curveperfectly matches the mean experimental dispersion data). In the studies noted above, from the many trialmodels searched, all models with a dispersion misfit value less than 1.0 (tens of thousands in most cases)were selected as a means to propagate the experimental dispersion data’s uncertainty into the resulting Vsprofiles. This approach has three shortcomings: (1) it may result in a highly variable number of acceptableprofiles, depending on the number of trial models attempted in the inversion and the subjective user-definedquality threshold (i.e., misfit less than 1.0 in this case), (2) for inversions with many models below thethreshold (e.g., tens-of-thousands) it may readily become computationally unmanageable to again propagatethe Vs uncertainty into subsequent engineering analyses (e.g., ground response analyses), and (3) the use ofan M0 inversion will likely cause the majority of the theoretical dispersion curves to be clustered around themean of the experimental dispersion data (i.e., a misfit of zero) rather than following the distribution of theexperimental dispersion data uncertainty.Foti et al. (2009) used a pure Monte-Carlo global-inversion algorithm (i.e., not Dinver) to develop suitesof acceptable models (more than 50,000 in the examples presented) whose theoretical dispersion curves fit theexperimental dispersion data. To select a more manageable subset of models, they sampled from these largesuites using a statistical test to select models which could be considered equivalent in term of their fit to theexperimental dispersion data, given its uncertainties. While the selection process reduced the computationalburden of having to consider ten-of-thousands of Vs profiles, it still maintained the disadvantage of theprevious technique that the number of profiles varied between application, ranging between 6 and 270 Vsprofiles for the cases presented. Furthermore, the theoretical dispersion curves from the selected modelscan be seen to underestimate the experimental dispersion data uncertainty, even for the most favorable casewhere 270 profiles were selected by the statistical test.Hollender et al. (2018), when developing statistics for the time-averaged shear-wave velocity in the upper30m (Vs30), used an “acceptable-misfit” approach [i.e., an inversion where the minimum misfit was not zero(M0+)] to develop a large suite of acceptable models (over 50,000) to represent the experimental dispersiondata uncertainty. From those models, they randomly extracted at least 6000 Vs profiles to develop statisticson Vs30. They then selected a smaller set of 33 representative Vs profiles, whose uncertainty in Vs30 matched5hat of the randomly selected 6000 profiles, for use in subsequent analyses. While this approach ensuredan equal and manageable number of profiles (i.e., 33) for each site, and permitted the propagation of theestimated uncertainty in Vs30, it does not guarantee that the initial large suite of acceptable models properlyaccounted for the experimental dispersion data’s uncertainty. Thereby, potentially allowing the subsequentuncertainties (i.e., the Vs30 uncertainty) upon which they are based to not fully capture the site-specificuncertainty in the experimental dispersion data.This final section groups a number of studies that are similar, in that they all use a M0 inversion todevelop their acceptable suites of models from which to select some fixed number of lowest misfit modelsto represent Vs uncertainty. Di Giulio et al. (2016) selected the 100 lowest misfit profiles out of over amillion trial models. Griffiths et al. (2016a; 2016b) and Teague and Cox (2016) chose 50 randomly selectedmodels from the 1000 lowest misfit models obtained from hundreds-of-thousands of trial models, as the50 randomly selected models were shown to have similar statistical properties as the 1000 best, but weremore computationally manageable. Cox and Teague (2016), Teague et al. (2018a) and Deschenes et al.(2018) presented the 1000 lowest misfit Vs profiles from hundreds-of-thousands of trial models. Vantasselet al. (2018), Cox and Vantassel (2018), and Yust et al. (2018) presented the 100 lowest misfit Vs profilesfrom hundreds-of-thousands of trial models considered to account for Vs uncertainty, as these 100 lowestmisfit Vs profiles were shown to be statistically similar to the 1000 lowest misfit Vs profiles, but weremore manageable. This approach (i.e., selecting some number of the lowest misfit/“best” models) resolvesthe issues of variability between analyses, avoids potentially large numbers of profiles, and is relativelysimple. However, this approach is still deficient in two specific ways: (1) it relies on an M0 inversion which,as mentioned previously, will likely produce profiles with theoretical dispersion that is clustered aroundthe experimental dispersion data’s mean rather than being properly distributed, and (2) the uncertaintyaccounted for in the 100 or 1000 best models is indirectly tied to the number of trial models attempted bythe analyst (i.e., different apparent uncertainty will results whether the analyst uses 1000, 10,000 or 1,000,000trial models).In summary, all of the presented approaches follow a basic two-step process. First, a large numberof acceptable models are obtained using a global-search algorithm, and second, some subset of profilesare selected from the acceptable models as a basis for assessing Vs uncertainty. In all of the presentedapproaches we observe a dependence on two basic assumptions: (1) that the large suite of acceptable models(i.e., from the first step) properly accounts for the experimental dispersion data’s uncertainty, and (2) thatthe selection process will guarantee a set of Vs profiles which rigorously propagates the uncertainty intosubsequent analyses. The remainder of this section is devoted to quantitatively assessing the veracity ofthese assumptions by examining the effectiveness of these approaches.To quantitatively assess the ability of the approaches from the literature to propagate experimentaldispersion data uncertainty into the resulting Vs profiles, four variations were applied to an experimentaldispersion dataset. The experimental dispersion data is taken from a large synthetic study focused on theperformance of surface wave inversion (Vantassel and Cox, 2020b). The data itself has been published as oneof twelve surface wave inversion benchmarks, which is publically available on the DesignSafe-CI (Vantasseland Cox, 2020a). While a full detailed discussion is provided in the previous references, in short, the syntheticdata was developed by taking the theoretical dispersion curve from an assumed ground model, resamplingit in log-wavelength, and assuming a normal distribution in Rayleigh wave velocity (Vr) with a coefficient ofvariation (COV) of 0.05. A COV of 0.05 was based on typical experimental dispersion data uncertainty valuesfrom several blind analyst studies (Cox et al., 2014; Garofalo et al., 2016). For reference, the experimentaldispersion data is the same as that shown in Figure 1c. The proposed methods to be considered were selectedto address two primary questions of interest, namely: (1) what effect, if any, does the type of inversion (pureoptimization, pure random, or some combination) have on the resulting acceptable models, and (2) howshould the analyst sample from the acceptable models. With regard to the type of inversion preformed, weconsider two inversion alternatives: (a) a pure-optimization (i.e., M0) inversion, as this approach is the mostcommonly used in the literature, and (b) a combined approach that is partly optimized and partly random.This second approach will be referred to as a minimum misfit of 1.0 (M1) inversion because the inversionalgorithm is forced to search randomly in those regions of the model space where the calculated misfit is6elow 1.0. A pure-random search was not explicitly considered here, as it is expected to produce similarresults to the M1 inversion at greater computational expense. With regard to how the representative Vsprofiles should be sampled from the acceptable inversion models, we consider two alternatives: (a) the 100lowest misfit/“best” models (b100), as this was the most popular approach in the literature, and (b) 100randomly selected models from all models with a misfit less than 1.0 (n100). Note that the n100 alternativewas a necessary adaptation to the approach of selecting all profiles below the misfit threshold of 1.0 to ensurea computationally manageable number of profiles and a fair comparison to the b100 profiles.For illustration purposes and to simplify initial discussions, the experimental dispersion data was invertedusing only a single, three-layer parameterization (i.e., a Layering by Number (LN) = 3 parameterization),as this layering parameterization was shown to perform the best when inverting the example dataset inthe previous study by Vantassel and Cox (2020b). The experimental dispersion data were inverted usingthe Neighborhood Algorithm (Sambridge, 1999) as implemented in the Dinver module of Geopsy (Watheletet al., 2020). For both the M0- and M1-style inversions, an initial 10,000 random trial models followedby 50,000 neighborhood-algorithm trial models were considered. Note that the number of neighborhood-algorithm models used here is less than those used in previous studies from the literature, which tended touse a hundred-thousand or more. We do this to mitigate bias in the resulting suite of acceptable modelscaused by the inclusion of many (potentially tens-of-thousands) very similar models with misfits close tozero. We believe the use of a smaller, but sufficient (Vantassel and Cox, 2020b), number of models providesthese techniques the best possible chance of successfully propagating the experimental dispersion data’suncertainty into the resulting models. The selected b100 and n100 models are shown in terms of their Vsprofiles and associated uncertainty in Figure 2. Note, to better illustrate their concentration, the Vs profileshave been discretized in terms of depth and Vs, binned into cells, and color mapped in terms of the numberof profiles in each cell. Vs uncertainty is expressed in the form of the lognormal standard deviation ofVs ( σ ln,V s ), which is commonly used in seismic hazard studies, and is similar to the COV for values lessthan approximately 0.3. The theoretical dispersion curves for the b100 and n100 models are presented withthe experimental dispersion data in Figure 3 so the reader can view how the apparent uncertainty in theinversion-derived Vs profiles relates to the apparent uncertainty in the theoretical dispersion curves.Figure 2 illustrates that the resulting Vs profiles are most sensitive to how they are sampled (b100 vsn100), and significantly less dependent on how the inversion is performed (M0 vs M1). This is shown bothqualitatively when viewing the range of Vs profiles [i.e., panels (a), (c), (e), and (g)] and quantitativelywhen viewing σ ln,V s [i.e., panels (b), (d), (f), and (h)]. For example, the Vs profiles obtained from the b100approach show minimal scatter/uncertainty, while those obtained from the n100 approach show significantlymore scatter/uncertainty. Note that the spikes in σ ln,V s are due to the uncertainty in the layer boundaries,reflecting a shortcoming in how σ ln,V s has been traditionally calculated, and are not the result of uncertaintyin Vs directly. The trends in Vs scatter/uncertainty shown in Figure 2 are tied directly to visible trends inthe scatter/uncertainty of the theoretical dispersion curves shown in Figure 3. The color scale in Figure 3indicates the range of dispersion misfit values, which varies widely between approaches. The b100 models areprimarily composed of models with low misfit dispersion misfit models ( < > four common approaches from the literature for accounting for Vs uncertainty in inversion. In order tobetter illustrate their concentration, the profiles have been discretized in terms of depth and Vs, binned intocells, and color mapped in terms of the number of profiles in each cell. The presented approaches include:(a) a minimum misfit of 0 (M0) inversion with the selection of the 100 lowest misfit/“best” models (b100),(c) a M0 inversion with 100 random profiles selected from all models with a misfit less than 1.0 (n100), (e)a minimum misfit of 1.0 (M1) inversion with the selection of the 100b models, and (g) an M1 inversion withthe selection of the n100 models. The Vs profiles from these approaches are shown alongside the solutionprofile in red. Adjacent to each suite of profiles is their corresponding lognormal standard deviation of Vs( σ ln,V s ). The sharp spikes in σ ln,V s are due to the uncertainty in the profile’s layer boundaries, reflecting ashortcoming in how σ ln,V s has been calculated historically, and are not the result of actual uncertainty inVs directly.at 4 Hz). Figure 5b illustrates that using b100 sampling, which was shown to fit the mean Vr so preciselyin Figure 5a, severely underestimates the true dispersion uncertainty at all frequencies. Recall that thedispersion uncertainty for this synthetic example was set using a COV of 0.05 at all frequencies; meaningthat the b100 sampling with a residual COV of nearly -0.05 underestimates the true dispersion uncertaintyby nearly 100%. In contrast, the n100 models that did not well-fit the mean trend tend to provide a betterrepresentation of the experimental dispersion uncertainty, although in this case they still underestimate themeasured uncertainty at both high and low frequencies by approximately 0.02, or 40%. These results clearlyindicate that the currently available methods of accounting for Vs uncertainty in surface wave inversionare unable to propagate the uncertainty in the experimental dispersion data into the resulting suites of Vsprofiles. 8 M0b100 (a)
M0n100 (b) M1b100 (c) M1n100 (d)
MeasuredSlice @ 4Hz0.00.20.40.60.81.0 D i s p e r s i o n M i s f i t Frequency, f (Hz) R a y l e i g h W a v e V e l o c i t y , V r ( m / s ) Figure 3: Comparison of the experimental dispersion data from a synthetic dataset with the individualtheoretical dispersion curves obtained from inversion using four common approaches from the litera-ture for accounting for Vs uncertainty. Each of the four panels (a)-(d) illustrate one of the four commonapproaches from the literature. They include: (a) a minimum misfit of 0 (M0) inversion with the selectionof the 100 lowest misfit/“best” models (b100), (b) a M0 inversion with 100 random profiles selected from allmodels with a misfit less than 1.0 (n100), (c) a minimum misfit of 1.0 (M1) inversion with the selection ofthe b100 models, and (d) an M1 inversion with the selection of the n100 models. The vertical dashed blueline at approximately 4 Hz, shown on all panels, indicates the location where a “slice” is shown in Figure 4for each of the four approaches.
M0b100 (a)
MeasuredInverted
M0n100 (b)
MeasuredInverted
200 250 300
M1b100 (c)
MeasuredInverted
200 250 300
M1n100 (d)
MeasuredInverted
Rayleigh Wave Velocity, Vr (m/s) P r o b a b ili t y D e n s i t y Figure 4: Distributions of the measured (i.e., experimental) and inverted (i.e., theoretical) dispersion data ata frequency of approximately 4 Hz using four common approaches from the literature to account forVs uncertainty in surface wave inversion. The four panels have the same caption as those shown in Figure 3.9 μ r e s / μ m e a s ( % ) Overestimate Mean VrUnderestimate Mean Vr (a) δ i n v − δ m e a s Overestimate Vr UncertaintyUnderestimate Vr Uncertainty (b)
M0, b100M0, n100M1, b100M1, n100
Figure 5: Quantitative assessment of four common approaches from the literature for propagatinguncertainty in experimental dispersion data to the inverted Vs profiles. The comparison is made on the basisof: (a) the dispersion residual mean ( µ res ) [i.e., the difference between the mean of the inverted theoreticaldispersion curves ( µ inv ) and the mean of the measured experimental dispersion data ( µ meas )] normalized by µ meas , expressed in percent, and (b) the dispersion residual coefficient of variation [i.e., the difference betweenthe coefficient of variation of the inverted theoretical dispersion curves ( δ inv ) and the coefficient of variationof the measured experimental dispersion data ( δ meas )]. The vertical dashed blue line at approximately 4Hz in panels (a) and (b) indicates the location where a slice was shown in Figure 4 for each of the fourapproaches. The four example approaches from the literature considered in the previous section had three specific short-comings: (1) the profiles tended to capture either the mean trend or the variance of the dispersion data,but not both (refer to Figure 5), (2) theoretical dispersion curves from the suites of inverted models didnot the follow the distribution of the synthetic experimental dispersion data (refer to Figure 4), and (3)when the results were found to be unacceptable, as was the case for the synthetic example discussed above,there was no clear procedure to remedy the inconsistency. This section presents an alternative procedureto account for Vs uncertainty in surface wave inversion. The new procedure is presented schematically inFigure 6. Figure 6a shows experimental dispersion data with site-specific measurements of uncertainty. Thisuncertainty is quantified in terms of surface wave phase velocity at each frequency using a mean, standarddeviation, and, very importantly, correlation coefficients between all frequency pairs. All of these statisticscan be easily quantified by the analyst using the approach to developing dispersion data with site-specificuncertainty discussed previously (recall Figure 1). Consider the experimental dispersion data in Figure 6a,which contains 13 frequency points. Its statistics are completely described using 13 mean values (shown witha circle), 13 standard deviations (shown with error bars), and a matrix of 13x13 correlation coefficients (notshown) relating the phase velocity at each frequency to the other 12 values. With the experimental dispersiondata described in terms of its statistics, it’s now possible to simulate a realization from that experimentaldispersion data, shown schematically in Figure 6b. It is during this simulation stage that the inclusion ofthe correlation information is so critical. Without such information one is left to assume the correlations, orworse, independence. However, doing so would likely result in simulated dispersion curves that were erraticfrom frequency-to-frequency and inconsistent with the dispersion data that was used to estimate the uncer-10 u r f a c e W a v e V e l o c i t y Disp. RealizationFrequency S u r f a c e W a v e V e l o c i t y Exp. Disp. DataImplied Disp. DataFits to Realizations (a) (b)(d)
Simulate Realization of Disp. DataRepeat (b) and (d) for N realizationsand Calcuate Implied Disp. Data I n v e r t R ea li z a t i on FrequencyDisp. Realization (c)
Fit to Realization C o m pa r e I m p li ed and M ea s u r ed Figure 6: A new procedure for developing uncertainty-consistent shear wave velocity (Vs) profiles fromsurface wave dispersion data. The procedure involves: (a) describing the site’s experimental dispersion datain terms of its statistics, (b) generating a realization of the experimental dispersion data, (d) inverting therealization to obtain a single, best theoretical fit, repeating (b) and (d) for N realizations, and (c) comparingthe statistics of the resulting theoretical dispersion curves (i.e., the implied dispersion data) to that of themeasured experimental dispersion data.tainty. In essence, the inclusion of the correlations between dispersion points imposes a relational constraintthat helps to encourage simulated dispersion curves that are consistent with the measured data (i.e., in thiscase, smooth and continuous). Once a curve has been simulated, it is inverted to obtain a single, “best”fit model, as shown in Figure 6d. Note that while all previous procedures to account for Vs uncertaintyrequired the use of a global-search algorithm to ultimately develop suites of Vs profiles, in contrast, thisnew procedure is inversion-algorithm agnostic, thereby permitting the analyst to use any algorithm of theirchoice (global or local). Furthermore, by selecting only the single, “best” model the issue of dependence onthe number of trial models is avoided altogether, provided of course that a reasonable/sufficient number ofmodels has been attempted to produce a good fit to the simulated dispersion curve, which is easily verifiable.Returning to Figure 6d, it is important to note that the fit may not be perfect, however, the proposed proce-dure does not require it. The issue of primary importance is that the simulated experimental dispersion dataand the theoretical fit are in good agreement over the experimental frequency range, though not necessarilystrictly identical. The procedure of simulating a realization of the dispersion data and fitting it throughinversion is repeated for some number of trials (N), which is the only input parameter to be defined by theanalyst. The best inversion-derived fit to each realization and the corresponding ground model is saved as apotential solution to the experimental data considering its uncertainties, refer to Figure 6c. The N fits (i.e.,N theoretical dispersion curves) can then be used to calculate implied dispersion data statistics which canbe compared directly against the experimental dispersion data statistics to assess the successfulness of theprocedure. Note that the value of N is expected to be problem dependent, however, a minimum acceptablevalue for N can be checked prior to performing any inversions by comparing the statistics of the simulateddispersion realizations and the experimental dispersion data. If the statistics of the simulated dispersiondata (i.e., the targets used for inversion) are unable to reproduce the experimental statistics, it can then be11nferred that the inversion-derived fits to those simulated curves will also not reproduce the experimentalstatistics. Thus, care must be taken to select N to be sufficiently large to ensure agreement between thesimulated and measured experimental data, although not too large to avoid excessive computational costwhen performing N different inversions. An N=250 was shown to perform well over the course of this study,however the choice of N is left to the discretion of the analyst. Importantly, following the inversion of theN simulated dispersion curves to obtain N theoretical dispersion curves, the statistics implied by the Ntheoretical dispersion curves (i.e., the implied dispersion data, refer to Figure 6c), must be compared quan-titatively with the measured experimental dispersion data in a manner similar to that shown in Figure 5. Ifthe agreement between the implied and measured experimental dispersion data is found to be unsatisfactory,additional simulations and/or increasing the value of N can be used to improve the results. The remainder ofthis paper will address the application of this new approach to two synthetic tests and a real-world example.
To demonstrate the effectiveness of the newly proposed procedure, we will test it using the same syntheticdataset as that was used previously to evaluate the literature-based approaches for accounting for Vs un-certainty. Recall, the experimental dispersion data (i.e., mean and standard deviation as a function offrequency) for this dataset is presented in Figure 1c. As the new approach requires the correlations betweenfrequencies, and the synthetic data, as developed, does not include such information, the correlations hadto be synthesized. The procedure for synthesizing the correlations involved first simulating ground modelswhose theoretical dispersion curves were consistent with the experimental dispersion’s uncertainty. Thestatistics used to inform the simulation of these ground models were developed based on the mean valuesfrom an M0-type inversion and the uncertainty from an M1-type inversion, and while these types of inversionshave been shown to be lacking in their ability to propagate dispersion uncertainty (refer to Figures 4 and5), they provided a reasonable correlation structure. The synthesized correlations were combined with theknown/assumed statistics (refer to Figure 1c) to define the experimental dispersion data. With this infor-mation, the procedure outlined in Figure 6 was performed using N=250 realizations. The choice of N=250was made prior to performing the inversions by checking that the statistics of the realizations consistentlyreproduced the statistics of the experimental dispersion data. Each realization was inverted using the same3-layer parameterization and 10,000 random plus 50,000 neighborhood algorithm trial models, as was usedin the previous section. Inversions were performed on the Texas Advanced Computing Center’s (TACCs)cluster Stampede2 using a single Skylake (SKX) node. The entire analysis (i.e., the inversion of all N=250simulated dispersion curves) took less than 2 hours to complete.The results from the inversion analyses are provided in Figures 7, 8, and 9. Figure 7a compares thetheoretical dispersion curves fit to the N=250 realizations of the experimental dispersion data with theoriginal experimental dispersion data. Note the theoretical dispersion curves have been colored in termsof their dispersion misfit relative to the original experimental data and not their respective realizations.In contrast to the previous methods (refer to Figure 3) we observed a reasonable mixture of theoreticaldispersion curves with low misfits (i.e., < .0 0.2 0.4 0.6 0.8 1.0Dispersion Misfit3 10 V r ( m / s ) (a) Exp. Disp. Data 275 325Vr (m/s) (b) ~3Hz 175 215Vr (m/s) (c) ~6Hz95 110Vr (m/s) (d) ~13Hz 85 100Vr (m/s) (e) ~30Hz P r o b a b ili t y D e n s i t y Figure 7: Qualitative assessment of the newly proposed procedure for propagating experimental disper-sion uncertainty into the inverted Vs profiles considering a single inversion parameterization for asynthetic example . Panel (a) shows the experimental dispersion data and 250 inverted theoretical disper-sion curves fit to the 250 realizations of the experimental dispersion data. The theoretical dispersion curveshave been colored according to their dispersion misfit values relative to the original experimental dispersiondata, not their respective realizations. The vertical dashed lines in panel (a) at approximately 3, 6, 13, and30 Hz denote the location of the “slices” shown in panels (b), (c), (d), and (e), respectively. These “slices”compare the measured experimental dispersion data’s distribution (solid black line) with the distribution ofRayleigh wave velocity (Vr) derived from inversion (histogram).illustrate excellent agreement between the experimental and inverted dispersion data’s mean and uncertainty(i.e., residuals approximately equal to zero at all frequencies), especially when compared to previous methods(refer to Figure 5).With the theoretical dispersion curves from the inverted ground models having been shown to be con-sistent with the uncertainty of the experimental dispersion data, we can now examine the effects on the Vsprofiles. Figure 9a presents the single-lowest misfit Vs profile from the inversion of each of the N=250 disper-sion realizations. To better illustrate their concentration about the true solution the Vs profiles have beendiscretized in terms of depth and Vs, binned into cells, and color mapped in terms of the number of profilesin each cell. The discretized lognormal median Vs profile and the 95% lognormal confidence interval (CI)profiles are also shown for reference. Figure 9a shows that the discretized median profile reasonably capturesthe true solution and that the 95% confidence interval qualitatively captures the variance in the Vs profiles,with the notable exception of the sharp spikes at layer boundaries. When we compare σ ln,V s from the newapproach (i.e., Figure 9b) directly with those from the literature (i.e., Figure 2b, d, f, and h) we observethat the new approach shows increased uncertainty from the (b100) alternatives, as should be expected,but surprisingly slightly less uncertainty than the (n100) alternatives. This observation indicates that (atleast for this example) the n100 models tend to over-estimate Vs uncertainty even when under-predictingthe dispersion data uncertainty (refer to Figure 5). Further examination of Figure 9b shows that σ ln,V s forthis example is approximately 0.06, which is quite low compared to values commonly assumed in practice(EPRI, 2012; Stewart et al., 2014; Toro, 1995). This is certainly due in part to the use of only a singleinversion layering parameterization, as others have shown that the variability within a parameterization isgenerally much less than that between parameterizations. Hence, if the true subsurface layering is unknowna priori, it is not simply enough to accurately represent the experimental dispersion data’s uncertainty usinga single, assumed layering parameterization. Rather, one must also incorporate the epistemic uncertaintyin the layering parameterization itself (Cox and Teague, 2016; Di Giulio et al., 2012; Vantassel and Cox,13 μ r e s / μ m e a s ( % ) (a) Overestimate Mean VrUnderestimate Mean Vr δ i n v − δ m e a s (b) Overestimate Vr UncertaintyUnderestimate Vr Uncertainty
Figure 8: Quantitative assessment of the newly proposed procedure for propagating experimental disper-sion uncertainty into the inverted Vs profiles considering a single inversion parameterization for asynthetic example . The comparison is made on the basis of: (a) the dispersion residual mean ( µ res ) [i.e.,the difference between the mean of the inverted theoretical dispersion curves ( µ inv ) and the mean of themeasured experimental dispersion data ( µ meas )] normalized by µ meas and expressed in percent, (b) the dis-persion residual coefficient of variation [i.e., the difference between the coefficient of variation of the invertedtheoretical dispersion curves ( δ inv ) and the coefficient of variation of the measured experimental dispersiondata ( δ meas )]. The vertical dashed lines in panels (a) and (b) denote the location of the “slices” shown inFigure 7.2020b). The incorporation of multiple inversion parameterizations into the procedure and their effects onthe Vs profiles are presented in the following section. Uncertainty within surface wave inversion is generally split into intra- and inter-parameterization variability,which address the uncertainty inside a single inversion parameterization (i.e., the previous example) andbetween various parameterizations (i.e., current example), respectively. Previous studies have shown thatthe inter-parameterization variability is generally more significant than the inter-parametrization variability.However, as the methods used in previous studies to account for intra-parameterization variability tendedto underestimate the experimental uncertainty (refer to Figures 4 and 5), it is of interest to reexamine thisconjecture using the new procedure.The new procedure was repeated using the same synthetic experimental dispersion data as the previousexample (refer to Figure 1c), except rather than using only a single layering parameterization consistingof 3 layers, five different layering parameterizations were used to account for epistemic uncertainty in theparameterization selection. The five parameterizations were based on the Layering by Number (LN) schemaand included LNs of 3, 5, 7, 9, and 14. These parameterizations were selected because they were all deemedacceptable when this same set of experimental dispersion data was inverted by Vantassel and Cox (2020b).In practice, the analyst must decide carefully which parameterizations to pursue, as this can strongly impactthe resulting inversion-derived Vs profiles. We recommend investigating a number of trial parameterizations,and the reader is referred to the previous study by Vantassel and Cox (2020b) for a more thorough discus-sion on the number and type of parameterizations to use in their inversions. The five parameterizationsnoted above were used to invert different sets of N=250 dispersion realizations (i.e., 1250 realizations in14igure 9: Uncertainty-consistent Vs profiles considering a single inversion parameterization for asynthetic example using the newly proposed procedure . In order to better illustrate their concentration,the profiles have been discretized in terms of depth and Vs, binned into cells, and color mapped in termsof the number of profiles in each cell. Panel (a) summarizes the Vs profiles resulting from the inversion ofthe N=250 realizations of the experimental dispersion data alongside the true solution, the discretized log-normal median, and the 95% confidence interval (CI). Panel (b) illustrates the lognormal standard deviationof Vs ( σ ln,V s ) of the 250 profiles shown in panel (a). The sharp spikes in σ ln,V s are due to the uncertaintyin the profile’s layer boundaries, reflecting a shortcoming in how σ ln,V s has been calculated historically, andare not the result of actual uncertainty in Vs directly.total). N=250 realizations per parameterization was selected to ensure that each parameterization wouldhave a sufficient number of realizations to each individually capture the uncertainty of the experimentaldispersion data, as it is important to not discount the intra-parameterization variability when investigatingthe inter-parameterization variability. This, of course, is not the only approach available to consider theinter-parameterization variability. An alternative method examined during the course of this study is toinvert a single set of N realizations with various parameterizations. This alternate approach was found toproduce similar results to those using the presented method, however, since the presented method is believedto be more robust, as it ensures each parameterization remains separate, it is the one selected here.The results from the newly proposed procedure when accounting for multiple layering parameterizationsare shown in Figures 10, 11, and 12. Figure 10a compares the theoretical dispersion curves fit to the1250 realizations of the experimental dispersion data with the original experimental dispersion data. Thetheoretical dispersion curves are colored in terms of their dispersion misfit relative to the original experimentaldispersion data (not their respective realizations) and are again shown to follow a reasonable distribution oflow misfit ( < .0 0.2 0.4 0.6 0.8 1.0Dispersion Misfit3 10 V r ( m / s ) (a) Exp. Disp. Data 275 325Vr (m/s) (b) ~3Hz 175 215Vr (m/s) (c) ~6Hz95 110Vr (m/s) (d) ~13Hz 85 100Vr (m/s) (e) ~30Hz P r o b a b ili t y D e n s i t y Figure 10: Qualitative assessment of the newly proposed procedure for propagating experimental dis-persion uncertainty into the inverted Vs profiles considering multiple inversion parameterizationsfor a synthetic example . Panel (a) shows the experimental dispersion data and the inverted theoreticaldispersion curves fit to the 1250 realizations (5 parameterizations with N=250 realizations each) of the ex-perimental dispersion data. The theoretical dispersion curves have been colored according to their dispersionmisfit values relative to the original experimental dispersion data, not their respective realizations. The ver-tical dashed lines in panel (a) at approximately 3, 6, 13, and 30 Hz denote the location of the slices shown inpanels (b), (c), (d), and (e), respectively. These slices compare the measured experimental dispersion data’sdistribution (solid black line) with the distribution of Rayleigh wave velocity (Vr) derived from inversion(histogram).dispersion realizations (i.e., N=250 profiles per parameterization * 5 parameterizations = 1250 profiles). Tobetter illustrate their concentration, the Vs profiles have been discretized in terms of depth and Vs, binnedinto cells, and color mapped in terms of the number of profiles in each cell. The inverted Vs profiles areshown alongside the true solution, discretized lognormal median Vs profile, and lognormal 95% CI profiles.As expected, the inverted profiles qualitatively show more uncertainty than when a single parameterizationis considered (refer to Figure 9a for comparison). This distinction is most apparent at the layer boundaries,as the locations of these boundaries are controlled predominantly by the assumed layering parameterization.Parameterizations with many layers will tend to result in smoother profiles with gradual changes in Vs,while parameterizations with only a few layers will tend to result in profiles with sharper contrasts (Coxand Teague, 2016). As the true site layering may not be known prior to surface wave inversion, there isoften a need to consider multiple parameterizations with different numbers of layers to properly address oneof surface wave inversion’s main sources of epistemic uncertainty. Figure 12b presents σ ln,V s of the 1250profiles. The uncertainty across multiple parameterizations is shown to increase in regards to what wasobserved previously for a single parameterization (refer to Figure 9b for comparison). Specifically, σ ln,V s is now closer to 0.1 within any given layer, thereby confirming quantitatively what was already observedqualitatively in Figure 12a. The increase in σ ln,V s lends confidence to the conjecture that the uncertaintyresulting from multiple parameterizations will tend to exceed that for a single parameterization, althoughperhaps not to the same extent as may have been inferred previously due the tendencies to underestimateintra-parameterization variability in previous studies. Regardless, at least for this synthetic example, σ ln,V s remains quite low compared to values commonly assumed in practice (EPRI, 2012; Stewart et al., 2014;Toro, 1995). 16 μ r e s / μ m e a s ( % ) (a) Overestimate Mean VrUnderestimate Mean Vr δ i n v − δ m e a s (b) Overestimate Vr UncertaintyUnderestimate Vr Uncertainty
Figure 11: Quantitative assessment of the newly proposed procedure for propagating experimental dis-persion uncertainty into the inverted Vs profiles considering multiple inversion parameterization fora synthetic example . The comparison is made on the basis of: (a) the dispersion residual mean ( µ res )[i.e., the difference between the mean of the inverted theoretical dispersion curves ( µ inv ) and the mean of themeasured experimental dispersion data ( µ meas )] normalized by µ meas and expressed in percent, (b) the dis-persion residual coefficient of variation [i.e., the difference between the coefficient of variation of the invertedtheoretical dispersion curves ( δ inv ) and the coefficient of variation of the measured experimental dispersiondata ( δ meas )]. The vertical dashed lines in panels (a) and (b) denote the location of the “slices” shown inFigure 10. To illustrate the application of the newly proposed method at a real site, we present the following casestudy. For reference, the site is located in the southern California and has a known geology consistingof alluvial soils (sands and clays) approximately 15 – 25 m thick, overlaying highly variable weatheredbedrock, over competent bedrock. While an extensive set of active-source and passive-wavefield surface wavemeasurements are available at this site, to keep the example simple, we will only describe how dispersion datawas extracted from some of the passive-wavefield microtremor array measurements (MAM). To incorporatealeatory variability across the 70m by 70m site into the estimates of dispersion uncertainty, four 35m diametercircular arrays each composed of nine sensors will be considered. In essence, each of the four MAM arrayswere responsible for characterizing one quarter of the overall site. While using only one size of MAMarray will result in somewhat bandlimited dispersion data, we do so herein to maintain clarity in how thecombined aleatory and epistemic uncertainties are estimated. MAM array sensors were buried and leftundisturbed to record ambient noise, seven hours of which were used to extract dispersion data from eacharray. Each of the four MAM arrays were processed separately to obtain experimental dispersion data usingthe Rayleigh three-component beamforming method developed by Wathelet et al. (2018) and implementedin the open-source software Geopsy (Wathelet et al., 2020). Constant length time blocks were preferredover frequency-dependent time blocks to facilitate statistical calculations on the experimental dispersiondata. Time blocks were selected to be 60 seconds long to ensure at least 100 cycles at the lowest processingfrequency. The recommended 4 time blocks per sensor per block set were used to compute the average cross-correlation matrices. As sufficiently long noise records were available, no overlap was permitted betweenblocks or block sets. No higher modes were apparent in preliminary dispersion processing, so only thespectral peaks with the single highest power at each frequency were selected from each block set to representthe experimental dispersion data. This resulted in 11 estimates of Vr at each frequency for each array. Or,17igure 12: Uncertainty-consistent Vs profiles considering multiple inversion parameterization for asynthetic example using the newly proposed procedure . In order to better illustrate their concentration,the profiles have been discretized in terms of depth and Vs, binned into cells, and color mapped in termsof the number of profiles in each cell. Panel (a) summarizes the Vs profiles resulting from the inversion of1250 realizations (5 parameterizations with N=250 realizations each) of the experimental dispersion dataalongside the true solution, the discretized log-normal median, and the 95% confidence interval (CI). Panel(b) illustrates the lognormal standard deviation of Vs ( σ ln,V s ) of the 250 profiles shown in panel (a). Thesharp spikes in σ ln,V s are due to the uncertainty in the profile’s layer boundaries, reflecting a shortcomingin how σ ln,V s has been calculated historically, and are not the result of actual uncertainty in Vs directly.in other words, 44 estimates of Vr per frequency across the entire site. The experimental data was binnedand resampled in terms of log-wavelength following the recommendations of Vantassel and Cox (2020b). Theraw experimental dispersion data and the resulting statistical representation is shown in Figure 13a. Thefrequency-dependent COVs tend to increase from approximately 0.02 at short wavelengths to approximately0.08 at long wavelengths. The minimum and maximum wavelength of the experimental dispersion data areapproximately 10m and 250m, respectively. The reader will note that the maximum wavelength (250m) issubstantially longer than what one would generally expect to be able to resolve based on the array resolutionlimits (typically 2 or 3 times the maximum array aperture). However, this was done only after confirming thatthe dispersion data from these arrays were consistent with data from larger aperture arrays also performedat the site (not shown). The minimum depth of profiling can be approximated by dividing the minimumwavelength by a depth factor (df) of 3 or 2 (Foti et al., 2018) to obtain an estimate of the minimum resolvablethickness of the near-surface layer between 3m and 5m. These relatively large near-surface layers are thedirect result of using only passive arrays of medium size for this example. The maximum depth of profilingcan be approximated by dividing the maximum wavelength by a df of 3 or 2. However, the use of even a df of3 may be optimistic when the dispersion curve has an “L”-shape which does not flatten at long wavelengths(Vantassel and Cox, 2020b). Nonetheless, we adopt a df of 3 for parameterizing the inversion’s maximumdepth based on the experimental dispersion data, and then limit the depth to which the inverted Vs profilesare presented based on the quality of the resulting Vs profiles (i.e., depth limited based on high σ ln,V s valueswhen resolution is poor).The frequency-dependent uncertainty in the experimental dispersion data was propagated through theinversion process using the newly proposed procedure outlined in Figure 6. To incorporate the additionalepistemic uncertainty from the inversion’s parameterization, five LN-type parameterizations using LN = 3,5, 7, 9, and 14 were considered following the same procedure as presented in the previous synthetic example.Inversions took 2 hours to complete using 5 SKX nodes on the Stampede2 cluster. After inversion, the18 Wavelength, λ (m)2505007501000 V r ( m / s ) (a) Raw Disp. DataExp. Disp. Data 3 10 Frequency, f (Hz) (b)
Exp. Disp. Data 0.00.20.40.60.81.0 D i s p e r s i o n M i s f i t Figure 13: Dispersion data used to study the application of the newly proposed procedure for propagatingexperimental dispersion uncertainty into the inverted Vs profiles considering multiple parameterizationsat a real site : (a) raw experimental dispersion data with its binned and resampled statistical representation,and (b) agreement between the experimental dispersion data and 1000 inverted theoretical dispersion curvesfit to the 1000 realizations (4 acceptable parameterizations with N=250 realizations each) of the experimentaldispersion data. The theoretical dispersion curves have been colored according to their dispersion misfitvalues relative to the original experimental dispersion data, not their respective realizations. The verticaldashed lines in panel (b) indicate the location of the “slices” shown in Figure 14.parameterization-quality criteria proposed by Vantassel and Cox (2020b) was used to assess the performanceof the five parameterizations. The criteria showed the LN=3 parameterization underperformed its coun-terparts and was therefore removed from further consideration to avoid biasing the results. The resulting1000 theoretical dispersion curves (i.e., 250 curves/parameterization * 4 acceptable parameterizations =1000 curves) that were fit to the 1000 experimental dispersion data realizations are shown alongside theexperimental dispersion data in Figure 13b. Note that the theoretical curves in Figure13b have been coloredin terms of their dispersion misfit relative to the original experimental data, and not their respective realiza-tions. Again, we observe that the majority of the curves have a misfit less than 1.0 (i.e., are on average within1 standard deviation of the mean), with quite a few below 0.2 that closely follow the mean experimentaldispersion trend. Importantly, however, a number of curves with misfits greater than 1.0 are also present, asshould be expected if the dispersion uncertainty is accurately being represented in the inversion results. Thisleads to a qualitative comparison that the experimental dispersion data’s uncertainty and the uncertaintyimplied by the suites of theoretical dispersion curves derived from the inversion are in excellent agreement.Figure 14 presents a quantitative comparison between the uncertainty of the experimental dispersion dataand inverted theoretical dispersion curves, highlighting four different “slices” using vertical colored dashedlines at approximately 3, 6, 13, and 30 Hz. The quantitative comparison confirms that the inversion-derivedmodels are accurately capturing the experimental dispersion uncertainty in terms of their mean (Figure14a), standard deviation (Figure 14d), and distribution shapes (Figures 14b, c, e, and f). It is importantto note that the minor discrepancies between the measured and inverted data’s mean value and uncertainty(refer to Figure 14a and 14d) are to be expected, because unlike synthetic data, real experimental datamay have some error in its mean and is likely to have variable frequency-dependent uncertainty, making theexact replication of some realizations difficult. For example, the dispersion COVs at this site varied between0.02 – 0.08, unlike the synthetic examples which had constant dispersion COVs of 0.05. However, despitethese complicating factors the new procedure is able to produce results that fit the experimental dispersiondata remarkably well, with good agreement being observed between the measured and inverted dispersionuncertainty across the entire bandwidth of interest.Figure 15 shows the inversion-derived, uncertainty-consistent Vs profiles at two depth scales; Figure19 μ r e s / μ m e a s ( % ) (a) O(ere%ti ate Mean VrUndere%ti ate Mean Vr Frequenc), f (Hz),0.050.000.05 δ i n v , δ m e a s (d) O(eresti ate Vr Uncertaint)Underesti ate Vr Uncertaint)
600 800Vr ( /s) (b) ~3Hz 255 265Vr ( /s) (c) ~6Hz190 200Vr ( /s) (e) ~9Hz 180 190Vr ( /s) (f) ~14Hz P r o b a b ili t ) D e n s i t ) Figure 14: Quantitative assessment of the newly proposed procedure for propagating experimental disper-sion uncertainty into the inverted Vs profiles considering multiple parameterizations at a real site .The comparison is made on the basis of: (a) the residual mean ( µ res ) [i.e., the difference between the meanof the inverted theoretical dispersion curves ( µ inv ) and the mean of the measured experimental dispersiondata ( µ meas )] normalized by µ meas in percent, (d) the residual coefficient of variation [i.e., the differencebetween the coefficient of variation of the inverted theoretical dispersion curves ( δ inv ) and the coefficientof variation of the measured experimental dispersion data ( δ meas )]. The vertical dashed lines in panels (a)and (d) at approximately 3, 6, 9, and 14 Hz denote the locations of the “slices” shown in panels (b), (c),(e), and (f), respectively. These “slices” compare the measured experimental dispersion data’s distribution(solid black line) with the distribution of Rayleigh wave velocity (Vr) derived from inversion (histogram).15a to a depth of 60m and Figure 15b to a depth of 30m. The Vs profiles reveal a three-layered systemconsistent with the anticipated geologic conditions and comprised of approximately 15m of soft soil (Vs ofapproximately 200 m/s), overlying at least 35m of weathered rock (Vs of approximately 600 m/s), overlyingstiffer material. The uncertainty in the Vs profiles is quantified using σ ln,V s and presented in Figures 15band 15d at the same depth scales as those for Vs (i.e., Figures 15a and 15c). σ ln,V s increases from closeto zero for the surface layer, where the high frequency/short wavelength surface-wave dispersion data hadlow COVs, to a σ ln,V s of 0.6 at depth, where due to the lack of very low frequency/long wavelength datacoupled with relatively high dispersion COVs the velocity of the material cannot be accurately resolved.We consider the depth where the Vs uncertainty becomes so large as to make meaningful inferences aboutthe site’s stiffness intractable as the maximum depth limitation of the presented Vs profiles (i.e., a depthof approximately 45m, where σ ln,V s starts to exceed about 0.2). To illustrate the accuracy of the surfacewave measurements, a P-S-suspension log from the site is shown in comparison to the Vs profiles in Figures15a and 15c. The P-S-suspension log is seen to reside solely within the 95% CI of the inverted Vs profiles,despite two complicating factors. First, the inversion-derived Vs profiles are the result of measurementsmade over a 70m by 70m area with known spatial variability, whereas the P-S-suspension log is more-or-lessa point measurement. And second, the P-S-suspension log was performed at the western edge of the 70m by70m site rather than at its center due to site constraints. When these two factors are combined, one wouldexpect some differences between the invasive and non-invasive Vs profiles. Nonetheless, the agreement isquite good and the Vs uncertainty associated with the surface wave profiles has been robustly quantifiedfor use in subsequent engineering analyses, while that of the single P-S log is unknown and would have to20igure 15: Uncertainty-consistent Vs profiles considering multiple inversion parameterizations at areal site using the newly proposed procedure . Panels (a) and (c), summarize the Vs profiles resultingfrom the inversion of the 1000 realized experimental dispersion curves (4 parameterizations with N=250realizations each) at two separate depth scales. In order to better illustrate their concentration, the profileshave been discretized in terms of depth and Vs, binned into cells, and color mapped in terms of the numberof profiles in each cell. The profiles are shown alongside their discretized log-normal median and 95%confidence interval (CI) and a P-S-suspension log from nearby. Panels (b) and (c), show the lognormalstandard deviation of Vs ( σ ln,V s ) of the 1000 profiles shown in panel (a) and (c), respectively. The sharpspike in σ ln,V s at approximately 15m is due to the uncertainty in the profile’s layer boundaries, reflecting ashortcoming in how σ ln,V s has been calculated historically, and is not the result of actual uncertainty in Vsdirectly.be assumed. This case study illustrates that the proposed procedure can rigorously propagate experimentaldispersion uncertainty through the surface wave inversion process and produce suites of Vs profiles that,while uncertain, can produce results that are consistent with invasive characterization methods. Surface wave methods are increasingly being preferred over traditional site characterization methods formeasuring a site’s shear-stiffness due to their non-invasive nature. The quantification and propagation ofuncertainty from non-invasive measurements into Vs profiles for use in subsequent engineering analyses hasbeen the focus of much work in recent years, and while significant progress has been made, no study has yetbeen able to show quantitative evidence of the successful propagation of experimental dispersion uncertaintyinto the resulting Vs profiles. Therefore, this study began by examining methods currently available from theliterature and showing their deficiencies in three specific categories. First, current methods are shown to behighly sensitive to their many user-defined inversion input parameters (e.g., number of global-inversion trialmodels), making it difficult/impossible for them to be performed repeatedly by different analysts. Second, thesuites of inverted Vs profiles derived from these methods, when viewed in terms of their implied theoreticaldispersion data, are shown to systematically under-estimate the uncertainty present in the experimentaldispersion data, though they may appear satisfactory when viewed purely qualitatively. And third, if theuncertainties in the implied theoretical dispersion data were to be examined quantitatively, which has notbeen done previously, there is no obvious remedy available to the analyst to resolve any inconsistency betweenthe measured and inverted dispersion uncertainty. Therefore, a new approach has been proposed that seeksto remedy these shortcomings. First, in addition to the necessary considerations for any inversion, the21ethod is governed by only one additional user-defined input parameter (N), which controls the numberof dispersion curve realizations. A value of N=250 has been shown to perform well throughout this study,however, as N is expected to be problem dependent the selection of an appropriate value of N is left to theanalyst. Second, the new method requires the quantitative comparison between the measured and inverteddispersion uncertainties to ensure suites of Vs profiles which quantitatively reproduce the uncertainties inthe experimental dispersion data. And third, should the measured and inverted dispersion uncertainty bein disagreement, specific guidance is provided on the actions necessary (e.g., increasing the value of N)to produce Vs profiles that better account for the experimental uncertainty. The application of the newprocedure has been demonstrated using two synthetic tests and a real-world example. Results show theprocedure’s ability to produce suites of Vs profiles which accurately capture the site’s Vs structure, whilerigorously propagating the measured, site-specific dispersion data’s uncertainty through the inversion process,yielding uncertainty-consistent Vs profiles that can be used in subsequent engineering analyses.
References
Ronald D. Andrus and Kenneth H. Stokoe. Liquefaction Resistance of Soils from Shear-Wave Velocity.
Journal of Geotechnical and Geoenvironmental Engineering , 126(11):1015–1025, November 2000. ISSN1090-0241, 1943-5606. doi: 10.1061/(ASCE)1090-0241(2000)126:11(1015).B. R. Cox and C. M. Wood. Surface Wave Benchmarking Exercise: Methodologies, Results, and Uncertain-ties. In
GeoRisk 2011 , pages 845–852, Atlanta, Georgia, United States, June 2011. American Society ofCivil Engineers. ISBN 978-0-7844-1183-4. doi: 10.1061/41183(418)89.B. R. Cox, C. M. Wood, and David P. Teague. Synthesis of the UTexas1 Surface Wave Dataset Blind-Analysis Study: Inter-Analyst Dispersion and Shear Wave Velocity Uncertainty. In
Geo-Congress 2014Technical Papers , pages 850–859, Atlanta, Georgia, February 2014. American Society of Civil Engineers.ISBN 978-0-7844-1327-2. doi: 10.1061/9780784413272.083.Brady Cox and Joseph Vantassel. Dynamic Characterization of Wellington, New Zealand.
DesignSafe-CI[publisher] , November 2018. doi: 10.17603/DS24M6J. Type: dataset.Brady R. Cox and David P. Teague. Layering ratios: a systematic approach to the inversion of surface wavedata in the absence of a priori information.
Geophysical Journal International , 207(1):422–438, October2016. ISSN 0956-540X, 1365-246X. doi: 10.1093/gji/ggw282.Michael R. Deschenes, Clinton M. Wood, Liam M. Wotherspoon, Brendon A. Bradley, and Ethan Thomson.Development of Deep Shear Wave Velocity Profiles in the Canterbury Plains, New Zealand.
EarthquakeSpectra , 34(3):1065–1089, August 2018. ISSN 8755-2930. doi: 10.1193/122717EQS267M.22iuseppe Di Giulio, Alexandros Savvaidis, Matthias Ohrnberger, Marc Wathelet, Cecile Cornou, BrigitteKnapmeyer-Endrun, Florence Renalier, Nikos Theodoulidis, and Pierre-Yves Bard. Exploring the modelspace and ranking a best class of models in surface-wave dispersion inversion: Application at Europeanstrong-motion sites.
Geophysics , 77(3):B147–B166, May 2012. ISSN 0016-8033, 1942-2156. doi: 10.1190/geo2011-0116.1.EPRI. Seismic Evaluation Guidance: Screening, Prioritization and Implementation Details (SPID) for theResolution of Fukushima Near-Term Task Force Recommendation 2.1: Seismic. Technical Report 1025287,Palo Alto, CA, November 2012.S. Foti, C. Comina, D. Boiero, and L.V. Socco. Non-uniqueness in surface-wave inversion and consequenceson seismic site response analyses.
Soil Dynamics and Earthquake Engineering , 29(6):982–993, June 2009.ISSN 02677261. doi: 10.1016/j.soildyn.2008.11.004.Sebastiano Foti, Carlo G. Lai, Glenn J. Rix, and Claudio Strobbia.
Surface Wave Methods for Near-SurfaceSite Charachterization . CRC Press, Boca Raton, FL, 1 edition, 2015.Sebastiano Foti, Fabrice Hollender, Flora Garofalo, Dario Albarello, Michael Asten, Pierre-Yves Bard, Ce-sare Comina, Ccile Cornou, Brady Cox, Giuseppe Di Giulio, Thomas Forbriger, Koichi Hayashi, EnricoLunedei, Antony Martin, Diego Mercerat, Matthias Ohrnberger, Valerio Poggi, Florence Renalier, Deb-orah Sicilia, and Valentina Socco. Guidelines for the good practice of surface wave analysis: a productof the InterPACIFIC project.
Bulletin of Earthquake Engineering , 16(6):2367–2420, June 2018. ISSN1570-761X, 1573-1456. doi: 10.1007/s10518-017-0206-7.F. Garofalo, S. Foti, F. Hollender, P.Y. Bard, C. Cornou, B.R. Cox, A. Dechamp, M. Ohrnberger, V. Perron,D. Sicilia, David P. Teague, and C. Vergniault. InterPACIFIC project: Comparison of invasive and non-invasive methods for seismic site characterization. Part II: Inter-comparison between surface-wave andborehole methods.
Soil Dynamics and Earthquake Engineering , 82:241–254, March 2016. ISSN 02677261.doi: 10.1016/j.soildyn.2015.12.009.Shawn C. Griffiths, Brady R. Cox, Ellen M. Rathje, and David P. Teague. Mapping Dispersion Misfitand Uncertainty in Vs Profiles to Variability in Site Response Estimates.
Journal of Geotechnical andGeoenvironmental Engineering , 142(11):04016062, November 2016a. ISSN 1090-0241, 1943-5606. doi:10.1061/(ASCE)GT.1943-5606.0001553.Shawn C. Griffiths, Brady R. Cox, Ellen M. Rathje, and David P. Teague. Surface-Wave Dispersion Ap-proach for Evaluating Statistical Models That Account for Shear-Wave Velocity Uncertainty.
Journal ofGeotechnical and Geoenvironmental Engineering , 142(11):04016061, November 2016b. ISSN 1090-0241,1943-5606. doi: 10.1061/(ASCE)GT.1943-5606.0001552.Fabrice Hollender, Ccile Cornou, Aline Dechamp, Kaveh Oghalaei, Florence Renalier, Emeline Maufroy,Clment Burnouf, Sylvette Thomassin, Marc Wathelet, Pierre-Yves Bard, Vincent Boutin, Clment Desbor-des, Isabelle Douste-Bacqu, Laetitia Foundotos, Cdric Guyonnet-Benaize, Vincent Perron, Julie Rgnier,Agathe Roull, Mickael Langlais, and Deborah Sicilia. Characterization of site conditions (soil class, VS30,velocity profiles) for 33 stations from the French permanent accelerometric network (RAP) using surface-wave methods.
Bulletin of Earthquake Engineering , 16(6):2337–2365, June 2018. ISSN 1570-761X, 1573-1456. doi: 10.1007/s10518-017-0135-5.J. D. Hunter. Matplotlib: A 2D Graphics Environment.
Computing in Science Engineering , 9(3):90–95, May2007. ISSN 1558-366X. doi: 10.1109/MCSE.2007.55.R. Kayen, R. E. S. Moss, E. M. Thompson, R. B. Seed, K. O. Cetin, A. Der Kiureghian, Y. Tanaka,and K. Tokimatsu. Shear-Wave VelocityBased Probabilistic and Deterministic Assessment of SeismicSoil Liquefaction Potential.
Journal of Geotechnical and Geoenvironmental Engineering , 139(3):407–419,March 2013. ISSN 1090-0241, 1943-5606. doi: 10.1061/(ASCE)GT.1943-5606.0000743.23arlo G. Lai, Sebastiano Foti, and Glenn J. Rix. Propagation of Data Uncertainty in Surface Wave Inversion.
Journal of Environmental and Engineering Geophysics , 10(2):219–228, June 2005. ISSN 1083-1363. doi:10.2113/JEEG10.2.219.Federico Passeri, Sebastiano Foti, Brady R. Cox, and Adrian Rodriguez-Marek. Influence of EpistemicUncertainty in Shear Wave Velocity on Seismic Ground Response Analyses.
Earthquake Spectra , 35(2):929–954, May 2019. ISSN 8755-2930, 1944-8201. doi: 10.1193/011018EQS005M.Federico Passeri, Sebastiano Foti, and Adrian Rodriguez-Marek. A new geostatistical model for shear wavevelocity profiles.
Soil Dynamics and Earthquake Engineering , 136:106247, 2020. ISSN 0267-7261. doi:https://doi.org/10.1016/j.soildyn.2020.106247.Ellen M. Rathje, Albert R. Kottke, and Whitney L. Trent. Influence of Input Motion and Site Property Vari-abilities on Seismic Site Response Analysis.
Journal of Geotechnical and Geoenvironmental Engineering ,136(4):607–619, April 2010. ISSN 1090-0241, 1943-5606. doi: 10.1061/(ASCE)GT.1943-5606.0000255.Ellen M. Rathje, Clint Dawson, Jamie E. Padgett, Jean-Paul Pinelli, Dan Stanzione, Ashley Adair, Pe-dro Arduino, Scott J. Brandenberg, Tim Cockerill, Charlie Dey, Maria Esteva, Fred L. Haan, MatthewHanlon, Ahsan Kareem, Laura Lowes, Stephen Mock, and Gilberto Mosqueda. DesignSafe: New Cyber-infrastructure for Natural Hazards Engineering.
Natural Hazards Review , 18(3):06017001, August 2017.ISSN 1527-6988, 1527-6996. doi: 10.1061/(ASCE)NH.1527-6996.0000246.Malcolm Sambridge. Geophysical inversion with a neighbourhood algorithm-I. Searching a parameter space.
Geophysical Journal International , 138(2):479–494, August 1999. ISSN 0956540X, 1365246X. doi: 10.1046/j.1365-246X.1999.00876.x.Jonathan P Stewart, Kioumars Afshari, and Youssef M A Hashash. Guidelines for Performing Hazard-Consistent One-Dimensional Ground Response Analysis for Ground Motion Prediction. Technical Report2014/16, PEER, Pacific Earthquake Engineering Research Center Headquarters, University of Californiaat Berkeley, October 2014.David P. Teague and Brady R. Cox. Site response implications associated with using non-unique Vs profilesfrom surface wave inversion in comparison with other commonly used methods of accounting for Vsuncertainty.
Soil Dynamics and Earthquake Engineering , 91:87–103, December 2016. ISSN 02677261. doi:10.1016/j.soildyn.2016.07.028.David P. Teague, Brady R. Cox, Brendon Bradley, and Liam Wotherspoon. Development of Deep Shear WaveVelocity Profiles with Estimates of Uncertainty in the Complex Interbedded Geology of Christchurch, NewZealand.
Earthquake Spectra , 34(2):639–672, May 2018a. ISSN 8755-2930. doi: 10.1193/041117EQS069M.David P. Teague, Brady R. Cox, and Ellen M. Rathje. Measured vs. predicted site response at the Gar-ner Valley Downhole Array considering shear wave velocity uncertainty from borehole and surface wavemethods.
Soil Dynamics and Earthquake Engineering , 113:339–355, October 2018b. ISSN 02677261. doi:10.1016/j.soildyn.2018.05.031.Gabriel Toro, R. Probabilistic models of site velocity profiles for generic and site-specific ground-motionamplification studies. Technical report, November 1995. Published as Appendix C of Silva et al 1996.Joseph Vantassel. jpvantassel/swprepost: latest (Concept), May 2020. URL http://doi.org/10.5281/zenodo.3839998 .Joseph Vantassel, Brady Cox, Liam Wotherspoon, and Andrew Stolte. Mapping Depth to Bedrock, ShearStiffness, and Fundamental Site Period at CentrePort, Wellington, Using SurfaceWave Methods: Implica-tions for Local Seismic Site Amplification Mapping Depth to Bedrock, Shear Stiffness, and FundamentalSite Period at CentrePort.
Bulletin of the Seismological Society of America , 108(3B):1709–1721, May2018. ISSN 0037-1106. doi: 10.1785/0120170287. 24oseph P. Vantassel and Brady R. Cox. Surface Wave Inversion Benchmarks.
DesignSafe-CI [publisher] ,May 2020a. doi: 10.17603/ds2-cpmr-v194. Type: dataset.Joseph P. Vantassel and Brady R. Cox. SWinvert: A workflow for rigorous surface wave inversions. 2020b.M. Wathelet. An improved neighborhood algorithm: Parameter conditions and dynamic scaling.
GeophysicalResearch Letters , 35(9):L09301, May 2008. ISSN 0094-8276. doi: 10.1029/2008GL033256.M. Wathelet, D. Jongmans, and M. Ohrnberger. Surface-wave inversion using a direct search algorithm andits application to ambient vibration measurements.
Near Surface Geophysics , 2(4):211–221, November2004. ISSN 15694445. doi: 10.3997/1873-0604.2004018.M Wathelet, B Guillier, P Roux, C Cornou, and M Ohrnberger. Rayleigh wave three-component beam-forming: signed ellipticity assessment from high-resolution frequency-wavenumber processing of ambientvibration arrays.
Geophysical Journal International , 215(1):507–523, October 2018. ISSN 0956-540X,1365-246X. doi: 10.1093/gji/ggy286.Marc Wathelet.
Array recordings of ambient vibrations: surface wave inversion . PhD thesis, Universite deLiege Faculte des Sciences Aplliquees, Liege, Belgium, February 2005.Marc Wathelet, Jean-Luc Chatelain, Ccile Cornou, Giuseppe Di Giulio, Bertrand Guillier, Matthias Ohrn-berger, and Alexandros Savvaidis. Geopsy: A User-Friendly Open-Source Tool Set for Ambient Vi-bration Processing.
Seismological Research Letters , April 2020. ISSN 0895-0695, 1938-2057. doi:10.1785/0220190360.Clinton M. Wood, Brady R. Cox, Russell A. Green, Liam M. Wotherspoon, Brendon A. Bradley, andMisko Cubrinovski. Vs-Based Evaluation of Select Liquefaction Case Histories from the 20102011 Canter-bury Earthquake Sequence.