Skip to main content
Log in

Frequency–Magnitude Distribution of −3.7 ≤ M W  ≤ 1 Mining-Induced Earthquakes Around a Mining Front and b Value Invariance with Post-Blast Time

  • Published:
Pure and Applied Geophysics Aims and scope Submit manuscript

Abstract

We investigated frequency-magnitude distribution (FMD) of acoustic emissions (AE) occurring near an active mining front in a South African gold mine, using a catalog developed from an AE network, which is capable of detecting AEs down to M W  −5. When records of blasts were removed, FMDs of AEs obeyed a Gutenberg−Richter law with similar b values, not depending on post-blasting time from the initial 1-min interval through more than 30 h. This result denies a suggestion in a previous study (Richardson and Jordan Bull Seismol Soc Am, 92:1766–1782, 2002) that new fractures generated by blasting disturb the size distribution of background events, which they interpreted as slip events on existing weak planes. Our AE catalog showed that the GR law with b ∼ 1.2 was valid between M W  −3.7 and 0 for AEs around the mining front. Further, using the mine’s seismic catalog, which covers a longer time period of the same area, we could extend the validity range of the GR law with the same b value up to M W 1.

This is a preview of subscription content, log in via an institution to check access.

Access this article

Price excludes VAT (USA)
Tax calculation will be finalised during checkout.

Instant access to the full article PDF.

Institutional subscriptions

Fig. 1
Fig. 2
Fig. 3
Fig. 4
Fig. 5
Fig. 6
Fig. 7

Similar content being viewed by others

References

  • Abercrombie, R. E. (1996), The magnitude-frequency distribution of earthquakes recorded with deep seismometers at Cajon Pass, southern California, Tectonophysics, 261, 1–7, doi:10.1016/0040-1951(96)00052-2.

  • Abercrombie, R. E., and J. N. Brune (1994), Evidence for a constant b-value above magnitude 0 in the southern San Andreas, San Jacinto and San Miguel Fault Zones, and at the Long Valley Caldera, California, Geophys. Res. Lett., 21, 1647–1650, doi: 10.1029/94GL01138.

  • Aki, K. (1965). Maximum likelihood estimate of b in the formula log N = a − bM and its confidence limits, Bull. Earthq. Res. Inst. Univ. Tokyo, 43, 237–239.

  • Aki, K. (1987), Magnitude-frequency relation for small earthquakes: A clue to the origin of \(f_{max}\) of large earthquakes, J. Geophys. Res., 92, 1349–1355, doi:10.1029/JB092iB02p01349.

  • Aki, K., and P. G. Richards, Quantitative Seismology, (W. H. Freeman and company, New York 1980).

  • Andrews, D. J. (1976), Rupture propagation with finite stress in antiplane strain, J. Geophys. Res., 81, 3575–3582, doi:10.1029/JB081i020p03575.

  • Boatwright, J. (1978), Detailed spectral analysis of two small New York state earthquakes, Bull. Seismol. Soc. Am., 68, 1117–1131.

  • Boettcher, M. S., A. Mcgarr, and M. Johnston (2009), Extension of Gutenberg-Richter distribution to \(M_W \, -\) 1.3, no lower limit in sight, Geophys. Res. Lett., 36, L10307, doi:10.1029/2009GL038080.

  • Dieterich, J. H., A model for the nucleation of earthquake slip, In Earthquake Source Mechanics, Geophys. Monogr. Ser., vol. 37, (ed. Das, S., J. Boatwright, and C. H. Scholz) (American Geophysical Union, 1986) pp. 37–47.

  • Finnie, G. J. (1999), Using neural networks to discriminate between genuine and spurious seismic events in mines, Pure Appl. Geophys., 154, 41–56, doi:10.1007/s000240050220.

  • Gutenberg, B., and C. F. Richter (1944), Frequency of earthquakes in California, Bull. Seismol. Soc. Am., 34, 185–188.

  • Hanks, T. C., and H. Kanamori (1979), A moment magnitude scale, J. Geophys. Res., 84, 2348–2350, doi:10.1029/JB084iB05p02348.

  • Horiuchi, S., Y. Horiuchi, Y. Iio, Y. Sawada, S. Sekine, H. Nakamura, T. Okada, M. Nakatani, and M. Naoi (2011), Automatic P and S wave arrival time picking compared to manually picking, IUGG General Assembly, 3430, Melbourne, Australia.

  • Ida, Y. (1972), Cohesive force across the tip of a longitudinal-shear crack and Griffith’s specific surface energy, J. Geophys. Res., 77, 3796–3805, doi:10.1029/JB077i020p03796.

    Google Scholar 

  • Ide, S., G. C. Beroza, S. G. Prejean, and W. L. Ellsworth (2003), Apparent break in earthquake scaling due to path and site effects on deep borehole recordings, J. Geophys. Res., 108, 2271, doi:10.1029/2001JB001617.

  • Ishimoto, M., and K. Iida (1939), Observations sur les séismes enregistrés par le microsismographe construit dernièrement (in Japanese with French abstract), Bull. Earthquake Res. Inst. Univ. Tokyo, 17, 443–478.

  • Kawakata, H., N. Yoshimitsu, M. Nakatani, J. Philipp, I. Doi, M. Naoi, T. Ward, V. Visser, G. Morema, S. Khambule, T. Masakale, A. Milev, R. J. Durrheim, L. Ribeiro, M. Ward, and H. Ogasawara (2011), Monitoring transmitted waves across a fault with a high potential for mining induced earthquakes—The Ezulwini gold mine in South Africa, AGU Fall Meeting, S31C-2257, San Francisco.

  • Kwiatek, G., K. Plenkers, G. Dresen, and JAGUARS Research Group (2011), Source parameters of picoseismicity recorded at Mponeng Deep Gold Mine, South Africa: Implications for scaling relations, Bull. Seismol. Soc. Am., 101, 2592–2608, doi:10.1785/0120110094.

  • Kwiatek, G., K. Plenkers, M. Nakatani, Y. Yabe, G. Dresen, and Jaguars-Group (2010), Frequency-magnitude characteristics down to magnitude −4.4 for induced seismicity recorded at Mponeng gold mine, South Africa, Bull. Seismol. Soc. Am., 100, 1165–1173, doi:10.1785/0120090277.

  • Malin, P. E., S. N. Blakeslee, M. G. Alvarez, and A. J. Martin (1989), Microearthquake imaging of the Parkfield Asperity, Science, 244, 557–559, doi:10.1126/science.244.4904.557.

  • Manthei, G. (2005), Characterization of acoustic emission sources in a rock salt specimen under triaxial compression, Bull. Seismol. Soc. Am., 95, 1674–1700, doi:10.1785/0120040076.

  • Mendecki, A. J., Seismic monitoring systems, In Seismic Monitoring in Mines (ed. Mendecki, A. J.) (Chapman and Hall, 1997).

  • Nakatani, M., Y. Yabe, J. Philipp, G. Morema, S. Stanchits, G. Dresen and Jaguars Group (2008), Acoustic emission measurements in a deep gold mine in South Africa—Project overview and some typical waveforms, Seismol. Res. Lett. 79, 311.

  • Naoi, M., M. nakatani, Y. Yabe, G. Kwiatek, T. Igarashi, and K. Plenkers (2011), Twenty thousand aftershocks of a very small (M 2) earthquake and their relation to the mainshock rupture and geological structures, Bull. Seismol. Soc. Am., 101, 2399–2407, doi:10.1785/0120100346.

  • Richardson, E., and T. H. Jordan (2002), Seismicity in deep gold mines of South Africa: Implications for tectonic earthquakes, Bull. Seismol. Soc. Am., 92, 1766–1782, doi:10.1785/0120000226.

  • Rydelek, P. A., and I. S. Sacks (1989), Testing the completeness of earthquake catalogues and the hypothesis of self-similarity, Nature, 337, 251–253, doi:10.1038/337251a0.

  • Scholz, C. H. (1988), The critical slip distance for seismic faulting, Nature, 336, 761–763, doi:10.1038/336761a0.

  • Shi, Y., and B. A. Bolt (1982). The standard error of the magnitude-frequency b value, Bull. Seismol. Soc. Am., 72, 1677–1687.

  • Spottiswoode, S. M. (2010), Mine seismicity: prediction or forecasting?, J.S.A. Inst. Min. Metall., 110, 11–20.

  • Spottiswoode, S. M., and L. M. Linzer (2004), Improved seismic locations and location techniques, SIM 020304 Final Project Report (Appendix 1), Mine Health and Safety Council, Johannesburg, South Africa.

  • Taylor, D. W. A., J. A. Snoke, I. S. Sacks, and T. Takanami (1990), Nonlinear frequency-magnitude relationships for the Hokkaido corner, Japan, Bull. Seismol. Soc. Am., 80, 340–353.

  • Umino, N., and I. S. Sacks (1993), Magnitude-frequency relations for northeastern Japan, Bull. Seismol. Soc. Am., 83, 1492–1506.

  • Utsu, T. (1965), A method for determining the value of b in a formula logn = a − bM showing the magnitude-frequency relation for earthquakes (in Japanese with English abstract), Geophys. Bull. Hokkaido Univ., 13, 99–103.

  • Wiemer, S., and M. Wyss (2000). Minimum magnitude of completeness in earthquake catalogs: Examples from Alaska, the western United States, and Japan, Bull. Seismol. Soc. Am., 90, 859–869, doi:10.1785/0119990114.

  • Woessner, J., and S. Wiemer (2005). Assessing the quality of earthquake catalogues: Estimating the magnitude of completeness and its uncertainty, Bull. Seismol. Soc. Am., 95, 684–698, doi:10.1785/0120040007.

    Google Scholar 

Download references

Acknowledgments

We thank Gold One International Ltd., operator of the Cooke 4 shaft, for approving our publication, and also personnel at the mine and at the Institute of Mine Seismology for facilitating our project. This project was funded by Grants-in-Aid from the Ministry of Education, Culture, Sports, Science and Technology (MEXT) of Japan (21224012, 21246134). It was also funded by JST/JICA, SATREPS, and the MEXT’s Observation and Research Program for Prediction of Earthquakes and Volcanic Eruptions.

Author information

Authors and Affiliations

Authors

Corresponding author

Correspondence to Makoto Naoi.

Appendices

Appendix A: Calibration of AE Sensor Frequency Responses

Here, we describe how we evaluated the frequency response of the AE sensors, largely following the method of Kwiatek et al (2011).

Let S i (f) be the amplitude spectrum of the ith-component ground acceleration incident on an accelerometer during a seismic event, A acci(f) be the ith-component accelerometer output spectrum, and I acc(f) be the accelerometer’s frequency response:

$$ A_{{\rm{acci}}}(f)=I_{{\rm{acc}}}(f)S_i(f). $$
(A1)

If we define the total spectrum A acc(f) of the three-component accelerometer output by

$$ A_{{\rm{acc}}}(f)=\sqrt{\textstyle\sum_{i=1}^3 A_{{\rm{acci}}}^2(f)}, $$
(A2)

we have

$$ A_{{\rm{acc}}}(f)=I_{{\rm{acc}}}(f)S_{{\rm{acc}}}(f), $$
(A3)

where \(S_{{\rm{acc}}}(f)=\sqrt{\sum_i S_i(f)^2}\) is the total acceleration spectrum.

Meanwhile, the AE sensor output spectrum A AE is given by

$$ A_{{\rm{AE}}}(f; \theta_1, \theta_2)=I_{{\rm{AE}}}(f; \theta_1, \theta_2)S_{{\rm{acc}}}(f), $$
(A4)

where I AE(fθ 1θ 2) is the frequency response of the AE sensor, θ 1 is the incidence angle of waves on the AE sensor with 0° pointing right to the sensor front and 180° pointing right to the sensor’s back, and θ 2 is the oscillation direction of the incident waves. The AE sensors are directional, and their frequency responses are thought to depend on θ 1 and θ 2. We assume their frequency responses are axially symmetric because of the sensors’ symmetry.

Let us also denote, for a trio of an AE sensor, an accelerometer and a seismic source, the source-AE sensor distance by R AE, the source-accelerometer distance by R acc, and the accelerometer-AE sensor distance by R. When RR acc and RR AE, the ground acceleration spectrum of waves incident on the accelerometer, S acc(f), is thought to agree with its AE sensor counterpart. The ratio of the frequency responses of the AE sensor and the accelerometer, T(f; θ 1, θ 2), can therefore be written as follows:

$$ T(f; \theta_1, \theta_2)=\frac{I_{{\rm{AE}}}(f; \theta_1, \theta_2)}{I_{{\rm{acc}}}(f)}=\frac{A_{{\rm{AE}}}(f; \theta_1, \theta_2)}{A_{{\rm{acc}}}(f)}, $$
(A5)

so that the total AE sensor output spectrum, A AE, can be converted to its accelerometer equivalent, A acc, by

$$ A_{{\rm{acc}}}(f)=A_{{\rm{AE}}}(f; \theta_1, \theta_2)/T(f; \theta_1, \theta_2). $$
(A6)

We used Eq. (A5) to calculate T(f; θ 1, θ 2). For P waves, we do not need to consider the θ 2 dependency because the angle is determined uniquely from θ 1. For S waves, although the θ 2 dependency should be considered independently of θ 1, we cannot correct for the θ 2 dependency for most AE sensors that do not have nearby accelerometers without knowing the focal mechanism of individual events, even if we know the θ 2 dependency of the sensor. We therefore ignore θ 2 dependency of sensor characteristics and assume I AE(f; θ 1,θ 2) = I AE(f; θ 1) and T(f; θ 1, θ 2) = T(f; θ 1).

We used records of two AE sensors (channels 25 and 42), only 2.0 m and 3.7 m apart from the nearest accelerometer, respectively, to determine T(fθ 1), whereby we took the average of the spectral ratios for both sensor pairs and applied the results to all AE sensors. We used only seismic events with hypocenters located within 200 m of the center of gravity of the AE network but at least 30 m apart from both the accelerometer and the AE sensor, and with incidence angles differing 10° or less for both types of sensors. For channel 25, we used P waves from 38,196 events and S waves from 33,413 events. For channel 42, the number of events used was 30,981 for P waves and 29,299 for S waves.

We calculated the spectra by using 2.048 ms (1,024 samples) portions of the accelerograms and the AE waveforms that include the direct P and S phases. After taking moving averages of the spectra with a window width of 2 kHz, we used equation (A5) to calculate T(fθ 1), and stacked the results separately in different bins of θ 1, the incidence angle with respect to the AE sensor. The θ 1 bins overlap partially, with a bin width of 20° and an angle increment of 5°.

Figure 8 shows the stacks of T(f; θ 1) in the different θ 1 bins. We stacked the spectral ratios only in the frequency bands where the S/N ratios, with respect to the pre-P noise spectra, were 10 dB or more. Further, we excluded the frequency bands where <10 T(f; θ 1) samples were available for stacking, labeling such bands as “unfit for analysis” and suppressing them as in Fig. 8.

Fig. 8
figure 8

Frequency responses of AE sensors in terms of AE sensor-to-accelerometer spectral ratios T(f; θ 1), estimated for AE sensors installed adjacent to accelerometers. The T(f; θ 1) was stacked separately in different bins of θ 1, the incidence angle with respect to the AE sensor. The θ 1 bins overlap partially, with a bin width of 20° and an angle increment of 5°. Indicated in color is the central angle of each bin. Frequency bands where <10 T(f; θ 1) samples were available for stacking, because of poor S/N ratios, are not shown. Stacked spectral ratios for an accelerometer-AE sensor pair near the AE network center (inter-sensor distance 2.0 m): a results for P waves. b Results for S waves. Stacked spectral ratios for an accelerometer-AE sensor pair on the north side of the AE network (inter-sensor distance 3.7 m): c results for P waves. d Results for S waves. Black curves average of all stacked spectral ratios shown in ad. We used this average spectral ratio to convert AE sensor spectra into equivalent acceleration spectra

In Fig. 8, we find the spectral ratios to behave slightly differently for the P and S waves and to vary systematically with the wave incidence angle θ 1, but both those variations are smaller than, or at most comparable to, the differences between the two AE sensors. We therefore decided to ignore all dependency on θ 1 and all differences between the P and S waves, and to use the average of all stacked results shown in Fig. 8a–d (solid black curves) to convert the AE sensor spectra into equivalent acceleration spectra. The frequency bands labeled “unfit for analysis” were not used in the averaging.

Appendix B: Method to Estimate Moment Magnitude of AEs

We estimated M o by fitting theoretical spectra to ground displacement spectra derived from accelerograms and AE sensor records. The AE sensor records were used after correction by T(f), an AE sensor-to-accelerometer frequency response ratio estimated in situ, estimated in "Appendix A".

Let S i (f) denote the ground acceleration spectrum of waves incident on an AE sensor or an accelerometer, and let S acc(f) be the total ground acceleration spectrum integrated from the three components: \(S_{{\rm{acc}}}(f)=\sqrt{\sum_i S_i(f)^2}. \) The total ground displacement spectrum D(f) can be obtained from the total acceleration spectrum S acc(f) as follows:

$$ D(f)=\frac{1}{(2\pi f)^2}S_{{\rm{acc}}}(f). $$
(B1)

We evaluated M o by fitting the theoretical total ground displacement spectrum D th(f) to the D(f) thus obtained. Subsequent analysis is similar in method to Ide et al. (2003). From Boatwright (1978),

$$ D^{{\rm{th}}}(f) = \frac{R^c}{4\pi \rho v^3 r} \frac{M_o}{\sqrt{1+(f/f_c)^4}} \exp\left( -\frac{\pi f r}{Qv} \right), $$
(B2)

where R c is the radiation pattern coefficient, ρ is the medium density, v is the elastic wave velocity, r is the source-receiver distance, f is the frequency, f c is the corner frequency, and Q is the quality factor. Taking the logarithms on both sides, we have

$$ \log D^{{\rm{th}}}(f) =\log M_o + \log \frac{R^c}{4\pi \rho v^3 r} -\frac{1}{2} \log \left[ 1+(f/f_c)^4 \right] -\frac{\pi f r}{Qv}. $$
(B3)

Because log D th(f) is linearly related to log M o ,   M o may be obtained by linear least squares fitting with the other parameters fixed. We defined the residual, Res, by

$$ {{\rm{Res}}} = \left[\log D(f)-\log D^{{\rm{th}}}(f)\right]^2 $$
(B4)

and used a least squares method to determine M o and a grid search method to determine f c that produced a D th(f) that minimizes the residual. We followed the definition of Hanks and Kanamori (1979) to convert M o to M W :

$$ M_W = \frac{2}{3} \left(\log_{10}M_o - 9.1\right). $$
(B5)

We analyzed the spectra of the individual sensor records near their P and S wave travel time readings (or theoretical travel times in their absence). We applied a cosine taper to the initial and final 5 % parts of the time-domain waveforms before putting them to spectral analysis.

The D(f) was obtained by resampling the FFTs of time-domain waveforms at logarithmically equal intervals and taking their moving averages with a window width of 0.05 on the logarithmic frequency scale. Fitting was conducted only on frequency bands where the S/N ratios, with respect to the noise spectra obtained from pre-P signals, were 10 dB or more. We labeled the M W estimates as unavailable when no continuous frequency band with S/N ≥ 10 dB was broader than a threshold width (0.4 on the logarithmic frequency scale for AE sensors and 0.5 for accelerometers).

We also labeled M o estimates as unavailable when logf c  − logf min < 0.2, where f c is the corner frequency estimate and f min is the minimum frequency used in the spectral fitting. That is because M o tends to be underestimated considerably when spectral fitting is done only in high frequency bands above the corner frequency.

We used a data length of 4,096 samples (∼8.2 ms) for the spectral analysis, but repeated the analysis with a doubled data length of 8,192 samples (∼16.4 ms) when the shorter data length produced an M W estimate of −3 or greater. If the doubled data length produced an M W estimate of −2 or greater, we repeated the analysis further with a tripled data length of 12,288 samples (∼24.6 ms). We also repeated the analysis with similarly extended data lengths when the M o estimate was labeled unavailable. We labeled it as “indefinite” for the sensor in question when it remained unavailable for a data length of 12,288 samples. When the SP time was short, we reduced the data lengths appropriately to prevent the P and S phases from being mixed up in the analysis.

We used the parameters V P = 5,700 m/s, V S = 3,600 m/s, ρ = 2,700 kg/m3,   Q P = 150, Q S = 200, \(R_P^c=2/\sqrt{15}\) and \(R_S^c=\sqrt2/\sqrt5\) in Eqs. (B2) and (B3). The V P and V S are the same values used to determine the hypocenters. The Q P and Q S were estimated by trial and error to allow the observed acceleration spectra to be well explained by Eq. (B3). Sound transmission tests have also produced the same Q P value (Kawakata et al. 2011). The R c P and R c S are rms values averaged over all azimuths (Aki and Richards 1980).

Figure 9 illustrates how theoretical spectra derived from Eq. (B3) were fit to the observed records. The figure shows that the theoretical spectra well explain the observed spectra for both the accelerograms (Fig. 9a, c) and the AE sensor records (Fig. 9b, d).

Fig. 9
figure 9

Examples showing how fitting was conducted on P and S wave spectra to evaluate the seismic moment M 0. a Waveform records of a three-component 25-kHz accelerometer for an event that occurred at 05:20:14.187 on 6 October 2011. The hypocenter was located 87.2 m from the accelerometer. The red, blue, and black bars beneath the waveforms show the time windows used to calculate the PS, and noise wave spectra, respectively. b Waveform records of an AE sensor for an event that occurred at 06:02:16.812 on 9 August 2011. The hypocenter distance was 39.6 m. c Total ground displacement spectra D(f) calculated from the waveforms shown in a. Red, blue, and black solid curves D(f) of the PS, and noise waves, respectively. The highlighted parts show the frequency bands with S/N ratios ≥ 10 dB. Dashed curves: theoretical spectra D th(f), which include the effect of Q, and provide the best fit for the highlighted parts of the D(f) curves. Dotted curves theoretical spectra without the effect of Q. The fitting led to a moment magnitude estimate of M WP   −1.9 based on the P waves and M WS   −2.0 based on the S waves, both before correction for the systematic differences between M WP and M WS . d Total ground displacement spectra D(f) calculated from the waveforms shown in b by using T(f) ("Appendix A"). See legend to c. M WP   −3.7,   M WS   −4.0

In this study, we used records of all six accelerometers and 11 of the 24 AE sensors in the analysis. We used only part of the available AE sensor records to save on computation time.

Of the 365,237 AE events, which were well located by the automatic picking and location software (Sect. 3.2), 17,050 had their hypocenters located within 150 m of AEnet’s center and their M o and M W determined for both the P and S wave parts of all six accelerograms. In those analysis results, the M W estimates from P waves (M WP ) exceeded their S wave counterparts (M WS ) by 0.24 on average. A similar tendency was confirmed for every single monitoring station, be it an accelerometer or an AE sensor. Possible causes include biases in radiation patterns and imprecision of the seismic wave velocities used.

In some cases, only one of the M WP and M WS was available. We therefore subtracted 0.12 from every single M WP and added 0.12 to every single M WS , including estimates from AE sensor records, to equilibrate the M W estimates from both P and S wave records on average. We took the mean of all available M WP and M WS and rounded the outcome to the nearest tenth to determine the final M W estimate. The standard deviation of M W among stations and phases for event i,   σ Mi , is shown in a histogram of Fig. 10. For 249,017 events which have more than ten available M W from different stations and phases, the mean of σ Mi was 0.136. In calculating σ Mi , the corrections of 0.12 on the M P W and M S W have been taken into account.

Fig. 10
figure 10

The histogram of σ Mi , the standard deviation of M W among stations and phases for event i, for 249,017 events which have more than ten available M W from different stations and phases

It should be noted that no low-frequency content below 1 kHz was available for estimating M o using AE sensor records, because they were all high-pass filtered at 1 kHz before analog-to-digital conversion. Accelerograms were available for use down to lower-frequency bands than the AE sensor records, because the accelerograms were high-pass filtered only at 50 Hz. The M o risks being underestimated when only limited low-frequency content is available, although we took care, as we stated above, to label the M o estimate as unavailable when only high-frequency content above the corner frequency was used in the fitting.

Figure 11 shows the relationship between the mean M W values determined using AE sensor records, M AE W , and those determined from accelerograms, M acc W , for 287,025 events, which were located within 150 m of the center of gravity of AEnet and for which records from at least one accelerometer and at least one AE sensor were used to evaluate M o . The corrections of 0.12 on the M P W and M S W have been taken into account. It shows that M AE W tends to fall short of M acc W and is likely underestimated for events above M W  −2. To avoid that bias, we used M acc W as the final M W value, without using analysis results from the AE sensor records, when M acc W was larger than −2.5.

Fig. 11
figure 11

Relationship between M AE W (the mean of M W values determined using AE sensor records) and M acc W (the mean of M W values determined using accelerograms) for the 287,025 events that were located within 150 m of the center of gravity of AEnet and for which records from at least one accelerometer and at least one AE sensor were used to evaluate M o

Appendix C: Consistency Between M W Estimates from the AE Network and the Mine’s Network

To compare size distributions from two different catalogs, the AEnet catalog and the Mnet catalog, we need to check the consistency between the M W estimates from the AEnet records (M AEnet W ) and the M W given in the Mnet catalog (M Mnet W ), using the events that are listed in both catalogs.

As we stated in Sect. 2, the clocks in both monitoring networks are not synchronized. Therefore, we first determine the time difference curve between both system clocks as described below, and used that result to correct them before we singled out what we thought to be pairs of identical event records in both catalogs, for which consistency of M W estimates were checked.

To estimate the time difference curve, we first searched for candidates of such pairs. For every single event in the Mnet catalog, we picked up events in the AEnet catalog that have a proximal origin time, location and M W estimate. We used 300 s ≤ dt ≤ 600 s, D ≤ 100 m and |dM W | ≤ 1 for our selection criteria, where dt is the difference between the Mnet and AEnet origin times, D is the distance between the Mnet and AEnet hypocenters, and dM W  = M AEnet W  − M Mnet W . We set the selection standard for dM W somewhat loose because our final objective is to verify the consistency of the M W estimates. The range of D was determined so as not to exclude pairs of identical event records for which mutually distant hypocenters were given by both catalogs, considering that location errors can sometimes be as large as 100 m in the Mnet catalog. As for dt, we started the actual analysis using broader ranges, such as −1,000 s ≤ dt ≤   1,000 s, and narrowed it by reference to the result. We show here only the final result using 300 s ≤ dt ≤ 600 s, which turned out to be the true range of the clock difference in the analysis term, but essentially the same results were obtained when we used broader dt range.

We used only 363,672 events in the AEnet catalog whose hypocenters were located within 200 m of the AEnet center. Figure 12 shows a time-series plot of dt for events in the AEnet catalog that satisfied the above condition, each being the only one event that could be paired with an event in the Mnet catalog. Figure 12 likely contains numerous plots of dt between what in fact are different events, but it does show a single, monotonously rising trend of densely populated data points. That trend likely represents the time difference between both system clocks.

Fig. 12
figure 12

A plot of time differences between the AEnet and the Mnet records, with the horizontal axis standing for the event occurrence times according to the clock of AEnet. Each data point represents an event pair satisfying 300 s ≤ dt ≤ 600 s, D ≤ 100 m and |dM W | ≤ 1, where dt is the difference in origin times between the AEnet catalog and the Mnet catalog, D is the distance between the two catalogued hypocenters, and dM W is the difference between the two catalogued M W values. Cases where two or more events in the AEnet catalog could be paired with an event in the Mnet catalog were not included in the plot. Open circles data points defining a single, monotonously rising trend. Gray triangles data points obviously at odds with the trend. Solid curve time correction curve determined from the open circles

We manually eliminated data that are obviously at odds with the time differences of preceding and following events (gray triangles in Fig. 12) and determined a time difference curve (solid curve in Fig. 12) by taking the moving averages of the remaining data points (open circles) with a window width of 2 days.

The origin times in the Mnet catalog were used only after correcting for the time differences according to the curve. We subsequently used stricter selection criteria to identify which pairs of events are believed to represent the same events with strong certainty and compared the M W estimates in those pairs.

Let dt′ be the residual difference between the catalogued origin times following the corrections. We singled out events in the Mnet catalog, with hypocenters located within 100 m of the AEnet center, each being associable with only one event in the AEnet catalog that satisfied −5 s ≤ dt′ ≤ 5 s, D ≤ 50 m and |dM W | ≤ 1. Figure 13 shows the relationship between the M Mnet W and the M AEnet W for the 117 events thus singled out. It demonstrates that the M AEnet W tends to be underestimated with respect to the M Mnet W for events with M Mnet W ≥ −1.3.

Fig. 13
figure 13

Differences between the M AEnet W (M W determined by AEnet) and the M Mnet W (M W determined by Mnet) for 117 pairs that are most certainly associated with the identical events. The sizes of the circles and the numeral figures superimposed on them indicate the number of events falling on each (M Mnet W ,   dM W ) grid point. The gray line approximating the trend was used to correct the underestimation of M AEnet W

As we stated in "Appendix B", we used the accelerograms alone in estimating the M AEnet W for events with M acc W ≥  −2.5 to avoid the influence of the 1 kHz high-pass filter on the AE sensor waveforms, but the accelerograms, for their own part, have been high-pass filtered at 50 Hz and their lower-frequency content has been eliminated. That presumably accounts for why the M W estimates from AEnet are underestimated for events upward of M Mnet W   −1.3. To correct those underestimations, we used the relationship between the M W estimates in the Mnet catalog and the dM W , indicated by the kinked gray line in Fig. 13.

Appendix D: b Value Analysis

Here, we describe how we evaluated b values, their standard errors, and the minimum validity magnitude M c of the Gutenberg–Richter law.

We calculated the b values with the maximum likelihood method (Aki 1965; UTSU 1965) by using the formula:

$$ b = \frac{1}{\overline{M_W} - (M_c-\Updelta M/2)} \log_{10} e, $$
(D1)

where \(\overline{M_W}\)is the mean magnitude of all M c and greater events, and \(\Updelta M\) is the sampling increment for M W in the catalog. Following Shi and Bolt (1982), we used the following formula to evaluate e b , or the standard error of a b value estimate:

$$ e_b = 2.30 b^2 \sqrt{\sum_{i=1}^{n} \left( M_i-\overline{M_W} \right)^2/n(n-1)}, $$
(D2)

where M i represents the magnitude of the ith event and n is the total number of events.

We used the maximum curvature method of Wiemer and Wyss (2000) to estimate M c , the minimum validity magnitude of the Gutenberg-Richter law in a given data set. For a given threshold magnitude M min, that method uses a “goodness of fit”, R, defined by

$$ R(a, b, M_i) = 100 - \frac{ \displaystyle \sum\nolimits_{M_i = M_{{\rm{min}}}}^{M_{{\rm{max}}}} | B_i - S_i | } {\displaystyle \sum\nolimits_{i}B_i} 100, $$
(D3)

to evaluate the discrepancy between the size distribution predicted by the Gutenberg–Richter law and the observed size distribution, where a and b are parameters of the Gutenberg–Richter law, M i is the M W of the ith magnitude bin, and B i and S i are the observed and predicted cumulative numbers of events, respectively, in the ith magnitude bin. We evaluated R(abM i ) for different choices of M min , and designated the smallest M min that gave R in excess of a certain threshold value (typically 90) as our estimate of M c , the magnitude where the frequency-magnitude distribution begins to deviate from the power law.

Woessner and Wiemer (2005) said this method tends to slightly underestimate M c , and recommended adding 0.2 to the M c estimates. Following their recommendation, we added 0.2 to the minimum M min with R > 90 to obtain our final M c estimate.

Rights and permissions

Reprints and permissions

About this article

Cite this article

Naoi, M., Nakatani, M., Horiuchi, S. et al. Frequency–Magnitude Distribution of −3.7 ≤ M W  ≤ 1 Mining-Induced Earthquakes Around a Mining Front and b Value Invariance with Post-Blast Time. Pure Appl. Geophys. 171, 2665–2684 (2014). https://doi.org/10.1007/s00024-013-0721-7

Download citation

  • Received:

  • Revised:

  • Accepted:

  • Published:

  • Issue Date:

  • DOI: https://doi.org/10.1007/s00024-013-0721-7

Keywords

Navigation