Research ArticleGEOPHYSICS

Revisiting evidence for widespread seismicity in the upper mantle under Los Angeles

See allHide authors and affiliations

Science Advances  22 Jan 2021:
Vol. 7, no. 4, eabf2862
DOI: 10.1126/sciadv.abf2862


We revisit the finding of widespread deep seismicity in the upper mantle imaged with a dense, temporary nodal seismic array in Long Beach, California using back-projection to detect candidate events and trace randomization to develop a reliable imaging threshold for candidate detections. We find that nearly all detections of small events at depths greater than 20 kilometers in the upper mantle fall below the reliability threshold. We find a modest number of small, shallower events in the crust that appear to align with the active Newport-Inglewood Fault. These events occur primarily at 15- to 20-kilometer depth near the base of the seismogenic zone. Localized seismicity under fault zones suggests that the deep extensions of active faults are localized and deforming, with stress concentration leading to a concentration of small events, near the seismic-aseismic transition.


Metropolitan Los Angeles is located within the actively deforming Pacific-North America plate boundary. The Newport-Inglewood Fault Zone (NIFZ) is among the most prominent active faults within and/or bounding the Los Angeles Basin. It is mapped through the Los Angeles Basin from a northwestern terminus near the Santa Monica–Hollywood Fault to Newport Beach (Fig. 1A), where it continues offshore. The NIFZ is an active right-lateral strike-slip fault system, with at least five ML (local magnitude) ≥ 4.9 earthquakes since 1920, including the 1933 Long Beach earthquake (ML = 6.3) (1), which is the second deadliest earthquake in California history.

Fig. 1 Map of Los Angeles Basin and Long Beach phase B deployment.

(A) Map of Los Angeles Basin showing the Newport-Inglewood Fault and other faults (LB, Long Beach; VF, Verdugo Fault; ERF, Eagle Rock Fault; EMF, East Montebello Fault; WHF, Workman Hill Fault). Red star shows epicenter of the 1933 Long Beach earthquake. Blue and green polygons outline the Long Beach phase A and B deployments. Blue dot is epicenter of a local M 2.1 earthquake that occurred during the deployment of phase B. OO′ is a profile across the epicenter. (B) Map of Long Beach phase B deployment. Stripe of missing sensors in the upper left is the Long Beach Airport runway. Narrow gap in the northern and eastern part of the deployment tracks the I-405 Freeway. Black lines show mapped surface trace of the Newport-Inglewood Fault. Red dashed rectangle is the surficial boundary for the 3D imaging cube we use in our analysis.

The shallow part of the NIFZ has been resolved to some extent (24), but its deep structure is less clear. Microearthquake monitoring is the primary tool for delineating the deep geometry of active faults (57). Microseismicity also provides crucial information on deformation processes, including potentially earthquake nucleation (8). Detection of small earthquakes, however, is hampered by the sparsity of seismic networks and high-level seismic noise in this setting. The level of small earthquake signals barely exceeds the low noise floor at Earth’s surface on individual seismograms, and in urban settings, such as that for the NIFZ, the noise floor is high. As a result, traditional earthquake detection methods based on single or a few station measurements may fail to detect small microseismic events. Template matching (9), data mining (10), and machine learning (11) methods have been developed for this purpose; however, in the presence of strong noise, these methods can be susceptible to a combination of lower sensitivity and more false detections.

Inbal et al. (12) applied a novel approach to detect deep microseismic events along the NIFZ using waveforms recorded by the Long Beach nodal deployment with dense, ~100 m, spacing, in urban settings. To suppress strong cultural noise, they performed downward continuation to back propagate the wave field from Earth’s surface to 5-km depth and then applied back-projection (BP) to image/locate seismicity at depth. They found widespread seismicity in the lowermost crust and upper mantle at depths greater than 20 km. These results paired with geochemical evidence from helium isotopes along the NIFZ (13) led Inbal et al. (14) to conclude that the NIFZ is seismically active down to the upper mantle. They proposed a conceptual framework of deep deformation that involves a combination of aseismic ductile flow with interspersed seismogenic deformation.

Inbal et al.’s finding (12, 14) challenges the widely held notion that seismicity on continental faults is confined to the crust, with deformation in the lower crust and upper mantle accommodated aseismically (1517). Li et al. (18) followed this work with a new detection approach that exploits local waveform similarity and applied it to the Long Beach dense array data to detect very weak microseismicity. They analyzed 1 week of data and found that most of the earthquakes occurred in the shallow crust. They suggested that the difference is because Inbal et al. used downward continuation to enhance earthquake detection at greater depths and suppress earthquake detection above 5 km. They were measured in their interpretation, however, because they only performed their grid search over horizontal space and not over depth. These seemingly contrary results on the seismogenic depth of the high-risk active fault zone led us to reexamine the seismic location/imaging results. Here, we study the spatial distribution of earthquakes imaged based with Long Beach dense array data, demonstrate that most of the candidate upper-mantle events are likely artifacts, and identify a subset of shallower, small events that are more likely to be real events with their own characteristic distribution pattern.

Challenges with the Long Beach dataset

The data were collected by Signal Hill Petroleum, initially for the purpose of oil and gas exploration, but because of its deployment straddling important faults in Southern California and unusually high spatial density, some of the data have been openly shared for seismological research. The Long Beach dense array with ~100-m spacing was deployed for two phases, in urban settings (Fig. 1A, phase A and phase B). Phase A covered an area of 10 km by 7 km and included approximately 5200 vertical 10-Hz geophones with sampling frequency at 500 Hz. It operated for the first 6 months of 2011. Phase B extended the original survey toward the east, using the identical sensors and sharing similar tectonic setting above a branch of the NIFZ. It covered an area of 8.5 km by 4.5 km with approximately 2500 sensors and operated for the first 3 months of 2012. We work on the phase B dataset, the volume of which exceeds 50 terabytes. Figure 1B shows the map of Long Beach phase B deployment. With time-frequency domain analysis of sections of seismograms, we find a variety of cultural noise sources (traffic, airplane, and helicopter) within a 220-s seismogram (fig. S1A), which are similar to those in Meng and Ben-Zion (19). We show other unidentified transient signals in fig. S1 (B and C). Figure S1 conveys some of the challenges with working at high noise levels. The conventional approach of detecting phases on individual seismograms, phase association with events, and event location is infeasible for weak signals in the presence of strong nonstationary noise from such sources.

Another factor that limits earthquake detection is the narrow aperture of the Long Beach deployment with respect to the target depth of microseismic sources. We perform synthetic tests to assess the imaging resolution of the BP method (Materials and Methods and fig. S2). For shallow sources, we find that BP has good resolution in both horizontal and vertical directions; however, for deep sources, the imaging resolution deteriorates, particularly in the vertical direction, due to small variability of travel time difference at these depths (20). With an approximate 10-km aperture for the detection at 24-km depth (fig. S2C), we cannot easily differentiate between sources ranging from 22- to 28-km depth. The real data will be yet more challenging, with stronger and more varied noise sources, inaccurate velocity model, etc. Thus, with the limited available aperture, we should be able to detect sources, but we expect it to be challenging to localize the depths, for events greater than 20 km.


Evaluation of the reliability in earthquake detections

We evaluate the significance of a detection by comparing the strength of BP images in the original versus trace-randomized Long Beach dense array data for the validation of earthquake location (Materials and Methods and fig. S3). Trace randomization (TR) involves randomly reassigning the seismic recordings to receivers at different locations rather than the original receiver location. We would not expect a true image to form, but the signals accurately represent the noise levels, bandwidth of the original signals, and geometry of the array. Maxima in our image space for trace-randomized BPs provide a lower bound on the strength of detections required for confidence that the results represent true earthquake detections rather than artifacts due to random and fortuitous alignment of local noise sources. Only if the maximum amplitude of the imaged source location is higher than that from randomized waveforms are the original waveforms likely to come from a real source.

Inbal et al. (12, 14) analyzed phase A of the Long Beach array deployment. We analyzed the adjacent phase B. Given the immediate proximity of these datasets and the massive number of detections reported previously, we expect that our results should generalize to phase A. We image earthquake sources below the Long Beach phase B deployment following the data processing workflow shown in Materials and Methods and fig. S4. We do not perform downward continuation, which is equivalent to back propagating the time-reversed seismic recording from surface and recording the new wave field at another depth (21). Figure S5 shows that after downward continuation, the signal-to-noise ratio (SNR) for traces near the source increases while that for traces farther away from the source decreases. Downward continuation does not improve the BP imaging quality; however, it suppresses the detection of near-surface sources.

Figure 2 (A and B) show the BP imaging results of Julian day 28, 2012 (during local night between 10 p.m. and 4 a.m.), based on original data and trace-randomized data, respectively. We adopt the same detection threshold value for trace-randomized data as that for original data. In the three-dimensional (3D) imaging cube, each dot represents one detection, with its color indicating the depth. From Fig. 2 (A and B), we find similar patterns of the earthquake depth distribution before and after TR. The statistics of the earthquakes at different depths and the BP energy before and after data TR is shown in Fig. 2D. In summary, we find the following:

Fig. 2 Earthquake imaging/location and statistical results.

(A) and (B) are the BP imaging results of Julian day 28, 2012 night time (local) on original and trace-randomized data, respectively; (C) is the BP relocation result on the same data in (A) using the new detection criterion; (D) and (E) are the statistical results for detections on 1 and 10 days, respectively. Energy in (D) is calculated by summing up the imaging amplitudes of all the detections.

1) Within the range of 0- to 20-km depth, the number of detections in both the original and trace-randomized data decreases with increasing depth. The BP imaging results with only random noise show a similar pattern (fig. S6);

2) At 20 to 25 km, the deepest 5-km interval, the number of detections in the original data is greater than in intervals of 5 to 10 km, 10 to 15 km, and 15 to 20 km. This is likely because the incident energy from the local/regional events outside the 3D imaging cube can result in detections at the boundary of the cube as we demonstrate below;

3) After TR, the number of the detections decreases by 23.73%, while the BP energy for these events decreases by 27.35%. The fact that most of the detection and energy remains before and after TR indicates that most of the detections could be artifacts rather than real earthquake sources.

The results for 10 days (Julian days 27 to 37, except 29) shown in Fig. 2E are similar to the 1-day result. The greater drop in detections at 0 to 5 km than those at greater depths after TR is probably caused by a combination of anthropogenic sources at the surface and better detection sensitivity at shallower depth.

Absence of upper mantle seismicity under new detection criterion

In addition to assessing the reliability of the prior detections, TR provides a useful criterion for the event detection from data with nonstationary noise. That is, only if the maximum amplitude for a potential detection decreases by an amount larger than the detection threshold after TR, can it be confirmed as a reliable detection. We apply this criterion to image microseismic events. Because the SNR for most microseismicity is low, we repeat BP on multiple sets of trace-randomized data to obtain a stable result. We remove from further consideration those detections that locate at the cube boundary from the earthquake relocation results (see Discussion).

Figure 2C shows the microseismic relocation result for Julian day 28, 2012 under the new criterion. We find a large reduction in the number of detections. Figure 3 (A and B) shows the slices for different depth ranges. We do not find evidence for widespread, deep earthquakes below 20 km as reported by Inbal et al. (12, 14). Although no clear earthquake trend is observed from the near surface to 15-km depth, the detections in the 15- to 20-km-depth range appear to correlate somewhat with the surface fault traces. The densest detection area (roughly x: 1293.8 to 1294.4 km; y: 1230 to 1231 km) matches the fault junction. When we slightly increase the detection threshold, the remaining dots in the depth slice of 0 to 5 km match the fault trace (Fig. 3B), and the missing dots in the upper-right part of the panel follow the traces of the freeway at the surface (narrow gaps within the deployment sensors in Fig. 1B). Horizontal slices for 10 days of earthquake detections are shown in fig. S7.

Fig. 3 Horizontal slices for earthquake detection on Julian day 28, 2012 (night time) under the new criterion.

(A) Horizontal slices for different depth ranges of 0 to 5 km, 5 to 10 km, 10 to 15 km, 15 to 20 km, and 20 to 25 km, respectively. (B) Horizontal slice for 0- to 5-km-depth range under a slightly higher detection threshold.


Our new detection criterion based on TR is effective in reducing false detections. Without this new criterion, the earthquake location result in Fig. 2A is under the same detection threshold as that in Inbal et al. (12), which is five times the median absolute deviation (MAD) of the noise distribution. Inbal et al. set this detection threshold based on the assumption of a Gaussian noise distribution. They produced synthetic tests by generating 1000 realizations (windows) of Gaussian noise with variance equal to the variance of the BP images, fitting a Gumbel distribution to the maximum value of each realization, and calculating the probability of exceeding 5-MAD at 3.22 × 10−7. We compare the Gumbel distribution to the real data by fitting the distribution of peak amplitude values from all time windows at each imaging grid point. We find that the real data do not follow the Gumbel distribution (fig. S8A) but are better described by a generalized extreme value (GEV) distribution with c = 0.168527 (fig. S8B). Assuming that the maxima of noise data follow the GEV distribution, we estimate that the probability exceeding 5-MAD is 0.001549, which leads to an expectation of 17 of 10,800 time windows (9 hours). We count a total number of 59 peak values from the time windows exceeding 5-MAD. The number of false detections (17) is nontrivial compared with the total detections (59).

Some false detections in the 3D imaging results by BP method (e.g., imaging result shown in Fig. 2A) may come from local/regional earthquakes originating outside of the 3D cube and localized nonearthquake sources at the surface. We show a local M 2.1 earthquake with hypocenter 2.5 km away from the western boundary of the 3D imaging cube (blue dot in Fig. 1A) accounting for one detection (red circle) at the boundary of the 3D cube in fig. S9A. This can explain the concentration of detections at the bottom and four vertical surfaces that bound the imaging volume, which should be removed for accurate location results. Movie S1 shows how the Long Beach dense array recorded the ground motion from this earthquake. BP imaging results also show an anomalously large number of events clustering at the northwest corner of the surface and far more detections through the entire imaging cube from top to bottom on Julian day 29, 2012, which were from vibroseis sources (fig. S10A). The vibroseis operating around 5.5 hours (fig. S10B) near the northwest corner of the imaging volume led to many apparent detections located far from the causative source (fig. S10A), indicating that localized nonearthquake sources at that surface can result in deep location artifacts in BP imaging results. These false detections can be excluded by the new detection criterion in our relocation results.

From the earthquake relocation results, we do not find widespread seismicity in the upper mantle under Los Angeles. The result showing that small events at 15 to 20 km follow the surficial fault traces (Fig. 3A) suggests a seismogenic model in which the base of the seismogenic zone accommodates more seismic activity than faults at shallower depth (Fig. 4). Assuming that the base of the seismogenic zone extends to 15-km depth and transitions from locked to slipping below that, the shear stress should increase at the boundary, much as stress/strain are focused on stuck sections of otherwise creeping faults (22). These, and other, results from dense-array monitoring suggest that there may be much we can learn from such data, particularly as dense nodal instrumentation and distributed acoustic sensing are more routinely deployed.

Fig. 4 Schematic plot for the strain rate at seismogenic depth.


BP method and synthetic test

The implementation of the BP algorithm mainly consists of time shift and stacking, which can be expressed as (23)Stacki(t)=1NΣk=1NSk(ttik)(1)where Sk(t) is the seismogram recorded at the kth station, tik is the predicted travel time from the ith grid point to the kth station calculated based on a known 3D velocity model, N is the number of stations, and Stacki(t) is the stacked seismogram for the ith grid-searching point.

We represent the imaging volume as a 3D cube and perform grid search on each potential source location within it. For each grid point, the seismograms are time-shifted on the basis of the travel time between the geophones and this point, and the aligned seismograms are stacked to a single time series just for this grid point. We thus obtain a 4D imaging result. At each time step, a grid point will only be confirmed as the real source if it has the largest stacked amplitude and it exceeds the detection threshold.

The imaging quality of the BP method is influenced by the density of deployment, array aperture, accuracy of velocity model, and SNR. The horizontal resolution is largely determined by the spacing between stations, while the vertical resolution depends on the array aperture. To understand better how BP works on Long Beach data, we perform 2D synthetic tests for BP imaging on a single point source. For the analogy with the Long Beach dataset, we design an imaging profile 10 km in horizontal extent and 30 km in depth extent. We use the vz component for imaging, because the Long Beach deployment only has vertical component.

Figure S2 shows the BP result with the earthquake source located at 2, 10, and 24 km, respectively. We observe that for shallow sources, the BP method has good resolution in both horizontal and vertical directions; however, when the source is deeper, the imaging result loses resolution in both directions, and the vertical resolution deteriorates more rapidly. For a source located at 24 km (fig. S2C), we cannot tell the exact depth from 22 to 28 km.

The maximum depth the deployment can resolve reliably is influenced by the array aperture, the velocity of the medium, and the sampling frequency. With fairly high sampling frequency at 500 Hz of Long Beach data, the main controlling factor is the array aperture. If the source is too deep for the array to image, which means the travel time differences from the source to different sensors at the surface are comparable to or smaller than the sampling interval, the array will not resolve the source depth (fig. S2D). The real data are more complicated with ample noise, inaccurate velocity model, etc. We conclude from the synthetic tests that although we should be able to detect deeper events, the limited array aperture < 10 km means that resolving the depth of earthquakes below 20 km is challenging.

BP imaging with original and trace-randomized data

We illustrate the strategy of BP imaging with original and trace-randomized data with fig. S3. The source is located at 10-km depth, and a 10-km-long seismic survey line containing 100 geophones is at the surface. Figure S3A shows the seismic profile with a hyperbolic move out, and the BP imaging result is shown in fig. S3D. Randomizing the traces of the seismic recordings leads to seismic profiles with an irregular alignment (fig. S3B). When we back-project the trace-randomized data, we get an imaging result with widely dispersed amplitude as shown in fig. S3E. The ratio between the maximum amplitude before (fig. S3D) and after (fig. S3E) TR is 1:0.252, indicating energy defocusing due to TR. If we perform TR on nonearthquake seismic recordings (here, we take the incoherent waveforms in fig. S3B as such an example), both the BP imaging results before (fig. S3E) and after another TR (fig. S3F) show dispersed energy and have similar maximum amplitudes (0.252:0.239).

Seismic data processing

To reproduce the imaging result by Inbal et al. (12), we process the seismic data with a similar workflow to theirs outlined in fig. S4. We convert the original data in SEG-D format to binary format, decimate from 500 to 50 Hz, and apply a band-pass filter at 5 to 15 Hz. The original continuous data are segmented to daily data, and we only use data during local night time (9 p.m. to 6 a.m.) for BP. To suppress the influence of strong space-dependent noise level and large noise bursts, we normalize the data with their 15-min maximum value. To reduce the sensitivity to the inaccuracy of the velocity model, we then calculate the envelope by smoothing the data with a three-point median window on the squared waveforms.

We perform BP in a 6 km by 4.4 km by 25 km 3D cube with grid spacing of 200 m in each dimension. The geographic boundary of the imaging cube is shown as red dashed rectangle in Fig. 1B. We calculate the travel time between each grid point and each geophone at the surface based on the Southern California Earthquake Center Community Velocity Model and store them in a travel time table for convenience. Then, for each grid point in the 3D cube, the seismograms are time-shifted and stacked to a single stacked seismogram. The stacked seismogram is segmented into 3-s time windows, and the maximum value within each time window is assigned to the BP value of this grid point. We thus obtain a 3D imaging cube for each 3-s time window. If the maximum BP value for a time window exceeds the detection threshold, the corresponding grid point will be marked as a detection. We set the same detection threshold as that in Inbal et al. (12), which is five times the MAD of the noise distribution.


Supplementary material for this article is available at

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.


Acknowledgments: We thank Signal Hill Petroleum for sharing the 2012 Long Beach 3D Seismic Survey Extension Continuous Data. We appreciate the comments from W. Ellsworth and Y. Ben-Zion that help improve the manuscript. We acknowledge A. Lellouch for useful discussions. We thank A. Inbal for sharing velocity models. Funding: This work is funded by US Geological Survey (USGS) External Grant award No. G20AP00015 and Southern California Earthquake Center (SCEC) funding No. 19101. This work is supported by SCEC contribution 10871. SCEC is funded by the NSF cooperative agreement EAR-1600087 and USGS cooperative agreement G17AC00047. Author contributions: G.C.B. conceived the study. L.Y. analyzed the seismic data and wrote the manuscript. X.L. converted the seismic data format and made statistical analysis. Competing interests: The authors declare that they have no competing interest. Data and materials availability: All results needed to evaluate the conclusions in the paper are provided in the paper and/or the Supplementary Materials. The raw seismic data can be requested from Signal Hill Petroleum.

Stay Connected to Science Advances

Navigate This Article