Production of a homogeneous seismic catalog based on machine learning for northeast Egypt

1 Introduction

Earthquakes are natural disasters that significantly impact culture, politics, and health in the regions affected by them [1]. The magnitude of an earthquake, which is the energy that an earthquake releases [2], is a crucial variable for assessing the human impact of an earthquake [3]. Large-magnitude earthquakes can cause considerable damage [1].

Located in the northeastern region of the African continent, Egypt is bound by three major tectonic elements (see panel (a) in Figure 1): the African-Eurasian plate margin, Red Sea plate margin, and Levant transforming fault. Tectonic distortion is related to contact and relative movement of the Egypt’s crust below Egypt’s territories and its remote effects [4,5]. Egypt’s seismicity is low to moderate relative to that of North-West Africa, and much more so relative to the extreme seismicity of neighboring Hellenic and Cyprus subduction arcs and the Dead Sea transformation zone (Figure 1). Local moderate earthquakes, however, pose a major threat to the population, as demonstrated by the 1992 case in Cairo that caused 561 fatalities and by the historical 1847 event [5]. Regarding the transtensional stress regime triggered by the expanding Red Sea margin, most earthquakes are characterized by a normal fault with variable strike-slip components, which increases towards the edge of the Sinai subplate. A few cases, which are mainly inland, have reverse focal mechanisms [4].

Figure 1

(a) Regional tectonic setting of Egypt. (b) Historical seismicity of Egypt based on the compiled earthquake catalog of ref. [53]. (c) Significant historical and instrumental events affecting the area (adopted from ref. [27]). (d) Instrumental seismicity of Egypt based on the compiled ENSN catalog.

Earthquake catalogs contain data on the time of occurrence, hypocentric coordinates, magnitude of the event, and other related information about the earthquake [6]. These details are extracted from observations of various types of seismic waves, usually from heterogeneous seismological networks in the space-time domain. The data obtained from seismic catalogs provide valuable information for the study of seismotectonics and seismic hazard distribution in space and time [6]. However, accurate results require the use of a homogenized seismicity catalog. A crucial point in the process of homogenization is the identification of all possible earthquake occurrences by using a single target magnitude scale [7].

Due to the occurrence of the catastrophic earthquake, in Cairo in October 1992, the Egyptian National Seismological Network (ENSN) project was implemented with the aim to cover and report the entire Egyptian territory (Figure 1) seismic activities. The ENSN started operating in August 1997, and since then, accurate parameters of relevant seismic activity have been determined. Despite the establishment of ENSN, rapid changes in the method used for estimating the magnitude of seismic incidents, such as seismometers, and data reduction procedures have caused variances in the reported magnitude; thus, a consistent magnitude scale is not always available [8]. There have been several scales of magnitude utilized, and each is specified in terms of amplitudes measured with a specific type of seismometer over a very narrow spectral band; thus, variable measurement errors are manifested. For this reason, different types of magnitudes have been reported for the same earthquake, such as local magnitude [9], surface-wave magnitude [10], body-wave magnitude [11], and moment magnitude [12]. This will have tremendous effects on the precision and homogeneity of seismic hazard, as well as on the calculation of seismic deformation [13,14]. Various empirical relations have been utilized to overcome that problem. Different types of magnitudes are transferred and unified to the same scale, prior to inclusion in the catalog for seismicity and seismic hazard studies. The majority of researchers prefer to utilize the moment magnitude (M _w), as it is no suffer saturation for earthquakes of higher magnitude, as opposed to other types of magnitude. This is because the size of an earthquake is related to the ruptured field, average slip volume, and force necessary to overcome friction, which hold the rock together from offset by failures that are directly derived from the seismic moment (M ₀) [12].

Traditionally, magnitude conversation relationships are utilized to produce homogenized earthquake catalogs that implement several regression procedures, including simple linear regression (SLR), inverted standard regression (ISR), orthogonal regression (OR), and general orthogonal regression (GOR), to produce various magnitude conversion relations [15,16,17,18,19,20,21]. SLR assumes that either the independent variable is error-free, or its error order is significantly small in comparison with the dependent variable measurement error. In this regression type, vertical residues are minimized. By comparison, as the locations of dependent and independent variables are reversed, ISR is accomplished by reducing horizontal offsets to the best-fit rows [22]. A comparison among regression methods applied to earthquake magnitude conversions is given in refs [23,24].

Empirical conversions between different types of magnitude are typically utilized in Egyptian catalogs to homogenize the data for statistical analysis and hazard assessment. The standard least-square linear regression method is the most common method used to determine the relationship between different magnitude scales [25,26,27]. These empirical relationships may be improved if we understand the spatial clustering of earthquakes, as earthquakes are known to cluster in space and time due to stress redistribution in the Earth’s crust [28]. Hence, current research presents a new clustered linear regression approach, which is utilized to overcome the shortages of multiple regression. The methodology we describe in this paper uses a multistage machine-learning approach to deal with traditional magnitude conversions methods’ limitations. Typical multistage learners combine two or more inferential approaches so that more complex learning problems can be solved. Our approach is a hybrid combination of unsupervised and supervised machine-learning methods. The three cluster algorithms involved in most studies are the partitional algorithm, hierarchical algorithm, and density-based algorithm [29,30]. Because of its applicability and reliability, the partitional algorithm is most commonly used. As with most common cluster techniques, k-means (hard clustering) and Fuzzy C-means (soft clustering) are categorized in terms of the partitional algorithm. k-means vary from Fuzzy C-means in that a point belongs only to one cluster in k-means, whereas the point belongs to a cluster to some extent depending on its Fuzzy Membership, in Fuzzy C-means.

In this research, k-means was utilized as a convenience to prevent the additional computational costs that are incurred by fuzzy inclusion in fuzzy C-means. k-means clustering is used before the application of a multiple regression method which is a supervised analysis procedure that can be described as a model of relationships between a target continuous variable and a collection of other input variables. This multistage approach has been adopted with the objective to improve the accuracy of learned models.

2 Methods and data

This section explains the proposed multistage seismological and machine-learning approaches that are used to predict the conversion relation between reported local magnitude M _L and estimated moment magnitude M _w with the intention of enhancing our learned models and predictive accuracy. In the first phase of our strategy, we acquired the frequency-independent, long-period spectral level Ω₀ below the corner frequency of a displacement spectrum [31,32]. The long-period spectral amplitude of the source spectrum can be related to seismic moment M ₀ [31] from which moment magnitude M _w can be computed [12].

For seismic moment estimation, the earthquake signals from the three component seismometers were inverted by assuming a Brune’s model for the source and connecting it to a rupture area with a slip value. According to Brune’s model [31,32], the displacement spectrum specifies that the frequency response of the seismic signal is considered to be flat until it reaches the corner frequency ω; then, the amplitude of the displacement spectrum drops to ω ⁻². By fitting such a model to the displacement spectrum, the low-frequency plateau W can be estimated and used directly in the calculation of the seismic moment because it relates to the area of seismic slip on the fracture plane by adapting the model to the displacement spectrum. Geometric spreading R, seismic wave velocity c, rock density ρ, and a radiation pattern imposed by the seismic source moment tensor F _c were included in the calculation to estimate the seismic moment as given in the following equation:

where μ is the rock’s shear modulus, d is the average fault displacement, and A is the fault area.

For P-waves W _c = W _p. According to ref. [33], F = 0.52 for P-waves could be adopted as an approximate value for F _c if the radiation pattern of the moment tensor is not known. With seismic moment, measured in Nm, the moment magnitude is estimated following the ref. [12] formula as:

In the second step of our strategy, unsupervised machine-learning data clustering utilizing k-means algorithm is carried out to find a data partition such that points in the same subgroup are close to the similar Euclidean distance ( d ( x , y ) = ∥ x − y ∥ 2 ) between two earthquakes x and y [34]. The k-means algorithm is a heuristic technique that often produces effective clustering outcomes. In that algorithm, the global minimizing partition S* is usually intractable because the group of partitions is very large. The algorithm specifies an initial cluster mean µ ⁽¹⁾ for each i = 1,…, N. Then, given a set of cluster means (t), the k-means algorithm searches for a partition S(t) of the earthquake epicenters such that [34]:

The algorithm then updates the group of epicenters for each i = 1, N and continues to iterate until the partition ceases shifts [35]. The within-cluster sum of squares (WCSS) is used as an indicator of the performance of clustering efficiency (lower value for WCSS is better).

For the final phase of our strategy, linear regressions between the local magnitude and moment magnitude ranges will be performed. Using the most consistent linear relationship, we convert the local magnitude M _L into moment magnitude M _w, for each available pair of magnitudes. In this approach, the relationship between the predictor variable and the response variable with SLR can be described as a one-to-one relationship, while multiple linear regression is a many-to-one relationship. For multiple linear regression, the relationship between X and Y is expressed as:

Similar to SLR, the regression coefficients are estimated using the least-square approach so that β ₀, β ₁, …, β _p is selected to minimize the residual sum of squares (RSS). The residual is the difference between the expected magnitude value and the true magnitude value; y – yˆ. yˆ can be replaced with β ₀ + β ₁ x and the residual sum of the squares can then be defined as:

We used earthquake data collected from the ENSN to test the efficiency of our proposed method and to establish a reliable relation between the reported local magnitude (M _L) and the moment magnitude (M _w). To estimate M _w for the region between 27–31°N and 31–36°E (see panel (a) in Figure 2), displacement spectra of 138 earthquakes that occurred from 1998 and 2020 were used. Such events occurred in various seismic areas of Egypt, including the Red Sea, the Aqaba, and Suez Gulfs, and localized seismicity in the northern part of the Eastern Desert. In addition, only stations with epicentral distances greater than 10 km were considered in order to eliminate near-field effects and obtain a point-source approximation [36]. The earthquake’s hypocentral parameters were taken from the annual ENSN bulletin. Selected waveforms were captured by broadband stations with good events-stations azimuth coverage and display a higher signal-to-noise ratio and clear P-wave phase. Low-pass filters with corner frequencies of 35 Hz were used to prevent an aliasing effect. The analyzed seismograms were corrected using an instrumental response for all stations. In addition, the local site effect was ignored since all stations were in hard rock. The source parameters were estimated by utilizing the Fourier spectra of the pure P-phase displacement spectra of the integrated velocity data and after applying instrumental and baseline corrections. Before marking the corner frequency and long-period spectral level, attenuation correction was applied, assuming the well-accepted Q-value for region as suggested by ref. [37].

Figure 2

(a) Selected events utilized to estimate M _w. (b–f) Examples of earthquake spectra for some ENSN stations.

Next, the source term of displacement was processed based on the refs [31,32] model, which utilizes circular faults for small seismic events. Rock densities and P-wave velocities were taken from the ref. [38] velocity model. Following ref. [39] correction of the radiation pattern for P-waves, a 0.52 value was utilized. Following the procedure of ref. [40], the corner frequency (f _c) and spectral level (Ω₀) were determined for each P-wave source spectra and used to calculate the seismic moment (M ₀) by fitting each P-wave displacement spectrum to the theoretical displacement spectra model as expressed by equation (1). The nonlinear Levenberg–Marquard fitting algorithm [41] was used to minimize the discrepancy from the theoretical to the observed range of displacement spectra. The approach will match the spectra effectively without having to compensate for the loss of the frequency of the corner due to the flat component. In Figures 2 and 3, the best matches for the observed displacement spectra to the theoretical spectra at a number of stations are plotted. Following the definition given in equation (6), the average seismic moment <M ₀> was calculated from the average logarithmic values obtained at different stations as proposed by ref. [42].

where N is the number of stations used, and M ₀ is the seismic moment determined from equation (1) for the ith event record. Then, the average moment of magnitude (M _w) is calculated using equation (2) for that event. Table 1 shows all variables in the description of the analysis.

Figure 3

(a–f) Another examples of earthquake spectra for some ENSN stations. (a) KOT spectrum, (b) NAT spectrum, (c) RYAN spectrum, (d) TAMR spectrum, (e) ZAF spectrum, and (f) ZNM spectrum.

Most studies developed scaling relationships between various magnitudes (see for example, refs [25,26]), but did not consider the spatial distribution of the events that were utilized in the analysis. This sometimes causes large fitting errors. Therefore, we need to cluster these earthquake data into active seismic regions as outlined in the second step of our proposed approach before establishing the scaling relationships. In this step, every earthquake epicenter, with coordinates being their latitude and longitude, is considered to be a point source in R ². However, since the Earth is not flat, that assumption is not accurate. It is important to view the latitude and longitude values in R ³ as a modification of spherical coordinates.

Based on the Euclidean distance of our k-means algorithm, we need to transform our database into a 3D Euclidean coordinate. A convenient way to achieve this transformation is to convert latitude and longitude values into spherical coordinates, then into Euclidean coordinates. Recall that a spherical coordinate in R ³ is a triple (r, θ, ϕ), where r is the distance from the origin, θ is the radial angle in the xy-plane from the x-axis, and ϕ is the angle from the z-axis. In our earthquake data, the longitude is already the appropriate θ value and the ϕ value (in degrees) is simply 90° minus the latitude.

Based on the data shown in Figure 2(a) and Table 1, the spatial distribution of the data utilized to estimate M _w is sparse and does not have enough events to obtain good clustering results. To solve the problem, we incorporated the data published in refs [17,18,43,44] in our database for utilization in the second and third phases of our hybrid clustered-regression approach. Spatial distribution of our collected data and incorporated data is given in Figure 4. A manual check for duplicates was performed. Table 2 provides descriptive statistics for our final database attributes.

Figure 4

Final spatial distribution of seismic events utilized for machine-learning clustering and regression analyses. M _w data are incorporated from the current study and the data published in refs [4,26,44,45].

The application of the k-means algorithm requires the number of clusters k to be known. Since our seismic data are an unlabeled dataset, finding the optimal value of k is a challenge. In this case, it is not certain how many clusters are required. Therefore, we used a heuristic approach called the Elbow method [18] to identify the optimal k value from the relation between the number of clusters k and the WCSS. For the elbow method, the distance from each earthquake data point and the cluster centroid to which it was assigned is determined by the within-cluster sum-of-squares (inertia) using the following equation:

where C _k represents cluster k, x _i is a data point, µ _k is the centroid for cluster k.

A plot of k vs the sum of squares reveals an elbow point in which the inertia starts descending much more slowly and can be selected as the best value for k. Figure 5(a) displays the approximate value of the sum of the k-means squared distance for different values of k. In the plot we can see a strong elbow at k = 5, which is the best value for k. However, the elbow method is limited, since the solution is heuristic, so we used the Silhouette method [18]. That method measures how well each data-point x _i “fits” into its assigned cluster and how poorly it fits into other clusters. This is a different way of looking at the same objective. Denote a _xi as the average distance from x _i to all other points within its own cluster k; the lower the value, the better. On the other hand, b _xi is the minimum average distance from x _i to points in a different cluster, which has been minimized for all clusters. That is, compute separately for each cluster the average distance from x _i to the points within that cluster, and then take the minimum. The silhouette s(x _i ) is defined as ref. [18]:

Figure 5

Comparison between WCSS (a) and Silhouette score (b) methods for finding the optimal number of clusters k-value for k-means algorithm.

Figure 5b shows the average silhouette score compared with results obtained by the elbow method for k-means clustering of the sample data, while Figure 6a and b show the results of the silhouette analysis with the corresponding scores computed for every earthquake data point in every cluster. The obtained silhouette scores depicted in Figure 6(a–d) and (e–h) confirmed the number of clusters that were obtained from the elbow method. Accordingly, we chose five as the optimal value of k to feed into the k-means clustering algorithm.

Figure 6

Silhouette analysis for k-means clustering on sample data with number of clusters ranging from two clusters to five clusters (a–d) and from six clusters to nine clusters (e–h).

In earthquake catalogs, M _w is not commonly considered the main magnitude, so understanding how certain magnitudes respond to M _w in each area is a matter of considerable concern due to seismic hazard estimation. Usually regression analysis is utilized for that purpose. However, Linear Regression Analysis consists of more than just fitting a linear line through a cloud of data points [18]. Before we begin supervised machine-learning regression analysis and modeling, it is recommended that we perform an exploratory data analysis and plot the estimated and collected magnitude data, as shown in Figure 7.

Figure 7

Comparison of the statistical distributions of local magnitude (M _L) with the moment magnitude (M _w) in terms of scatter plot (a), histogram plot (b), and box-plot (c).

By examining these initial plots, we can quickly determine that the data have a strong linear relationship as shown in Figure 7(a). From the histogram plots shown in Figure 7(b), we can clearly observe that both local and moment magnitudes are normally distributed. An example of this pattern is that the data are unimodal, or the data have a single mode, which is marked by the curve’s peak. On the other hand, the local magnitude values are skewed to the right, as the histogram has a middle peak that steadily tapers toward the right side of the graph. This dataset is unimodal, as the mode closest to the left of the graph is smaller than the mean or the median. The mean for the local magnitude data is located on the right side of the diagram and is higher than the median or the norm. This graph shows that there are few data points that are larger than the mode, which may be representative of outliers in our database. Finally, box-plots are used to test the presence of outliers in our data. Box-plots (Figure 7c) are used to show the overall response patterns for each magnitude type. They provide a useful way to visualize the magnitude value ranges.

Figure 7(c) shows a variety of different box plot shapes and positions were obtained for each magnitude type. It is clear that both M _L and M _w suffer from the existence of outliers. The observed outliers exerted a large influence on the overall outcome of the regression model parameter’s estimate.

3 Results and discussions

Through an analysis of earthquake source parameters, local tectonics and basic earthquake mechanics can be understood. To improve scaling between moment magnitude and local magnitude for earthquakes in the northeastern Egypt, the moment magnitude (M _w) estimates were extended to lower magnitudes by using the spectral analysis of P-waves from vertical component seismograms. This dataset has been processed, analyzed, and interpreted to study the source parameters of local earthquakes occurring around northeastern Egypt. Only events from at least three ENSN stations at different distances with reliable trace amplitude readings were selected. In this way, a total of 138 events with reliable M _L values (Figure 2a) were chosen for M _w determination. The long-period spectral amplitude of the source spectrum, which could be related to the seismic moment, was utilized to compute M _w following the definition of ref. [16]. Figures 2 and 3 display examples of P-wave Fast Fourier Transform spectra fitted by the smoothed spectral [17] model at several ENSN stations.

Displacement of the P-wave portion is considered to be an approximate source time function and could be integrated to obtain scalar moments at each station. Multiple records are combined to give a good estimation of the moment magnitude without correcting for the source mechanism [18]. Table 1 summarizes the average value of the obtained source parameters for small earthquakes in northeast Egypt.

Values of the M _w calculated using the method described in this paper were compared with those obtained from automatic moment tensor solutions and published in refs [43,44,45]. The average deviation between the determined M _w and that of the aforementioned studies is ±0.085, meaning that, on average, our moment magnitude was slightly overestimated.

After finalizing the first phase and estimating the moment magnitude M _w, the datasets were utilized to compare the partitions with the distribution of recognized active sources and were applied to areas that may not be adequately captured by existing source models. In addition, k-means cluster analysis was applied to the hypocentral distribution of observed earthquakes in the mapped region. Since the moment magnitudes of the utilized events were lower than six (M _w < 6.0), we can consider them to be point-source attributes. The obtained clusters comprise sets of earthquake data that are similar to other data in the same group, but dissimilar from data in other groups, according to the Euclidean distance. Based on results of the Elbow and Silhouette methods shown in Figures 5 and 6, we picked five (k = 5) to be the final number of clusters, as it was considered to be the optimum number of clusters for which the highest average silhouette score was taken from the value of the silhouette of each cluster.

Our implementation of the k-means algorithm began by assigning each seismic occurrence randomly to a cluster, and then the mean center of each cluster was determined. At this point, the Euclidean distance between each occurrence and the initial clusters was determined and reassigned to a new cluster centered on the closest mean center. Then, the mean centers were recalculated, which proceeds until the cluster elements stop changing. Using the optimal number (k) as suggested by the Elbow and Silhouette score methods (Figures 5 and 6), we were able to group the earthquakes from the dataset into five different clusters. We grouped the seismic activities according to their geographic location. The latitude and longitude values were converted into Universal Transverse Mercator coordinates for correct estimation of the Euclidean distance and cluster centers or cluster centroid identification. Distribution of delineated clusters is shown in Figure 8(a) and the frequency of each cluster is given in Figure 8(b).

Figure 8

Northeastern Egypt cluster zones (a), frequency percentage derived from the k-means clustering algorithm (b), and Davies-Bouldin score verses number of clusters in dataset (c).

The Davies-Bouldin index [46] was used to test the clustering performance and to evaluate the k-means model, where a lower Davies-Bouldin (DB) index was correlated to a model with a better separation between the clusters. This index signifies the average “similarity” between clusters, where the similarity is a measure that compares the distance between clusters with the size of the clusters. As shown in Figure 8(c), the lowest (DB) index obtained corresponded to the five clusters, which confirmed the appropriateness of the utilized value of number of clusters for the k-means model and indicated a better partition of seismic activity zones affecting the mapped area.

The obtained clusters are in accord with various seismotectonic studies of Egypt (see for example, refs [27,38,47,48]). The clusters that were found may be linked to three main seismotectonic areas. The first region is the Gulf Suez, which comprises the second (cluster ₁) and fourth (cluster ₃) clusters. The area is the location of the March 31, 1969, Shadwan Island event with a surface-wave magnitude M _S = 6.8 in the mouth of the Gulf of Suez. The majority of the focal mechanism solutions of reported events indicate that the current tectonic activity had a normal fault mechanism in this region. Events of this region represent 23% of our database. The second region, which is known as the Gulf Aqaba, is an area of extreme activity situated along the main boundaries of the tectonic plate. Its major earthquakes were triggered by the displacement of this transform boundary. In this area, the largest recorded earthquake (M _w = 7.2) occurred on November 22, 1995. This area’s events represent 42% of the data in our database, and they were part of the first (cluster ₀) and fifth clusters (cluster ₄). The last area contains the largest percentage of events in our database and is known as the Northern Red Sea source region. The fault plane solutions have a normal fault mechanism with a minor strike-slip component that is compatible with the Gulf’s extensional rhomb-shaped grabens system. The number of seismic activities available in each region varies considerably. The minimum number is observed in the second clustered zone with 41 delineated events, while the maximum number of 136 delineated events was obtained for the fourth zone. Table 3 provides concise statistics for each cluster.

The last stage of our multistrategy approach consists of a supervised learning process following previous clustering-based partitioning of the magnitude data. The M _w versus M _L relationships for the delineated five cluster zones were determined.

The most common approach to perform homogenization of different magnitude scales into a single magnitude scale and to compare various measures of magnitude is through a regression analysis that utilizes the least-squares method [49]. However, its basic assumptions are rarely satisfied in practice [15]. Hence, refs [15,22] recommend that if the absolute variance of the data values is not understood, which is always the case in magnitude estimates, various regression forms should be performed and compared. Therefore, the M _L magnitude scale can be converted into an M _w magnitude scale by building empirical relationships between the correlated earthquake events.

Uncertainties related to various forms of magnitude play a critical role during transformation of magnitudes into a single magnitude scale. The utilized dataset also suffers from the presence of outliers (see box-plots in Figure 7c). Possible causes of observed outliers in magnitude values were not fully identified. The estimated local magnitude values that were extracted from the ENSN bulletin essentially depend on amplitude pickings, which occur over a long period of time at varying levels of seismic noise at a given station, with different levels of noise occurring at different stations. Moreover, the moment magnitude estimate also relies upon the impulse response of the medium, which was never perfectly evaluated and may be highly sparse from one seismic event to another, depending on the event’s location versus the seismic stations and earthquake source depth.

Although there is no clear description or definition of the outlier, outliers are often considered as observations that do not fit the trend of the other observations. It is not usually a concern if the outlier is merely an unusual finding made from the tail of a regular distribution; however, if the outlier comes from a non-usual calculation error or any other breach of typical ordinary least-square rules, then it undermines the integrity of the regression if a non-robust regression method is used [50].

To study the effect of the observed outliers and to develop a robust regression method, the RANdom SAmple Consensust (RANSAC) method [51] was applied to our study. RANSAC is a non-deterministic algorithm that partitions the whole input sample data into a collection of noise-prone inliers and outliers, e.g., due to incorrect measurements or invalid data hypotheses. Therefore, the resulting model is derived from established inliers only. The algorithm creates and validates a linear model based on the Minimum Least-Square (LMS) approach by extracting noisy samples (i.e., outliers), choosing the appropriate features (i.e., descriptors), deriving an LMS model from the training set samples, and finally, predicting the behavior of the test set samples while using the principle of applicability domain [51].

To develop the M _w versus M _L empirical relation, the RANSAC algorithm is utilized to estimate regression coefficients of the linear model by randomly sampling the observed earthquake magnitude data in each cluster to find the optimal fit result and implementing a voting scheme [51]. The implementation of this voting scheme is based on two assumptions: the noisy features will not vote consistently for any single model (few outliers), and there are enough features to agree on a good model (few missing data). Regression fitting is performed in two stepwise manners that are iteratively repeated. In the first step, a subset containing minimal magnitude data is randomly selected from the input earthquake dataset. A linear LMS fitting model and the corresponding model coefficient parameters are computed using only the elements of this magnitude subset. In the second step, the algorithm checks which elements of the entire magnitude dataset are consistent with the model instantiated by the estimated model parameters that were obtained from the first step. A data element will be considered an outlier if it does not fit the fitting model instantiated by the set of estimated model parameters within some error threshold that defines the maximum deviation attributable to the effect of noise [52].

The RANSAC regression relationships established between M _w and M _L for the various cluster zones are shown as blue regression lines in Figure 9a–f and Table 4. To check the rationality of the derived conversion relations between magnitude types, we applied ordinary least squares (OLS) and the GOR (see the red and green regression lines in Figure 9a–f). We assumed that the ratio between the variable variance (η), which was not known for our dataset, was equal to one, based on the regulations described in ref. [22], assuming that the error variances of different magnitudes were approximately equal.

Figure 9

(a) M _w M _L relationships for northeastern Egypt cluster zones derived from SLR, GOR, and RANSAC regression algorithms for all clusters, cluster 0 and cluster 1. (b) Same M _w M _L relationships for cluster 2, cluster 3, and cluster 4.

When comparing the different regression metrics, we can see that the RANSAC algorithm can boost our model scores when the dataset includes a significant number of outliers. In this specific case, all the regression metrics scores improved significantly. The difference is particularly observable in the R ² score, which is a statistical indicator reflecting the percentage of deviation for a dependent variable described by an independent variable. It reached almost 99% with the RANSAC regression model. The good agreement between the observed and estimated magnitudes provides greater confirmation of the results and could be easily used to validate the accuracy of the proposed relations for the earthquake early warning system in the northeastern portion of Egypt. It was also observed that the relations derived from this RANSAC technique had significantly lower errors in their regression parameters than those of the corresponding relations established from other techniques, as shown in Table 4, for all of the computed regressions.

However, the R ² score is not always enough to determine which model is a better fit for the results. Additional error analyses are required, given the distribution of the residuals along the fitted line. We added more measurements, including the scattering residual analysis, to complete the predictive analysis on the goodness-of-fit of our regression models. Residuals will show a random scatter, as was already mentioned. Prediction error analysis for the three regression models is illustrated in Figure 10a–c.

Figure 10

(a) Residuals for each M _w M _L relationship for northeastern Egypt cluster zones derived from SLR, GOR, and RANSAC regression algorithms for all clusters and first zone. (b) Same Residuals for each M _w M _L relationship for second and third cluster zones. (c) Same Residuals for each M _w M _L relationship for fourth and fifth cluster zones.

It is very clear that most of the points in the residual plot are distributed randomly along the horizontal axis (shown with a red line) for RANSAC and GOR; hence, those two regression models are suitable for predicting the conversion relations between the two types of earthquake magnitude. This also suggests that for any possible regression analysis, grouping the used collection by considering the spatial distribution of earthquakes will vastly enhance the creation of a suitable scaling relationship that can be applied to homogenize the magnitude of the earthquake.

Finally, the proposed relationships were compared to those provided by other authors in order to analyze the validity of the established relations. The results of this study are generally consistent with the observations reported by other authors, although some differences in value were revealed. Meanwhile, there were some differences in the values among the relationships. As geographic classification of seismic events was not discussed by other researchers, this could preclude the direct use of existing scaling relationships by other writers in a similar area such as the Gulf of Aqaba. This should be used to highlight the fact that the scaling relationships in region cannot be used explicitly with that of other regions and need to be investigated further.

4 Conclusion

The purpose behind the conversion of different magnitude scales into a single moment magnitude M _w lies in the fact that it remains unsaturated throughout the entire magnitude scale and achieves an accuracy that is two to three times higher than that of other magnitude types. Linear regression models are heavily impacted by the presence of outliers, which can come from extreme values of noise or erroneous measurements. The RANSAC algorithm can be used to remove sets of points that do not follow the dominant pattern of the data in the models, and as a robust machine-learning algorithm, RANSAC improves the model’s performance by estimating the parameters with a high degree of accuracy, even when a significant number of outliers are present in the dataset. RANSAC usually performs poorly when the number of inliers in the dataset is less than 50%. In addition, when the noise threshold is too small, the estimated parameters tend to be unstable. Finally, this algorithm removes data from the model, but loss of data should be avoided when developing models.