Introduction

Accurate seismic velocity models are essential for making realistic synthetic seismograms that are appropriate for predicting ground-motion parameters and waveforms for future large earthquakes. In particular, accurate three-dimensional (3D) models of deep sedimentary basins are critical for assessing the seismic hazard to tall buildings in those basins (for example, Marafi and others, 2019).

The purpose of this study was to determine if a deployment of broadband seismometers for only 3 weeks in the Seattle basin, Washington, was sufficient to produce noise correlations that could be used to constrain the shallow shear-wave velocities (V_S) in the Seattle basin. We found that noise recordings for a duration as short as 9 days were indeed useful for determining the V_S for depths less than 2.0 kilometers (km).

Computer simulations of ground shaking in future large earthquakes rely on an accurate 3D model of seismic velocities and densities in the crust and upper mantle. Extensive 3D simulations were recently completed for moment magnitude (M) 9 earthquakes on the Cascadia subduction zone (Frankel and others, 2018; Wirth and others, 2018). Previously, 3D simulations for the Seattle basin were conducted for the M6.8 2001 Nisqually earthquake, other recorded local earthquakes, and hypothetical M6.5–7 earthquakes on the Seattle Fault (Frankel and Stephenson, 2000; Frankel and others, 2007, 2009; Wirth and others, 2019; Thompson and others, 2020).

Ground motions from the simulations exhibited substantial amplification of the Seattle and Tacoma basins, similar to that observed in recordings of local earthquakes, for periods of 1 to 10 seconds (s) (Frankel and others, 2018). It is especially important to have accurate models of these sedimentary basins at depths of 10 km and shallower to accurately predict their basin amplification (Frankel and others, 2018).

Stephenson and others (2017) constructed a 3D model for the Cascadia region used in the M9 simulations (denoted herein as SAO17) from a variety of seismological, geological, and geophysical data. The model for the Seattle basin consists of as much as 1.0 km of Quaternary sediments overlying as much as 6 km of Tertiary sedimentary rock, above volcanic basement rocks. The SAO17 model used the map of the thickness of the Quaternary sediments from Johnson and others (1999), which was determined from marine seismic reflection data and the Jones (1996) map of sediment thickness derived from limited borehole data.

The V_S in the Quaternary sediments of the Seattle basin in the velocity model are not well constrained by previous data for depths greater than about 30 meters (m). The SAO17 model has the same V_S with depth within the Quaternary sediments for all locations in the Seattle basin. The V_S profiles from the surface to 30-m depth at several locations in the Seattle basin have been determined from seismic refraction studies (Williams and others, 1999; Stephenson, 2007) and surface-wave dispersion (Wong and others, 2011). The V_S to 100-m depth were measured in two boreholes (Odum and others, 2004). As stated in Stephenson and others (2017), primary-wave velocities (V_P) within the deep sediments of the Seattle basin were constrained by observations from marine seismic data in Puget Sound (Calvert and others, 2003). The SAO17 model was also informed by V_S values for depths greater than 1 km determined from active seismic sources by Snelson and others (2007) for an east-west line in the Seattle basin about 2 km north of the study area of this paper. The SAO17 model used a V_P/V_S ratio that ranged from 2.5 to 2.2 from 0- to 1-km depth, respectively, to estimate the V_S in the Quaternary sediments from the V_P values.

Over the past 15 years or so, correlation of seismic noise between sites has been applied to determine V_S structure on regional and local scales. This approach is based on the idea that the correlation of seismic noise between two sites represents the Green’s function or impulse response between the two sites (Lobkis and Weaver, 2001; Shapiro and Campillo, 2004). Delorey and Vidale (2011) used noise correlation and seismic tomography to determine the V_S of the Seattle basin at depths less than 3.5 km. That study used a recording duration of a few months and found somewhat lower V_S at depths of 1 to 3 km compared to the Stephenson (2007) model. However, the Delorey and Vidale (2011) study had limited resolution of V_S at depths less than 1 km, and the tomographic inversion smoothed the V_S contrast across the southern edge of the basin that is inferred from the surficial geology and observations of basin-edge generated surface waves (Frankel and others, 2002). More recently, Stephenson and others (2019) applied the spatially averaged coherency (SPAC) method to determine V_S profiles to about 1.5-km depth in the Seattle basin. This method is based on fitting the shape of the coherency of seismic noise across a seismic array as a function of wavenumber. They found that increases in V_S corresponding to the top of bedrock occurred at depths roughly comparable to those in the SAO17 model.

Empirical ground-motion models (GMM) that predict ground-shaking levels as a function of magnitude and distance often parameterize the amplification from deep sedimentary basins using the depths where the V_S reaches 1.0 kilometer per second (km/s) (Z_1.0) (for example, Chiou and Youngs, 2014) or 2.5 km/s (Z_2.5). Z_1.0 is generally assumed to be the depth to bedrock (Chiou and Youngs, 2014) and Z_2.5 is taken to be the depth to crystalline basement rocks (Campbell and Bozorgnia, 2014). Thus, the depth to the top of bedrock is important for the prediction of ground shaking using empirical GMMs, as well as to 3D ground motion simulations.

In this study, useable Green’s functions were determined from noise correlations of stations separated by about 5 km after only 9 days of recording. These Green’s functions were analyzed to determine the group velocity of Love waves as a function of period from 0.5 to 6 s. We found that the SAO17 model, with a depth to bedrock of 0.9 km, was generally consistent with these group velocities. We developed an inversion program to find other V_S profiles that matched the observed group velocities. These other models had depths to bedrock ranging from 0.9 to 1.3 km. Thus, the seismic noise data provided useful constraints on the near-surface V_S and Z_1.0 in the Seattle basin.

Data and Cross-Correlation Procedure

Three dense arrays totaling 20 broadband seismometers with digital recorders were deployed in Seattle, Wash., by personnel from the Pacific Northwest Seismic Network (University of Washington) and the U.S. Geological Survey (USGS). The instruments were operational for about 3 weeks in April and May 2017. Figure 1 is a map of the deployment with the arrays labeled NW (northwest), NE (northeast), and SE (southeast). Two of the arrays were located in the Seattle basin (NE, NW) and one (SE) was just to the south of the southern edge of the basin, as defined by geology and gravity (Blakely and others, 2002). The Seattle fault zone forms the southern boundary of the Seattle basin. Within each array, the stations had spacings of 0.5 to 1 km. The NE and SE arrays used compact Trillium sensors and the NW array consisted of Trillium 40 sensors. We removed instrument response before correlating recordings from different types of sensors. The original sampling interval was 0.01 s.

Cross correlations were calculated in the frequency domain for a set of consecutive time windows. We first decimated the seismograms to 10 samples per second. The horizontal seismograms were rotated into the transverse and radial components, relative to the path between the two stations. We initially tried one-bit normalization of the seismograms before cross correlation (see Bensen and others, 2007; Lin and others, 2008) but found it yielded similar cross correlations as the original non-normalized seismograms. The results shown here are from the original seismograms without normalization. We also tried pre-whitening the spectra of the waveforms but determined that this did not improve the results. For the inter-array correlations, we bandpass filtered the waveforms at 0.05 to 1 hertz (Hz) before cross correlation. For the intra-array correlations, a bandpass filter from 0.1 to 4 Hz was applied before correlating. Three stations were not used in the analysis because of a timing problem or intermittent recording (unlabeled symbols in fig. 1). Seismograms from 17 stations were analyzed in this study.

Cross-Correlation Results

Figure 2 shows that a 10-hour window of noise produces a correlation time series for two stations (RAD2 and AB09) 5 km apart with discernable pulses at 5- to 15-s positive and negative lag times. This is representative of the correlations of the transverse component between stations in the NW and NE arrays with paths entirely within the Seattle basin (fig. 1). We found consistent delay times of the peak of the correlated pulse for 10-hour windows over the duration of the deployment. The time series derived from stacking the correlation time series from 20 consecutive 10-hour windows (about 9 days) more clearly displays the pulses at 5- to 15-s positive and negative lag times (fig. 2). Positive lag time corresponds to waves propagating from northwest to southeast and negative lag time corresponds to waves traveling in the opposite direction. The results in this paper are based on stacking 20 consecutive 10-hour windows. Stacking additional windows did not substantially change the results. We did not include the period in the later part of the deployment when an M_L 3.3 earthquake occurred on May 3, 2017, near Bremerton, Wash., about 25 km southwest of the NW array, so that signals from this earthquake did not affect the cross correlations. In general, we found different amplitudes for the pulses at positive and negative lag times. For any given cross correlation, we only used the pulse (positive or negative lag) with the largest amplitude for the dispersion analysis. We typically observed larger correlation pulses for paths travelling west to east across the basin, consistent with the seismic noise microseisms being preferentially generated to the west near the coast. Using the smaller amplitude pulses would not make a substantial difference to the results.

Figure 3 documents the moveout of the correlation pulse across stations in the NE array, indicating that the noise correlation represents the Green’s function of the path. These correlations use the same station in the NW array and three stations in the NE array. As the inter-station distance increases by about 1 km, the peak of the cross correlation is delayed by about 2 s. We conclude that these pulses on the correlation time series of the transverse time series represent the Love-wave impulse response between the two stations.

Cross correlation of the transverse components produced more distinct pulses compared to the correlations of the vertical and radial components. Figure 4 shows the clear pulses on the transverse correlation at 4- to 10-s positive lag and 3- to 12-s negative lag. Identifying pulses on the radial and vertical components was generally problematic for the waveforms we examined. This indicates that Love waves had higher correlations between sites than Rayleigh waves at the same period. Therefore, we restricted our analysis in this paper to the transverse component.

We determined group velocities of Love waves by bandpass filtering each cross-correlation waveform using zero-phase Butterworth filters, enveloping the filtered waveform, and picking the time of the peak amplitude of the envelope (Herrmann, 1973). Figure 5 exhibits the results of bandpass filtering the noise correlation waveform at 3- and 5-s periods, for two sites separated by about 5 km (RAD2 and AB08). The dispersion of the cross-correlation pulses is clearly indicated by the earlier arrival of the 5-s pulse compared to the 3-s pulse. Group velocities were determined by dividing the distance between the stations by the time delay of the peak of the envelope of the bandpass-filtered cross-correlation waveform.

There were several challenges in the cross-correlation analysis. For periods shorter than 2 s, there were generally no distinct pulses in the cross correlation for inter-array paths of about 5 km. At periods longer than about 6 s, the peak of the filtered waveform for the inter-array paths was too close to zero lag to be resolved in time. Theoretical studies recommend that the path be long enough such that the travel time is longer than the period (Tsai, 2009). The ratio of the delay time and maximum period of the waveforms in this study is comparable to that of the laboratory ultrasound study of Lobkis and Weaver (2001), who initially proposed the use of ambient noise correlation to determine Green’s functions.

Bandpass filtering the cross-correlation time series around a center period does not always produce a waveform that is predominantly at that center period. We chose center periods of 0.5, 1.0, 2.0, 3.0, 4.0, 5.0, and 6.0 s. The frequency width of the bandpass filters was 10 percent of the center frequency. The band-limited nature of the cross correlation caused by the dominant period of microseisms can make the dominant period of the filtered waveform lower than the center period specified in the filtering (Bensen and others, 2007). We calculated the dominant period of the filtered waveforms by measuring the time between zero crossings near the peak amplitude.

We calculated group velocities for 12 paths between the NW and NE arrays from 3 to 6 s periods and for 4 paths at 2 s (fig. 6). The station pairs have separation distances of about 5 km and are listed in table 1 and shown in figure 1. The average group velocity values for these inter-array measurements are plotted in figure 6 at the average periods of the bandpass-filtered cross-correlation waveforms. The error bars in figure 6 denote one standard deviation of the measurements for different pairs of stations. The other pairs of stations we examined yielded poor cross correlations and were not used in this study.

We also found usable cross correlations at 0.5- and 1.0-s periods for sites separated by about 1 to 2 km within the NE array. The group velocities at 0.5- and 1.0-s periods for five pairs of NE stations (intra-array paths) are also plotted in figure 6. We did not find usable cross correlations at these periods for stations within the NW and SE arrays.

The group velocities increase from about 0.42 km/s at a 0.5-s period to 0.72 km/s at about a 6-s period (fig. 6). In general, the standard deviations were larger for the shorter periods. The standard deviation in group velocities for a given period is likely caused by lateral variations in the V_S depth profiles in the basin, as well as by measurement uncertainties.

We also cross correlated six paths crossing the southern edge of the Seattle basin (see table 1; fig. 1). These correlations were between stations in the NW and SE arrays and between stations in the NE and SE arrays. Figure 6 shows the average results found at 4- to 5-s periods. Shorter periods did not exhibit useful correlations. The group velocities for paths crossing the basin edge were about 1 km/s faster than those for the paths within the basin. This is expected because of the higher V_S in the sedimentary rocks south of the basin edge compared to the glacial sediments in the basin at depths less than 1 km. We did not find useful cross correlations for inter-array paths with endpoints at AB11 or AB03 (fig. 1). These stations were located in an area of shallow bedrock (Tertiary sedimentary rock) which may have different noise characteristics. Cross correlations could be used for paths entirely outside the Seattle basin to quantify the V_S in the rocks south of the basin.

Table 1.    

Station pairs used to determine group velocities.

[s, seconds]

Table 1.    Station pairs used to determine group velocities.

Paths and periods	Station pairs
Inter-array basin paths, 3- to 6-s periods	RAD2–AB08, RAD2–AB09, RAD2–AB10, RAD4–AB04, RAD4–AB09, RAD4–AB10, RAD4–AB14, RAD4–AB15, RAD6–AB04, RAD3–AB13, RAD3–AB16, RAD3–AB17
Inter-array basin paths, 2-s period	RAD2–AB08, RAD2–AB09, RAD4–AB14, RAD4–AB09
Intra-array, 0.5- to 1.0-s periods	AB14–AB04, AB16–AB10, AB17–AB14, AB08–AB09, AB04–AB10
Inter-array, crossing edge of basin, 4- to 5-s periods	AB10–AB05, AB10–AB06, AB04–AB06, AB16–AB06, RAD4–AB06, RAD2–AB06

Analysis of Group Velocities to Determine Shear-Wave Velocity Profiles

We compared the group velocities derived from the noise correlation for paths within the basin with those calculated for a flat-layered velocity model based on the SAO17 model, as well as other profiles. Love-wave group velocities for flat-layered velocity models were determined using the computer programs of Herrmann (2013). The layer thicknesses used correspond to those in the SAO17 model: 0.05-km thickness of the top layer, 0.1-km thickness to 1.3-km depth, and 0.3-km thickness from 1.3- to 10-km depth. The SAO17-model V_S profile has a jump in V_S at the base of the Quaternary sediments at 0.9-km depth (fig. 7). Tertiary sedimentary rock underlies the sediments. There is another distinct jump in V_S at the top of the volcanic rock (basement) at about 7-km depth. The original SAO17 model had a minimum V_S of 600 meters per second (m/s) for the top layer, based on the limitations of the 3D simulations. For the group velocity calculations, we lowered the V_S in the top layer of the SAO17 model to 0.45 km/s, to be consistent with observations of near-surface V_S on the glacial soils in the Seattle basin (Williams and others, 1999). We calculated densities from the V_P in the SAO17 model using the empirical relationship given in Brocher (2005).

The SAO17 model at the midpoint of the arrays (MP in fig. 1) fits the observed group velocities except for at a 1-s period, where it was slightly larger than one standard deviation above the mean observed value (fig. 6). We found that introducing a V_S gradient from 0.05- to 0.1-km depth improved the fit to the group velocity at a 1-s period. Thus, a minor modification of the SAO17 model matched the observed group velocities.

A key issue is the non-uniqueness of velocity profiles that fit the observed group velocities. To assess this non-uniqueness, we wrote a computer program to perform a nonlinear inversion of the observed group velocities to determine the V_S for depths less than 2 km. Each iteration of the inversion determines the adjustments to the input V_S profile that improve the fit to the observed group velocities. The program first finds the change in group velocity at each period given a perturbation in the V_S for each layer. The group velocity at each period is calculated from the computer codes of Herrmann (2013). Damping and smoothing is also applied to stabilize the inversion. The inversion uses singular value decomposition to invert the matrix containing the derivatives of group velocity with respect to V_S in each layer. The inversion program iterates until the root-mean-square error between the observed and predicted group velocities is no longer reduced significantly. The velocity in the top layer is constrained to be close to 0.45 km/s, based on the observed group velocity at a 0.5-s period.

Various V_S profiles were found that adequately fit the observed group velocities from 0.5- to 6-s periods. Our suite of starting models had a constant velocity in the layers that were inverted for V_S. Figure 7 contains the V_S profiles from the SAO17 model and those from the inversions that fit the observed group velocities. The group velocities derived from the inversions are shown in figure 6. The starting V_S profiles and the profiles determined from the inversions are displayed in figure 8. All of these profiles adequately fit the group velocity observations (fig. 6), with values within one standard deviation of the mean at each period.

For the first starting model, we solved for the V_S in 14 layers at depths less than 2.0 km (fig. 8A). The layers have a thickness of 0.1 km to 1.3-km depth and 0.3 km for depths greater than 1.3 km. The top layer has a thickness of 0.05 km. The starting model was identical to the SAO17 model at depths greater than 2.0 km and had a V_S of 1.38 km/s at depths less than 2.2 km. This inversion yielded a velocity model with a jump in V_S at 1.3-km depth, which we assume represents the top of bedrock. The velocity gradient at depths less than 1.3 km was greater than in the SAO17 model.

We also tried a starting model using the V_S values of the SAO17 model for depths greater than 0.9 km and a V_S of 1.2 km/s at depths less than 0.9 km (fig. 8B). The inversion solved for V_S at depths less than 0.9 km. This resulted in a model that also fit the data well. The inversion produced a model with lower velocities than 1.0 km/s at depths less than 0.9 km and also a small velocity decrease with increasing depth from 0.5 to 0.9 km. Therefore, we interpret the depth to bedrock (Z_1.0) in this model to be 0.9 km.

We found that V_S profiles from the inversions with Z_1.0 values of 0.5 and 0.7 km did not adequately fit the observed group velocities. The starting models and inversion results are shown in figure 9. We tried a starting model with V_S of 1.2 km/s in depths less than 0.5 km. In this case, the inversion solved for V_S in depths less than 0.5 km to see if models with a Z_1.0 of 0.5 km could fit the observations. This model could not produce the observed group velocities (fig. 6). A starting model with a V_S of 1.2 km/s in depths less than 0.7 km (fig. 9) yielded an inversion result with a poor fit to the observed group velocity at about a 6 s period (fig. 6). This indicates that models with a Z_1.0 of 0.7 km do not fit the observed dispersion.

We conclude that depths to top of bedrock of 0.7 km or less are inconsistent with the observed group velocities. Thus, the noise observations provide an important constraint on the velocity model.

We also tried two starting models that followed the SAO17 model deeper than 2.0 km but with V_S perturbations of ±0.2 km/s. We consider this to be a reasonable representation of the uncertainty of V_S at depths greater than 2.0 km. The starting velocities at depths less than 2.0 km were 1.18 and 1.58 km/s (fig. 8C, D). The inversions found two models that fit the observed group velocities. The faster of the starting models produced a V_S jump at 1.3-km depth in the inversion result. This jump occurred at the same depth as in the SAO17 model but had a larger increase. The slower starting model produced a final model without a distinct jump in V_S at any given depth less than 2.0 km. However, the Z_1.0 determined from the inversions using both of these starting models was 1.3 km (fig. 8C, D).

The velocity profiles from the surface to 0.5-km depth were quite similar between all of the V_S profiles that fit the observed group velocities (fig. 7). The velocity profiles diverged for depths of 0.5 to 0.9 km, indicating that the observed group velocities provide less resolution of the V_S in this depth range.

Better constraints on the V_S at depths greater than 2 km could be achieved with longer period group velocities. Determining these group velocities would require greater station separation, which would be difficult to achieve for sites in the Seattle area without encountering lateral variations in V_S.

The results of this report are consistent with the results of Stephenson and others (2019). Using SPAC analysis of seismic noise recorded on triangular arrays, they found an increase of V_S at 0.8±0.2 km depth for sites in the Seattle basin. This depth range is comparable to the 0.9- to 1.3-km depth for the jump in V_S determined from cross correlation of seismic noise from the wider spaced arrays used here.

Summary

This study demonstrated that the correlation of seismic noise with a short deployment time of 9 days is useful for evaluating the shear-wave velocity structure at depths to 2.4 km in the Seattle basin. We found several velocity profiles with depths to bedrock of 0.9 to 1.3 km that matched the observed group velocities. The results indicate that these noise data provide important information on the near-surface shear-wave velocity profiles and the depth to bedrock (Z_1.0). These parameters are important in three-dimensional simulations, as well as empirical ground-motion models. Thus, relatively short deployments of instruments provide noise correlations that constrain the depth to bedrock in the Seattle basin, which is important for estimating ground motions from future earthquakes and, therefore, seismic hazard.