Anthropogenic Seismicity Rates and Operational Parameters at the Salton Sea Geothermal Field

See allHide authors and affiliations

Science  02 Aug 2013:
Vol. 341, Issue 6145, pp. 543-546
DOI: 10.1126/science.1239213

This article has a correction. Please see:

Hot and Bothered

One of the most uncertain aspects of energy production at geothermal fields is the potential for the induction of earthquakes from extraction of hot fluids and injection of wastewater back into the subsurface. Brodsky and Lajoie (p. 543, published online 11 July) compared the historical rates of seismicity—after correcting for the occurrence of aftershocks—to operational parameters at the Salton Sea Geothermal Field in southern California. The net volume of fluid extracted and injected was correlated with recent seismicity. However, since the 1980s, when the first plant came online, the rate of earthquakes induced with additional increments of volume injected has decreased.


Geothermal power is a growing energy source; however, efforts to increase production are tempered by concern over induced earthquakes. Although increased seismicity commonly accompanies geothermal production, induced earthquake rate cannot currently be forecast on the basis of fluid injection volumes or any other operational parameters. We show that at the Salton Sea Geothermal Field, the total volume of fluid extracted or injected tracks the long-term evolution of seismicity. After correcting for the aftershock rate, the net fluid volume (extracted-injected) provides the best correlation with seismicity in recent years. We model the background earthquake rate with a linear combination of injection and net production rates that allows us to track the secular development of the field as the number of earthquakes per fluid volume injected decreases over time.

The exploitation of geothermal resources is rapidly expanding as society increases its reliance on renewable energy sources. Production of geothermal power induces seismicity as water is pumped both into and out of a reservoir (1, 2). Fluid injection, a major component of most geothermal operations, has been shown to induce seismicity in a variety of settings, and the earthquakes are often attributed to a decrease in the effective stress on faults due to increased pore pressure (37). Earthquakes can also be induced by fluid extraction through more complex processes, such as thermal contraction or subsidence-driven increases in shear stress (8, 9). Seismic consequences of geothermal production are of particular concern for facilities neighboring large tectonic systems. Therefore, it is important to understand the controls on geothermal production–related seismicity in a setting such as the Salton Sea Geothermal Field, which borders the diffuse terminus of the San Andreas Fault system in southern California.

The Salton Trough is home to four operating geothermal fields (Salton Sea, Brawley, Heber, and East Mesa) that generate a total of over 650 MW of electricity (Fig. 1). The Salton Sea Geothermal Field exploits a hot, geothermal brine, with temperatures in excess of 320°C at 2 km depth (10, 11). The field includes 10 operating geothermal power plants. Plants use flash technologies that extract fluid from depth and flash a portion to steam in order to power generators, while the remaining brine is either injected via separate wells back into the reservoir or subjected to further flashing. Steam condensate is also recaptured for injection. The volume of fluid loss during operations depends on a variety of conditions, including salinity in the flashing chambers and air temperature. Since 1992, the average reported monthly injection volume is 81% of the reported production volume, with a standard deviation of 7%. Production at the plant cycles annually in response to demand and local environmental conditions.

Fig. 1 Earthquake and geothermal facility locations and activity.

(A) Regional map with faults and location of the Salton Sea Geothermal Field. (B) Drill year of the wells in the field for 1960–2012. (C) Earthquakes (blue circles) and injection wells (red stars) in map view. (D) East-west cross-sectional view. Earthquake hypocenters cluster around and beneath injection wells. Depth is relatively poorly constrained in the sedimentary basin (27) and is not used for any subsequent statistics in this study.

The first plant came online at the Salton Sea Geothermal Field in 1982 (10). Operations expanded steadily from 1991, with new wells regularly added and fluid volumes increasing through 2005 (Figs. 1B and 2). For the highest production month in 2012, 11.2 × 109 kg (~107 cubic meters) of geothermal fluid was extracted from the reservoir at depths of ~1 to 2.5 km, and 9.2 × 109 kg was injected at similar depths (12, 13).

Fig. 2 Background seismicity rate (μ) over time compared with operational fluid volumes at the Salton Sea Geothermal Field.

The seismicity rate curve is identical for each panel (right axis, green curve), and the operational rate (left axis, blue curve) in each case is (A) production rate, (B) injection rate, and (C) net production rate. Seismicity rate estimation is on 2-year overlapping intervals centered on each month for which there is operational data. Confidence levels on μ are mapped in fig. S3.

We quantified the relationship between fluid volumes within the Salton Sea Geothermal Field and the local rate of seismicity by combining publicly available data sets. By law, in California the monthly fieldwide total production and injection volumes are released to the California Department of Conservation. For earthquake locations and times, we used the largest high-precision catalog available for the region (14) for January 1981 to June 2011 supplemented by the Southern California Seismic Network catalog for July 2011 to December 2012. Station configurations have changed over time, with a notable increase in data after the release of the geothermal field’s local seismic network in 2007. Therefore, we restricted the data to the local magnitude of completeness of 1.75 (supplementary text) to ensure that we only analyze events that are large enough to have been detectable throughout the study period. For this study, the Salton Sea Geothermal Field seismicity is defined as earthquakes shallower than 15 km in the region bounded by 33.1 to 33.25° N and 115.7 to 115.45° W (Fig. 1).

The seismicity rate has mirrored overall activity in the field (Fig. 1 and fig. S2). As noted by (15), earthquakes cluster around injection wells both at the surface and at depth. The seismicity rate was initially low during the period of low-level geothermal operations before 1986. As operations expanded, so did the seismicity. The maximum net volume production rate occurred in July 2005, which is a month before the largest earthquake rate increase. This August 2005 swarm has been linked to a creep transient but had not been compared previously with the production data (16, 17).

However, the relationship between seismicity and plant operations is not simple. Earthquakes are highly clustered owing to local aftershock sequences, and so it is difficult to untangle the direct influence of human activities from secondary earthquake triggering. In addition, seismicity rate varies over orders of magnitude, whereas pumping conditions evolve more smoothly. Operations are continually changing at the plants in response to both economic and natural factors. These issues complicate the data so that a simple correlation between the raw seismicity rate and operational parameters is suggested but unclear and requires further analysis (fig. S2).

Because earthquakes commonly have aftershocks, any statistical method that assumes independence of events is problematic. A better approach is to measure the background rate over time, separate from the secondary triggering of aftershocks. In this context, the term “background rate” means the primary earthquakes directly related to the driving stress from both tectonic and anthropogenic sources, and therefore the background rate can vary in time.

To separate background and aftershock seismicity, we used the Epidemic Type Aftershock Sequence (ETAS) model, in which background and aftershock rates are parameterized by using standard empirical relationships, and the resulting parameter set is simultaneously fit from the observed catalog (18). This strategy builds on previous work on identifying fluid-modulated signals in natural and induced seismicity (1821). The seismicity is modeled as a Poissonian process with history-dependent rate RETAS at time tE that is a combination of the modified Omori’s law, which describes the temporal decay of aftershocks; the Gutenberg-Richter relationship, which describes the magnitude distribution; and the aftershock productivity relationshipEmbedded Image (1)where μ is background seismicity rate, and the term inside the summation describes the component of seismicity due to aftershock sequences. In this formulation, KE is the aftershock productivity of a mainshock, α describes the efficiency of earthquakes of a given magnitude at generating aftershocks, c and p are Omori’s law decay parameters, Mi and ti are the magnitude and time of the ith event in the catalog, and Mc is the magnitude completeness threshold of the catalog (supplementary materials, materials and methods) (18, 22).

We performed maximum likelihood fits on overlapping 2-year windows and tracked background seismicity rate over time (Fig. 2 and fig. S3). We assumed that the background rate is stationary for relatively short times—that is, over the 2-year interval—but varies over the long duration of the full catalog. The 2-year interval allows sufficient events to have well-resolved parameters (fig. S4), while still capturing the high-frequency fluctuations (fig. S5).

The increasing trend of the background rate μ in Fig. 2 imitates all the metrics of fluid volume (injection, production, and net production, defined here as production minus injection), particularly in the earlier years of injection (until ~1991). In later years (2006 to 2012), μ tracks net production more closely (Table 1). In between, seismicity may track net production with a baseline shift (Fig. 2C), but the correlations are much less clear. Both total and net production seem important.

Table 1 Cross-correlation of operational parameters with seismicity rate.

Reported values are maxima of the normalized cross-correlation between background seismicity μ reported in Fig. 3 and the given operational rate, with means removed. Normalization is by the geometric mean of the autocorrelations. Lags are time shifts corresponding to the maximum cross-correlation and are reported in months. Lags are restricted to positive (seismicity lagging behind operation) so as to limit cases to only physical scenarios. Time series run from 1 January to 1 January of the reported years. Alternative model assumption results are in table S1 (supplementary text).

View this table:

The time variable behavior is captured by the best-fit coefficients of a linear modelΜ = c1I + c2N (2)where I and N are the injection and net production rate, respectively, and c1 and c2 are the coefficients. Because the correlation between total injection and total production is extremely high, only one of those variables is used in Eq. 2 so as to ensure a well-constrained solution. We measured time variation by fitting Eq. 2 over a moving data window that is longer in duration than the window used to fit the ETAS model and much shorter than the full study interval. A linear least-squares fit in a 6-year window with 0.5-year increments (~90% overlap) captures the essential features (Fig. 3).

Figure 3 Results of linear model of seismicity based on a combination of injection and net production.

(A) Sample seismicity rate and model prediction of seismicity rate using the observed fluid data and the best-fit linear model of Eq. 2. (B) Number of earthquakes per day triggered per rate of net volume of fluid extracted or total fluid injection. Symbols are best-fit coefficients for Eq. 2. The injection values are coefficient c1 in Eq. 2, and net production values are c2. Error bars are 2 SD of model estimates based on the linear regression.

The combined model of Eq. 2 results in well-constrained model parameters after the initial growth of operations in the mid-1980s. An F test rejects the null hypothesis of insignificant fit at the 95% level for all time periods except during the rapid growth of the field in 1993. The early years have substantial uncertainty in the fit coefficients, as might be expected during the highly nonsteady initial phase of operations. During the well-fit period, the seismicity rate generated per monthly volume of fluid injected steadily decreases over time. Net production generated more earthquakes per fluid volume than did injection both early and late in the study period.

In addition to the difference in correlations over time, there is a difference in phase lags between seismicity and the various operational parameters. The phase lag corresponding to the maximum correlation of seismicity relative to net production is usually 0 months (Table 1). However, the intervals with the strongest correlations between injection or production and seismicity can have several-month phase lags. Over the full data set, the maximum correlation between seismicity and injection has a lag of 8 months (Table 1). In intervals such as 1991–2006, unphysical, large phase lags accompany low correlations and are another indicator of the lack of predictive power of a single operational parameter during these periods.

We conclude from the observed correlations and the F test that net production volume combined with injection information is a good predictor of the seismic response in the short term for a fully developed field. Much previous work (3, 57, 23) has focused on the increase in pore pressure (decrease in effective stress) from injection as the primary driver of seismicity, and the importance of net production (volume removed) suggests that seismicity is responding to a more complex process. Earthquakes responding to net volume loss with no phase lag may imply that seismicity responds to elastic compaction and subsidence and not simply diffusion of high pore pressures at injection sites (8, 9).

An important issue for any induced seismicity study is the possibility of triggering a damaging earthquake. Like most earthquake sequences, the Salton Sea Geothermal Field seismicity is dominated by small earthquakes, and the magnitude distribution follows the Gutenberg-Richter relationship: The number of earthquakes of magnitude greater than or equal to M is proportional to 10bM, where b is nearly 1 (the maximum likelihood estimate of b for 1982–2012 is 0.99) (fig. S6). Static and dynamic stresses have been observed to trigger earthquakes on disconnected fault networks (24, 25), and the Gutenberg-Richter relationship generally holds for the aggregated sequences (18, 26). The major limitation in applying Gutenberg-Richter in a particular region is estimating the maximum size of an earthquake that can be hosted by the local faults. The largest earthquake in the Salton Sea Geothermal field region during the study period was of 5.1 magnitude, and the neighboring, highly strained San Andreas fault can have earthquakes of magnitude at least 8.

Supplementary Materials


Figs. S1 to S7

Table S1

References (2836)

References and Notes

  1. Acknowledgments: This work was funded in part by the Southern California Earthquake Center (SCEC contribution 1752). SCEC is funded by NSF Cooperative Agreement EAR-0106924 and U.S. Geological Survey Cooperative Agreement 02HQAG0008. Geothermal operation data are archived and distributed by the California Department of Conservation ( and seismic data are from the Southern California Seismic Network ( Comments from T. Lay, A. Shuler, P. Fulton, and M. Clapham improved the manuscript, and P. Fulton’s assistance with figures is greatly appreciated.
View Abstract

Navigate This Article