Strong Aftershock Study Based on Coulomb Stress Triggering—A Case Study on the 2016 Ecuador M_w 7.8 Earthquake

Abstract

The 2016 Ecuador M 7.8 earthquake ruptured the subduction zone boundary between the Nazca plate and the South America plate. This M 7.8 earthquake may have promoted failure in the surrounding crust, where six M ≥ 6 aftershocks occurred following this mainshock. These crustal ruptures were triggered by the high coulomb stress changes produced by the M 7.8 mainshock. Here, we investigate whether the six M ≥ 6 aftershocks are consistent with the positive coulomb stress region due to the mainshock. To explore the correlation between the mainshock and the aftershocks, we adopt a recently published high-quality finite fault model and focal mechanisms to study the coulomb stress triggers during the M 7.8 earthquake sequence. We compute the coulomb failure stress changes (ΔCFS) on both of the focal mechanism nodal planes. We compare the ΔCFS imparted by the M 7.8 mainshock on the subsequent aftershocks with the epicenter location of each aftershock. In addition, the shear stress, normal stress, and coulomb stress changes in the focal sources of each aftershock are also computed. Coulomb stress changes in the focal source for the six M ≥ 6 aftershocks are in the range of −2.17–7.564 bar. Only one computational result for the M 6.9 aftershock is negative; other results are positive. We found that the vast majority of the six M ≥ 6 aftershocks occurred in positive coulomb stress areas triggered by the M 7.8 mainshock. Our results suggest that the coulomb stress changes contributed to the development of the Ecuador M 7.8 earthquake sequence.

1. Introduction

The 2016 M 7.8 earthquake occurred offshore of the west coast of Northern Ecuador, where the plate boundary between the Nazca and the South America plates lies (Figure 1). During this earthquake, at least 668 people were killed, and 27,732 people were injured. More than 7000 buildings were destroyed, including most of the town of Pedernales and its surrounding urban areas [1]. The urban earthquake disaster is growing throughout the world and has become a serious problem in developing countries. Unless earthquake monitoring is improved, subsequent aftershocks will have more disastrous social and economic consequences [2]. Due to its unique geodynamic environment, seven M > 7 earthquakes have occurred within 250 km of this event since 1900.

The M 7.8 mainshock was followed by six M ≥ 6 aftershocks, the largest being an M 6.9 aftershock on 18 May about one month after the mainshock. The epicenter locations and focal mechanisms of these aftershocks are consistent with the plate boundary subduction zone between the Nazca and the South America plates.

In recent years, many seismology scientists worldwide have focused on coulomb stress triggering and the correlation between the mainshock and subsequent aftershocks in an earthquake sequence. Much research on large earthquake sequences has concluded that stress changes from the mainshock affect the locations of subsequent aftershocks [4,5,6]. The coulomb stress triggering theory has been proposed for evaluating aftershock hazards after great earthquakes. This theory implies that coulomb stress changes due to a mainshock can promote nearby faults closer to failure stress and trigger the occurrence of some aftershocks. It was thought that small coulomb stress changes can alter the likelihood of earthquakes on nearby faults [7,8]. The occurrence of earthquake activity can be promoted when the coulomb stress increases as little as 0.1 bar on a seismogenic fault [9,10]. Small positive changes in coulomb stress due to an earlier mainshock can trigger subsequent earthquakes.

The coulomb stress change triggered by a mainshock has been widely studied to explore its interaction with aftershock activities. Coulomb stress triggering of moderate-to-large earthquakes has been observed in the US [4,11,12,13], in Europe [14,15,16], on the Sunda Trench [17,18,19], and in Japan [20,21,22,23,24,25]. However, in South America and Ecuador, the stress triggering theory has yet to be extensively studied [26,27,28]. The 2016 Ecuador earthquake sequence provides a unique opportunity to test the interaction between a large subduction earthquake that imparted coulomb stress changes and the subsequent seismicity. In this study, based on the Coulomb 3.3 software [29] (3.3, USGS, Menlo Park, CA, USA), which implements the elastic half space of Okada [30], we investigate the coulomb stress triggering theory with the 2016 Ecuador M 7.8 earthquake sequence. Firstly, we collect the focal mechanism solutions of six M ≥ 6 aftershocks from the USGS (US Geological Survey). Secondly, we compute the coulomb stress changes on both nodal planes of these events. We determine the regional stress direction from the principal axes. Thirdly, we compute the coulomb stress, normal stress, and shear stress changes in the focal source of each M ≥ 6 aftershock. Our primary goal is to find the correlation between the coulomb stress changes and the spatial distribution of the following aftershocks. In addition, we want to determine whether the coulomb stress triggering theory can explain the Ecuador M 7.8 earthquake sequence.

2. The Coulomb Stress Triggering Hypothesis

2.1. Coulomb Failure Criterion

The coulomb failure criterion is commonly adopted to quantify the effect of static stress change on a fault plane. The trend of rocks to rupture in a brittle manner is considered to be a function of both shear and confining stresses, commonly formulated as the coulomb failure criterion [8]. A seismogenic fault rupture permanently deforms the crustal medium and modifies the tectonic stress field, increasing or decreasing the stress on nearby faults as a function of the location, geometry, and slip direction. The coulomb failure stress (CFS) is defined as follows [31]:

$$\mathrm{CFS}=\mathsf{\tau }-s-\mathsf{\mu }({\mathsf{\sigma }}_{\mathrm{n}}-p)$$

(1)

where τ is the shear stress along the fault plane (positive in the interred direction of slip). s is cohesion, and μ is apparent friction coefficient, which ranges from 0.6 to 0.8 for most intact rocks [31]. σ_n is the normal stress on the fault plane, and p is the fluid pressure, respectively. The coulomb failure stress change can be defined as follows:

$$\Delta \mathrm{CFS}=\Delta \mathsf{\tau }+{\mathsf{\mu }}^{\prime }\Delta {\mathsf{\sigma }}_{\mathrm{n}}$$

(2)

where Δτ is the shear stress change on a given failure plane (positive in the fault slip direction). Δσ_n is the normal stress change (positive for fault unclamping or compression). μ′ = μ(1 − β) is the effective coefficient of friction, which includes the effects of pore pressure changes and generally ranges from 0 to 0.8. β is Skempton’s coefficient, and the range is from 0 to 1. The previous study showed that the μ′ is typically found to be around 0.4 (μ = 0.75, β = 0.47) for subduction zones [28]. Previous studies have employed values of 0.2 ≤ μ′ ≤ 0.75 to compute coulomb stress changes and have discussed the robustness of the obtained conclusion for different settings of apparent coefficients of frictions [32,33,34]. In this paper, we investigated μ′ = 0.2, 0.4, 0.6, and 0.8, respectively. The computed results show that there is no significant differences, which is consistent with previous conclusions. A priori, the models with μ′ = 0.4 can minimize the uncertainty, as discussed by King et al [8]. Here, we take the computation results with μ′ = 0.4 for the following analysis and discussion. Positive values of the coulomb stress change promote rupture and negative values inhibit rupture. It seems that the coulomb stress increases as low as 0.1 bar can trigger the occurrence of subsequent aftershocks [7,8].

The coulomb stress change relies on the source model of an earthquake, the geometry, the slip direction of the receiver fault (fault receiving stress from a mainshock), and the effective coefficient of friction. Here, we utilize this model to estimate how the Ecuador M 7.8 mainshock transferred stress in the epicentral region. By assuming a Young's modulus of 8.0 × 10⁵ bar and a Poisson’s ratio of 0.25, we compute the coulomb stress change in an elastic half-space. The effective friction coefficient is 0.4. We make two complementary computations. First, we resolve coulomb stress changes on both nodal planes of each M ≥ 6 aftershock. Second, we compute the normal, shear, and coulomb stress changes in the focal source of each M ≥ 6 aftershock. The coulomb failure stress changes, produced by the Ecuador M 7.8 mainshock, is important for us to investigate the cause and mechanism of the larger aftershocks.

2.2. Finite Fault Model

A more realistic finite fault failure model is critical for the following calculation of coulomb stress changes. Here, we choose a variable finite fault model that was inverted from Global Seismic Network (GSN) broadband waveforms by Gavin Hayes (Figure 2) [35]. This model can be downloaded from https://earthquake.usgs.gov/archive/product/finite-fault/us20005j32/us/1460878 551639/web/20005j32_coulomb.inp. Distribution of the amplitude and slip direction for subfault elements of the fault rupture model are determined by the inversion of teleseismic body waveforms and long period surface waves [36,37]. Gavin Hayes used GSN broadband waveforms downloaded from the National Earthquake Information Center (NEIC) waveform server and analyzed 38 teleseismic broadband P waveforms, 9 broadband SH waveforms, and 53 long period surface waves, selected based on data quality and azimuthal distribution. Waveforms were first converted to displacement by removing the instrument response and then used to constrain the slip history using a finite fault inverse algorithm [38].

As shown in Figure 2, the different colors indicate the slip amount. The deeper the color, the greater the amplitude of the slip. Arrows indicate the slip direction (of the hanging wall with respect to the foot wall). For the M 7.8 mainshock, the strike of the fault rupture plane is 29° and the dip is 15°. The rupture surface is approximately 70 km along the strike and 20 km in the dip direction. The seismic moment release based upon this plane is 7.1e + 27 dyne.cm.

3. Coulomb Stress Changes Computed on Nodal Planes

We present here the high-quality focal mechanism solutions of six M ≥ 6 aftershocks from USGS (

Table 1. The basic parameters of the M ≥ 6 aftershocks and stress change calculation results.

Date	Aftershock	Latitude Longitude	Depth (km)	Strike Dip Rake (°)	Stress Change (Bar)
					Shear	Normal	Coulomb
20 April	M 6.2	0.639 80.210	11.5	181,79,82	6.458	2.766	7.564
				38,13,127	6.477	−0.260	6.373
20 April	M 6.0	0.708 80.035	11.0	199,79,97	3.006	−2.104	2.164
				346,13,58	3.015	0.119	3.062
22 April	M 6.0	0.292 80.504	11.5	187,74,86	2.022	3.688	3.497
				22,16,104	2.041	−0.485	1.847
18 May	M 6.7	0.426 79.790	16.0	168,77,79	−2.249	13.871	3.299
				29,16,130	−2.175	5.837	0.159
18 May	M 6.9	0.495 79.616	35.5	183,76,85	−1.711	−1.149	−2.170
				23,15,110	−1.739	6.251	0.762
11 July	M 6.3	0.581 79.638	25.5	178,71,79	−1.465	1.030	1.053
				29,22,119	−1.449	−1.123	1.898

). All of the aftershocks are within the depth range of 0–50 km, bounded by 0.3°–0.9° N latitude and 79°–81° W longitude. Based on the finite fault model of the M 7.8 mainshock inverted by Gavin Hayes [35], we compute the coulomb stress changes on both nodal planes of each subsequent aftershock in their rake directions (Figure 3). We simulate the M 7.8 mainshock as dislocations in an elastic half-space and compute the resulting static stress changes. Our computational results are precise in the near field where the details of the slip distribution are critical. Correspondingly, the details of the slip distribution in the far field become less important than the total moment of the earthquake and the average orientation of the rupture plane [39]. Coulomb stress changes exhibit a lobed pattern of alternating positive and negative stress changes on the receiver fault planes or nodal planes. The computation of coulomb stress changes on both of the nodal planes are a direct test of the coulomb hypothesis.

Figure 3 illustrates the coulomb stress changes triggered by the M 7.8 mainshock on both nodal planes of the six M ≥ 6 aftershocks. As illustrated in Figure 3, the stress field figures for each aftershock are completely different because two different nodal planes are represented. Each figure shows the characteristic areas of coulomb stress increases and areas of coulomb stress decreases caused by the M 7.8 mainshock. The distributions of the six M ≥ 6 aftershocks are consistent with these stress change patterns. Rises in coulomb stress appear to be sufficient to trigger aftershock events, while coulomb stress reductions effectively inhibit them. Among the six M ≥ 6 aftershocks, five of the hypocenters plotted fall in regions of positive stress changes. Conversely, only one hypocenter falls within the stress shadow of the mainshock, as shown in Figure 3i. Overall, we observe that the vast majority of the aftershocks undergo positive coulomb stress change, which would promote or trigger rupture. These findings suggest that the M 7.8 mainshock may have triggered the six M ≥ 6 aftershocks.

It should be noted that the focal mechanism solutions display geometry similar to the six M ≥ 6 aftershocks. Small uncertainties could remarkably affect the computation results. Therefore, we explore the coulomb stress change on the focal mechanism solutions given by USGS. We assume that each M ≥ 6 aftershock has a focal mechanism solution consistent with the regional tectonic stress field.

4. Coulomb Stress Changes Computed in the Focal Source

In this section, to further test whether the M 7.8 mainshock triggered subsequent aftershocks, we also examine the shear, normal, and coulomb stress changes in the focal source for each M ≥ 6 aftershock. The coulomb stress changes in the focal source more accurately test the stress triggering effect. Here, we present the USGS data to compute the coulomb stress change caused by the M 7.8 mainshock. As shown in Table 1, there is no obvious feature for the shear stress and normal stress calculation results. However, for the coulomb stress, the majority of the calculation results was positive, except for the result of the M 6.9 aftershock on the nodal plane: 183,76,85, which agrees well with the result in Section 3. For the other five M ≥ 6 aftershocks, the coulomb stress changes on both fault planes increased at varying degrees. We computed that the M 7.8 mainshock increased the coulomb stress in the focal source of five M ≥ 6 aftershocks by at least 0.1 bar. Although some of the coulomb stress changes are probably too small to cause subsequent aftershocks, these coulomb stress changes can promote failure if the fault is already near failure or an advance future earthquake occurrence if it is not. The results indicate that the majority of the M ≥ 6 aftershocks were triggered by the coulomb stress caused by the M 7.8 mainshock slip. Few M ≥ 6 aftershock occurred where the stress is computed to have dropped.

5. Discussion and Conclusions

Based on the seismic stress triggering theory and elastic dislocation theory, we compute the coulomb stress changes on two inferred rupture planes according to the focal mechanism solutions of each M ≥ 6 aftershock. The computation results show that the vast majority of the M ≥ 6 aftershocks occurred within areas of coulomb stress increases (>0.1 bar). Based on the focal mechanism solutions from USGS, we also computed the coulomb stress changes in the focal source of each aftershock. These possible coulomb stress changes are greater than the threshold of coulomb stress triggering (0.1 bar). As shown in Table 1, the majority of the coulomb stress calculation results is positive, except for the result for the M 6.9 aftershock on the nodal plane: 183,76,85. For the other five M ≥ 6 aftershocks, the coulomb stress changes on both fault planes increased at varying degrees. Among them, the maximum increase reached 7.564 bar and the minimum increase was 0.159 bar. Thus, the M ≥ 6 aftershock distribution is well explained by the seismic coulomb stress changes caused by the Ecuador M 7.8 mainshock.

In this study, the coulomb stress triggering among the Ecuador M 7.8 mainshock and subsequent M ≥ 6 aftershocks was investigated. It should be noted that there is much randomness and uncertainty as regards the time and location of an earthquake. Data and model uncertainties are still the biggest obstacle in determining quantitatively whether coulomb stress triggered the subsequent aftershocks. These uncertainties include the finite fault model, earthquake distribution, focal depth, and focal mechanism solutions. The nodal plane uncertainties add random noise to the stress change estimates and systematically reduce a potentially large coulomb index [11]. It is difficult to perfectly explain the complex geological phenomenon by only the coulomb stress triggering hypothesis. The original aim of this study is to provide a more reliable seismic hazard assessment in the Ecuador earthquake zone by considering the interaction of the earthquake sequences.

In this paper, it should be noted that we only adopted one finite fault model inverted from GSN broadband waveforms to compute the coulomb stress changes. However, the source model has much inherent uncertainty due to input data errors, model assumptions, and model estimation procedures [40]. Slip models using rectangle sub-faults artificially cause anomalous stress concentration along the edge of each sub-fault, which will affect the coulomb stress computation results. These uncertainties from the source model have yet to be fully studied. We attempted to approach source model uncertainties by considering multiple source models for the 2015 Nepal M 7.8 earthquake [41]. With the set of finite fault source models of the June 2000 M 5.9 Kleifarvatn earthquake, Woessner propagated correlated fault model uncertainties to determine uncertainties of coulomb failure stress change calculations [40]. In terms of the Ecuador M 7.8 earthquake, some limitations apply to our uncertainty and robustness analysis. The most important factor is the lack of reliable source models for this earthquake. Second, the fault model uncertainties should be regarded as a lower-bound estimate of the true uncertainties, as uncertainties in the modeling procedure are included via the global optimization procedure by Sudhaus [42].