# Stochastic Modeling of the Al Hoseima (Morocco) Aftershock Sequences of 1994, 2004 and 2016
## 5 Magnitude-Frequency Relationships  
The Gutenberg-Richter recurrence relationship [54; 55] is a frequently used approach for quantifying seismic activity in each region and has been shown to apply aftershock frequency-magnitude data. The equation is a reasonable approximation of the frequency-magnitude statistics that describe the correlation between earthquake occurrence frequency and magnitude  
\[log_{10}\ N\ (\geq m)=a-bm\ \ \ \ \ \ \ \ \ m\geq m_{c} \tag{1}\]  
Figure 7: Stress regime computed from focal mechanisms. **(A)** For the 1994 mainshock, **(B)** for the 1994 aftershock sequence, **(C)** for the 2004 aftershock sequence, **(D)** for the ‘eastern cluster’ considered by [34], **(E)** for the 2016 ‘main cluster’ and **(F)** for the 2016 ‘secondary cluster’. where \(N\) (\(\geq\)\(m\)) is the number of events with magnitudes greater than or equal to m. For the estimation of both the a and \(b\)-values, it is widely recommended to use a complete dataset for all ranges of magnitude. The threshold magnitude parameter \(m_{c}\), is typically determined using one of two methods: a network-based [56; 57] or catalog-based approach [56; 57]. The first group uses the day-to-night ratio to calculate the earthquake frequency, if noise reduces the detection threshold at night [58; 59]. The second set of approaches assumes that earthquake production is self-similar, allowing us to use a power law or the Gutenberg-Richter relationship to construct earthquake frequency-magnitude distributions. The most applied methods in this group reviewed by [60; 61] are the maximum curvature (Maxc) method [62], entire magnitude range (EMR) method [60; 63], median-based analysis of the segment slope (Mbase) [64], determination of \(b\)-value instability [65] and goodness-of-fit test (Gft) by [62], which was later modified by [66]. The threshold magnitude for each of the three sequences under consideration was thoroughly estimated in the current study. Although both the maximum curvature (Maxc), Gft, Gft at levels of 5% and 10% (Gft5% and Gft10%), and modified goodness of fit (mGft) methods were first investigated, the maximum curvature method produced better results, i.e., a better fit (Table 1). Threshold magnitude values of 2.8, 3.4 and 2.0 were obtained for the 1994, 2004 and 2016 aftershock sequences, respectively. The threshold value fluctuated over time, which was also studied. This temporal fluctuation is depicted in Figure 8 and was computed using a sliding-window method. Prior to each new computation, a window of 20 events was shifted by five events [67]. The middle time of the associated window is supplied to each new threshold magnitude value. The window length was chosen as the balance between the need for temporal resolution and smoothness. Multiple tests were conducted first, changing the number of occurrences per window and the shift; however, neither aspect had a significant influence on the definition. The maximum likelihood approach was used and is considered one of the most reliable approaches for computing the \(b\)-value among the available methods. The estimator by [68], given below, was used to calculate the \(b\)-value:  
\[b=\frac{log_{10}e}{\langle m\rangle-m_{c}} \tag{2}\]  
\(\langle m\rangle\) being the average value of the magnitude. An estimate of the standard deviation \(\sigma_{b}\) is obtained using the [69] relation, given as follows:  
\[\sigma_{b}=2.30\,b^{2}\sqrt{\frac{\sum_{i=1}^{N}\left(m_{i}-m\right)^{2}}{N(N -1)}} \tag{3}\]  
According to [69; 70], even when the \(b\)-value varies in time and/or space, this relationship provides a reliable approximation of \(b\)-uncertainty. When utilizing rounded magnitudes, the estimator by [68] is inaccurate but agrees with the modified distribution of [71].