Introduction

Earthquakes that occur offshore and along the eastern shore of the Red Sea are caused by tectonic processes acting in the vicinity of the Arabian Shield. Several of these seismic source zones are potential sources of seismic hazard for metropolitan areas. Because of the geodynamic processes in the region, the territory of the Arabian Plate is being dragged northward in intensive movements, resulting in complex local tectonic settings and causing this territory to trigger moderate seismicity and form complicated geological structures in many areas (Abdelfattah et al. 2017, 2020). This seismotectonic setting strongly indicates that Saudi Arabia has experienced considerable earthquakes in the past and will experience them in the future, in particular along the coastline both onshore and offshore the Red Sea. Historical documents indicate that an earthquake of magnitude 7 occurred in 1068 in the northwestern Arabian Peninsula (Ambraseys et al. 2005). During the last century, earthquakes of Mw 5.1 in 2004 (Aldamegh et al. 2009), Mw 5.4 in 2009 (International Seismological Centre), Mw 3.9 in 2009 (Aldamegh et al. 2010), Mw 4.6 in 2014 (Abdelfattah et al. 2017), and Mw 4.0 in 2017 (Abdelfattah et al. 2020) have occurred.

The seismic hazard assessment over Saudi Arabia Peninsula is limited and focuses mainly on the probabilistic seismic hazard analysis (PSHA). The PSHA was quantified by Al-Amri (2013), and Zahran et al. (2015, 2016) using GMPEs of different areas outside Saudi Arabia. Moreover, the deterministic seismic hazard (DSH) was analyzed in some local areas located along the eastern coastline of the Red Sea (Almadani et al. 2015). However, the DSH does not well cover the whole area of the western Saudi Arabia peninsula. Owing to no ground motion distance attenuation relationships developed specifically for Saudi Arabia, several different ground motion prediction equations (GMPEs) were used in regions where no strong ground motion records were available (e.g., Thenhaus et al. 1989; Fukushima and Tanaka 1990; Sadigh et al. 1997; Skarlatoudis et al. 2004; Pankow and Pechmann 2004; Ambraseys et al. 2005; Atkinson and Boore 2006; Zhao et al. 2006; Akkar et al. 2014). Recently, a regional GMPE was established by Kiuchi et al. (2019) for western Saudi Arabia. Abdelwahed et al. (2020) used predictive attenuation relationships from different regions to identify appropriate relationships capable of predicting ground motions in the region. Besides, Abdelfattah et al. (2020) developed an application based on an empirical Green’s function summation technique to simulate seismic ground motion records in the Jazan Region of Saudi Arabia using the waveform data and source parameters derived from the spectral amplitudes of an Mw 3.4 event. The estimated peak ground acceleration (PGA) values are consistent with the GMPEs derived for California, Japan, and western Saudi Arabia. The results appear to be capable of predicting ground motions for future earthquakes; therefore, GMPEs in the region can be optimized by incorporating the effects of local geological conditions.

The southernmost part of the Arabian Shield is located near a triple junction that connects the axial rift in the southern Red Sea, the rifting system in the Afar Region, and the rifting process in the Gulf of Aden. The deformation in this region is attributed to the relative movement between the Arabian and African plates and the active rifting processes in the region. Even though this region is not considered to have a high rate or size of seismicity, earthquake hazard prevention and mitigation need to be seriously considered; this region is of high interest because the government plans to erect some dams in the area. Also, this region is home to one of the major cities in the territory of Saudi Arabia where most of the mega projects exploit green technologies to produce electricity. Together with the population, this has motivated many researchers to assess the earthquake hazards and the related earthquake engineering functions in this region. For more precise seismic hazard assessments, a local predictive equation of earthquake ground motions needs to consider the local characteristics of the source effects, receiver-to-path effects, and local geological conditions in a specific region. Accordingly, it is crucial to define the attenuation relationships to more accurately predict GMPEs, specifically for urban development plans and earthquake hazard assessments in specific regions, which will be essential to assess earthquake hazards for fast-economic growth regions such as the Jazan Region. The Saudi Geological Survey records provide an opportunity to access records at sites that will provide an important reference to establish a regional GMPE.

Data

Waveform data were collected at broadband seismic stations, shown in Fig. 1, that belong to the seismographic network operated by the Saudi Geological Survey. This network consists of two types of velocity seismometers, STS-2 and Trillium 40, which have flat responses approximately ranging from 0.008 to 50 Hz and from 0.025 to 50 Hz, respectively.

Fig. 1
figure 1

Map showing the tectonic plate boundaries near the epicentral area with an inset showing the spatial distribution of the epicenters in the period from 2002 to 2017 in the epicentral area and its vicinities (light brown polygons indicate the wide distribution of ancient Cenozoic basaltic fields called Harrats)

Full size image

Figure 1 shows the geographical locations of the seismic stations and the distribution of the earthquake epicenters used in the present analysis. A plot showing the magnitude as a function of the distance is shown in Fig. 2. Note that the local magnitude is the only magnitude estimate available for most events, except for a few earthquakes. Following the method of Gupta et al. (2017), the signal-to-noise ratio between the S waves and the noise window prior to the P wave arrivals was examined to select the appropriate records for this study. The instrumental response was deconvolved from the selected seismograms, and then the records were differentiated into acceleration records. A total of 769 records with a good signal-to-noise ratio, related to 99 earthquakes and 150 stations, were selected to update the GMPEs for the region over the magnitude range from 4.1 to 6.9. Figure 3 shows example velocity records of an ML 4.9 earthquake that occurred on August 03, 2014, at 20:21:11.26 UTC. Figure 4 shows the distribution of the focal mechanism solutions that were available for a few events in the region, which indicate normal to strike-slip focal mechanisms. The focal depths of the earthquakes were generally shallower than 10 km. An extensional tectonic regime is associated with the earthquakes that occurred along the axial rift of the Red Sea, while a shear tectonic regime is associated with the earthquakes that occurred along the reactivation of the inherited conjugate fault system in the region (Abdelfattah et al. 2020).

Fig. 2
figure 2

Distribution of earthquake magnitudes as a function of the epicentral distance (left) and histogram of the records and corresponding magnitude classes for the data set used in the present analysis (right). The data set comprises 72 earthquakes having a local magnitude range from 2.0 to 5.1 recorded by seismic stations distributed over a distance range of 5–200 km

Full size image
Fig. 3
figure 3

Velocity seismograms of an event recorded at different stations for an earthquake that occurred on August 08, 2014, at 20:21:11.26 UTC (red lines show the picking marks for the P and S wave arrival times)

Full size image
Fig. 4
figure 4

Distribution of earthquake ground motions (open blue squares) and the corresponding predictive attenuation relations (solid red lines) for the mean values of PGV (a) and PGA (b) recorded by the horizontal components of the seismic stations listed in Table 1

Full size image

Method

The decay of ground motion amplitudes can be expressed as an exponential function that represents the effects of the geometric spreading and anisotropic attenuation coefficients in the form

$$ A\sim {r}^n{\exp}^{\upgamma r}, $$
(1)

where A is the amplitude of the ground motion, r is the hypocenter or epicenter distance, n is the geometrical spreading coefficient, and γ is the anelastic attenuation coefficient. Taking the logarithm of both sides of Eq. (1), we find

$$ {\log}_{10}(A)=c-\upgamma r-n{\log}_{10}(r), $$
(2)

where c is a constant. The geometric spreading coefficient is equal to 1 for the propagation properties that affect the transmission of seismic waves along the source-to-receiver paths. Taking into consideration the variation in the ground motion amplitudes with magnitude, the attenuation relationship to predict the earthquake ground motion can be written as

$$ {\log}_{10}(A)=a+b\cdotp {M}_L-c\cdotp {\log}_{10}(r)-d\cdotp r, $$
(3)

where A is either PGA or the peak ground velocity (PGV) and ML is the local magnitude. To solve Eq. (3), the scaling coefficients, a, b, c, and d, can be determined using regression analysis. Accordingly, Eq. (3) can be rewritten in matrix form as

$$ \left|\begin{array}{c}{\mathit{\log}}_{10}\left({A}_{11}\right)\\ {}{\mathit{\log}}_{10}\left({A}_{12}\right)\\ {}\vdots \\ {}{\mathit{\log}}_{10}\left({A}_{mn}\right)\\ {}{\mathit{\log}}_{10}\left({A}_{21}\right)\\ {}{\mathit{\log}}_{10}\left({A}_{22}\right)\\ {}\vdots \\ {}{\mathit{\log}}_{10}\left({A}_{mn}\right)\\ {}\vdots \\ {}\vdots \end{array}\right|=\left|\begin{array}{cccc}1& {M}_{L_1}& {\mathit{\log}}_{10}\left({r}_{11}\right)& {r}_{11}\\ {}1& {M}_{L_1}& {\mathit{\log}}_{10}\left({r}_{12}\right)& {r}_{12}\\ {}\vdots & \vdots & \vdots & \vdots \\ {}1& {M}_{L_m}& {\mathit{\log}}_{10}\left({r}_{mn}\right)& {r}_{mn}\\ {}1& {M}_{L_2}& {\mathit{\log}}_{10}\left({r}_{21}\right)& {r}_{21}\\ {}1& {M}_{L_2}& {\mathit{\log}}_{10}\left({r}_{22}\right)& {r}_{22}\\ {}1& \vdots & \vdots & \vdots \\ {}1& {M}_{L_m}& {\mathit{\log}}_{10}\left({r}_{mn}\right)& {r}_{mn}\\ {}\vdots & \vdots & \vdots & \vdots \\ {}\vdots & \vdots & \vdots & \vdots \end{array}\right|\cdotp $$
(4)

The two-step approach has been proposed by Joyner and Boore (1981) to find space model parameters that can be used to model predictive relationships of the earthquake ground motions. The approach aims to reduce the effect of the distance–magnitude correlation on observed strong motion data sets. In addition, Joyner and Boore (1993) revealed that the methods of the stagewise regression yield comparable results with that obtained from a single-step analysis using the maximum likelihood technique. However, the methods based on the stagewise regression provide results with less accuracy than a single-step maximum likelihood analysis (Draper and Smith 1981); accordingly, the method based on the singular value decomposition of (Lanczos 1961) was used to solve the least-squares problem in this analysis. The applied method controls the stability to obtain the optimum estimates and uncertainties in the final solution (Press et al. 1992). Furthermore, the present analysis overcomes and quantifies the influence of the distance–magnitude correlation on the observed data set.

The site-effect term can be estimated using the residual between the observed and predicted PGA values:

$$ {res}_i=\ln \left(\frac{A_i}{Ac_i}\right), $$
(5)

where Aci indicates the predicted PGA or PGV values. Taking the average of the residual at each station, a coefficient representing the average site-effect term can be computed using the equation of Wu et al. (2001):

$$ S=\exp \left(\frac{1}{n}\sum \limits_{i=1}^N{res}_i\right)\cdotp $$
(6)

Results and discussion

It is worth stressing that seismograms recorded at good sites are an important reference when establishing the predictive GMPE in a region. Using the data set composed of 638 records from 72 earthquakes over the local magnitude range of 2.0–5.1 and the hypocentral distance range of 4–200 km, the empirical predictive relationships of the earthquake ground motions were defined in the Jazan Region. Simultaneous non-linear regression was applied to solve Eq. (4), and the empirical predictive relationships were defined for the PGA and PGV values based on the mean values estimated by taking the average values recorded for the two horizontal components. The results are summarized in the following equations.

$$ {\log}_{10}(PGA)=-1.36+0.85\cdotp {M}_L-0.85\cdotp {\log}_{10}(r)-0.005\cdotp r $$
(7)
$$ {\log}_{10}(PGV)=-1.05+0.65\cdotp {M}_L-0.66\cdotp {\log}_{10}(r)-0.04\cdotp r $$
(8)

The influence of the focal mechanisms on the predictive attenuation relationships of the ground motion has been mentioned and discussed in many studies (e.g., Campbell 1984, 1997; Sadigh et al. 1997; Boore et al. 1997; Anooshehpoor and Brune 2002; Kiuchi et al. 2019). These studies revealed that PGAs released from trust source mechanisms tend to be higher than other source mechanisms. Unfortunately, the number of fault plane solutions in this study, shown in Fig. 1, was insufficient to explore the effect of the source mechanisms. The station locations are listed in Table 1 with their corresponding site classes and site correction factors.

Table 1 List of seismic stations and their parameters (site correction factors, instrumental types, and site classes) used in this study
Full size table

The geological setting of the study area, as a part of the Arabian Shield, is primarily composed of exposed Precambrian basement, ranging in age from 1000 Myr to 545 Myr, having formed from thick sequences of late Oligocene to Holocene volcanic rocks developed with the spreading of the Red Sea (Camp and Roobol 1989). Figure 5 shows the effect of the coefficient d in Eq. (3), which represents the effect of the anelastic attenuation in western Saudi Arabia, indicating low anelastic attenuation in this region, which is in good agreement with the stable Precambrian geology of the Arabian Shield. The attenuation characteristics were estimated using the following form:

$$ \frac{1}{Q}=-\frac{\beta d}{\uppi f}, $$
(9)

where β is the shear wave velocity in the study region set to 3.5 km/s as taken from Rodgers et al. (1999) and f is the frequency in Hz. It is obvious from Fig. 7 that the results are comparable, to some degree, to the low attenuation characteristics that were determined by Pasyanos et al. (2009) and Abdel-Fattah et al. (2014) in the Arabian Shield. The frequency-dependent attenuation frequency generally indicates an inversely proportional relationship.

Fig. 5
figure 5

Frequency-dependent attenuation frequency in the studied region, where the relationship represents the effect of the calculated anelastic attenuation using the value of the coefficient d represented in Eq. (3). Dashed lines represent the 16th and 84th percentile limits of the standard deviation curves

Full size image

On January 23, 2014, a moderate-sized earthquake of ML 4.9 occurred in the southernmost part of the Arabian Shield in the Jizan Region. This earthquake had a focal depth of 10 km and exhibited a dextral strike-slip focal mechanism over the ENE fault trend, implying the reactivation of a vertical fault in the Precambrian basement. A maximum modified Mercalli intensity of V–VI was widely felt in the epicentral area, as reported by El-Hadidy (2015). Fortunately, this event was recorded by a total of 10 strong ground motion stations, providing unique data that were previously unavailable. Figure 6 shows a comparison of the derived PGA attenuation relationships (this study and that of Kiuchi et al. (2019)) with the observed values scaled to ML 4.9; these are plotted with the 16th and 84th percentile limits of the standard deviation curves for the attenuation relationship derived in this study. The derived relationships account for the site correction coefficient, which was estimated based on the residuals between the observed and predicted peak ground motion values. However, further studies will be required to calibrate the source terms when a large number of focal mechanism solutions are available. It is noteworthy that the relationship Kiuchi et al. (2019) indicated lower predictive ground acceleration than the currently derived equation. There are some factors that may probably cause the difference between two predictive relations due to the multiple different steps followed in the processing and regression analysis, the different data sets used in the two studies, and the influence of source mechanisms that were not accounted for in this study. However, our predictive relationship shows good agreement with the observed PGA values of the 2014 earthquake, as opposed to that derived by Kiuchi et al. (2019), which means that our relationship is capable of predicting the ground motions of future earthquakes in this region.

Fig. 6
figure 6

PGA attenuation curves derived in this study (solid blue line) and in that of Kiuchi et al. (2019) for strike-slip focal mechanisms (solid red line), with the strong ground motion (open blue squares) of an ML 4.9 earthquake that occurred in the Jizan Region in 2014 clearly showing the effect of the local site conditions at two stations located at hypocentral distances of approximately 200 km (encircled in red). Dashed lines represent the 16th and 84th percentile limits of the standard deviation curves for the attenuation relationship derived in this study

Full size image

The relationship that was developed in this study can be used not only to predict the ground accelerations from future earthquakes but also to generate earthquake shaking maps within a few minutes of the event based on routine estimates for magnitudes and hypocentral distances by incorporating the effects of local geological conditions. Figure 7 shows a comparison between the derived GPME in this study and seven different GMPE models derived in the Eastern United States and the Mississippi Valley (Thenhaus et al. 1989), Japan (Fukushima and Tanaka 1990), California (Sadigh et al. 1997), Greece (Skarlatoudis et al. 2004), Europe and the Middle East (Ambraseys et al. 2005), Europe and the Middle East (Akkar and Bommer 2010), and western Saudi Arabia (Kiuchi et al. 2019). The seven models are consistent with our model for distances ranging from 3 to 300 km, revealing that our model is capable of predicting ground motions from future large earthquakes.

Fig. 7
figure 7

A comparison between several GMPEs from different areas and the GMPE derived in the present study, revealing that the present GMPE is capable of predicting ground motions from future large earthquakes. The GMPE models are abbreviated as THEN for the Eastern United States and the Mississippi Valley (Thenhaus et al. 1989), FT for Japan (Fukushima and Tanaka 1990), SAD for California (Sadigh et al. 1997), SKA for Greece (Skarlatoudis et al. 2004), AM for Europe and the Middle East (Ambraseys et al. 2005), AB for Europe and the Middle East (Akkar and Bommer 2010), and KIK for western Saudi Arabia (Kiuchi et al. 2019)

Full size image

Conclusions

Predictive attenuation relationships for earthquake ground motions are of great importance when assessing seismic hazards in a specific region. These relationships are used by engineers to overcome design problems; an overestimation of ground motions may lead to additional costs in structural designs, while underestimations result in disasters in the event of large ground motions. The predictive attenuation relationships account for the uncertainty in PGA and PGV as a function of the earthquake size, source-to-site paths, and local geological conditions.

Based on our analysis of the means of 638 PGA and PGV values calculated from the horizontal component records, attenuation relationships were derived to predict the earthquake ground motions in the southernmost part of the Arabian Shield for earthquakes with magnitudes ranging from 2.0 to 5.1 and over hypocentral distances ranging from 4 to 200 km. The predictive attenuation relationships obtained in this study revealed good agreement with the published GMPE for western Saudi Arabia. Moreover, the predictive attenuation relationship of the PGA values compared well with the observed ground motions recorded from the 2014 earthquake. The attenuation relationships obtained in this study can therefore be used to assess the probabilistic seismic hazards and earthquake shaking maps in the region within a few minutes of the earthquake occurrences based on the prior information on magnitudes and hypocentral distances by incorporating local site responses.