The question of whether large and small earthquakes exhibit different characteristics during early rupture is central to understanding the evolution of rupture over time. For point source EEWSs that rely on information from the first few seconds of source rupture to rapidly estimate the source parameters, the ability to quickly and accurately calculate magnitude is crucial. In the past decade, researchers have made efforts to investigate whether the initial rupture is deterministic, but the results remain inconclusive.

The concept of deterministic assumption refers to the variation in the final magnitude of earthquakes due to different nucleation processes, and the waveform of the initial few seconds can predict the final magnitude. Some scholars have held that the frequency or amplitude of the initial P wave is proportional to the final magnitude of the earthquake, which allows for an estimate of the magnitude before the completion of the rupture process [25,35]. Colombelli et al. [27] discovered that the peak displacement of small earthquakes increases rapidly in the initial phase, while the peak displacement of large earthquakes grows slowly, which led to the conclusion that the P wave peak displacement evolves differently with time for different earthquakes in the early stages of the rupture process. On the other hand, Rydelek and Horiuchi [29] questioned the idea of an deterministic assumption, contending that the earthquake nucleation process is universal and independent of the final magnitude, and that the rupture process is ultimately unpredictable. Rupture unpredictability suggests that it is impossible to predict the final magnitude of an earthquake based on early rupture behavior alone. Moreover, Meier et al. [30] analyzed moderate-to-large earthquake events and found no evidence of differences between small and large earthquakes in terms of the onset of rupture. Trugman et al. [31] compiled oversize seismic datasets by measuring the peak ground displacement (P_d) of a progressively longer time window and assuming a time-dependent saturation of the linear relationship between log10 (P_d) and magnitude (Mw). The results indicated a universal growth pattern in the evolution of the P_d with time after the initial rupture of the fault, which was inconsistent with the deterministic model of seismic rupture.

In recent years, the focus of earthquake research has shifted from deterministic rupture nucleation to a concept known as rupture weak determinism, which concerns the ability to infer the final earthquake size after the nucleation. Melgar and Hayes [32] proposed a weakly determined model of rupture evolution, in which they found that large earthquakes rupture by slip pulses with self-similarity. They analyzed the variation in the average seismic moment rate over time by constructing an average source–time function with 0.5 magnitude units within the magnitude range of Mw 7.0 to 8.5. The results indicated that the average seismic moment rate was significantly different in the first 10 s, whereby it was much smaller than it was during the duration of the earthquake. This suggests that a self-similar slip pulse is formed soon after the rupture starts. So, purely deterministic rupture mechanisms may be ruled out. Nonetheless, in some cases, weak or probabilistic forms of deterministic assumption can still be observed through a detailed analysis of seismic or geodetic data. Goldberg et al. [33] used seismic and geodetic data to study early rupture behavior and concluded that it was insufficient to infer the final earthquake magnitude (EM) from the first few seconds of the initial rupture. However, accurate estimates were possible in the tens of seconds before rupture completion, which indicated a weak certainty. Meier et al. [34] observed plenty of shallow crustal seismic data records and showed that inferring the final EM from the initial rupture characteristics was impossible. The rupture characteristics of small and large earthquakes can be distinguished only when the rupture development reaches a certain degree. Hutchison et al. [36] suggested that when the earthquake length reaches 20% of the total length, the final length of the fault and the EM can be predicted with accuracy, provided that the earthquake slip is accurately known and the fault structural maturity is considered.

Most of the regional warning methods based on seismic networks use point source algorithms (PSAs), which treat an earthquake as a point source, calculate the source parameters of the earthquake and predict future IMs based on the GMM. The implicit assumption of these EM estimates that are based on P wave signals forms the basis of deterministic assumption. When an earthquake starts, most of the rupture sliding of the fault is limited to a small area near the ES. Seismic energy is released in a brief period, and small-to-moderate earthquake events align with this assumption. However, predicting the final magnitude is difficult when the rupture extends hundreds of kilometers across the fault before it is complete, as the time available to reach the P wave window to predict the final magnitude is limited. The EM estimate is likely not yet stabilized after the arrival of the S wave at the station. Furthermore, high-pass filtering is carried out for strong motion records to minimize long-period drift during velocity and displacement integration, which diminishes the low-frequency content in seismic records [37]. Hence, one drawback of using PSAs for regional warnings is that for earthquakes larger than Mw 7.0, the EM prediction suffers from a saturation of the magnitude estimate [35]. Another disadvantage is the uncertainty of using the GMM to predict IMs. The GMM is grounded on a substantial amount of seismic event statistics and uses epicenter distance and magnitude to estimate IMs at varying distances. The calculated IMs values represent only the expected results of using the standard regression equation and are not indicative of the actual IMs values at the target location. Thus, even if the exact earthquake location and magnitude are known when the GMM has a substantial error or the EM is imprecise, the estimated IMs will be unreliable. So, PSAs face limitations in three main aspects: the uncertainty of the ES parameter estimation, the saturation of the EM estimation and the uncertainty of the IMs prediction. Improvements in the PSA should be implemented from these perspectives.

When it comes to regional warnings, to predict the IMs parameters, the seismic magnitude needs to be calculated first. The traditional methods for EM estimation are mainly carried out by using amplitude algorithms, period algorithms, multiparameter combination methods and Bayesian-based methods. The amplitude algorithm estimates the magnitude according to the amplitude parameter within a few seconds of the initial P wave. Wu and Zhao [38] established the correlation between the source distance (R), magnitude (M) and peak displacement (P_d) based on the P_d of the initial P wave in the first 3 s. The period algorithm was first proposed by Nakamura [39], who argued that the EM was proportional to the frequency of the seismic wave and that the EM could be derived from the characteristic period τ_p of the seismic wave. While the single parameter-based magnitude estimation had a high uncertainty, the multiparameter combination method combined multiple parameters for magnitude estimation, which effectively improved the accuracy of the magnitude prediction. Huang et al. [40] found that the joint estimation of the two parameters τ_c, P_d was more effective than the single τ_c method and could be used as a new EM estimation method. The virtual seismologist method proposed by Cua and Heaton [41] was an EM estimation method based on Bayesian conditional probability distribution theory. Real-time estimates and updates could be made based on the attenuation relationship between the observed IMs values, prior information and the distance from the earthquake epicenter. In addition, machine learning algorithms for magnitude estimation are also available. Mousavi et al. [42] used a deep learning approach to predict magnitude by directly extracting multiple P wave feature parameters from waveforms. Zhang et al. [43] proposed a fully convolutional neural network (CNN) model based on earthquake early rupture information, which could be used for real-time earthquake detection, earthquake localization and magnitude estimation, etc.

GMMs are mathematical models that predict how the intensity of IMs varies with earthquake size, distance and other factors, and they play an important role in earthquake engineering and seismic damage assessments. When researching the attenuation models of IMs, models generally consider the effects of three aspects: source characteristics, propagation medium and site conditions. The source characteristics cover the magnitude, fault type and plate location of the earthquake; the propagation medium focuses on the geometric dispersion and energy dissipation and absorption of seismic waves, which are usually expressed by the epicenter distance and fault distance; and the site conditions focus on the influence of the site type such as bedrock or soil fields on the IMs. In areas with abundant earthquake observational records, attenuation relationships are empirical formulas obtained from statistically analyzing strong motion observation records. A very well-known GMMs project is the Next Generation Attenuation (NGA) project led by the Pacific Earthquake Engineering Research Center (PEER), which shows the trends of the next generation of attenuation relationships in the global digital network with abundant strong motion data [44]. The first phase (NGA-West) started in 2003 and ended in 2008 and focused on the attenuation patterns in California, the U.S.A. and other seismically active regions around the world. Five working groups developed five NGA models based on different research objectives and the selection of different databases [45]. The second phase of the NGA project started in 2010 and is divided into two parts, NGA-West2 and NGA East. NGA-West2 is a continuation of NGA-West and aims to improve the models of NGA-West, such as in terms of their directional effects, basin effects and topographic effects. NGA East focuses on areas of stable seismicity in Central and Eastern North America [46]. In addition, NGA-sub is a derivative of the NGA project, and this is a model of attenuation based on subduction zone seismic records [47,48]. The NGA project promotes the development of GMM research, which has led to significant progress.

Among the EEWSs currently in operation, the ShakeAlert warning system in the U.S. is a typical regional warning system. The ElarmS algorithms in the ShakeAlert system calculate the seismic position, size and other source parameters by using the limited number of parameters extracted from the first few seconds of the initial P wave observed by seismic stations. ElarmS estimates the EM by using data from at least the first four stations of the seismic wave arrival network; it takes the two parameters τ_c and P_d of the initial P wave and then predicts the IMs values by using the GMM. This information is later integrated and updated in real time to generate IM prediction maps [26,49]. The second-generation ElarmS-2 algorithm has been recorded and modularized, with improvements made to the station network configuration. Its processing speed reduced the alert time by 6 s, which resulted in an overall enhancement in the warning performance [50]. Since ElarmS uses the short-term average/long-term average (STA/LTA) method to trigger earthquakes, this triggering method is very sensitive and prone to generating false seismic events. The third generation of ElarmS-3 uses a new teleseismic filter and trigger filter to reduce false alarms. The teleseismic filter distinguishes teleseismic signals with a filter bank by using the fact that the high-frequency components are more attenuated and the low-frequency components are less attenuated during the long-distance propagation of the teleseismic seismic waves. The peak ground velocity (PGV) values of the narrow bandpass filter traces with nine frequency bands ranging from low to high frequencies are calculated to distinguish whether they are teleseismic, which is determined by assessing whether the frequencies are in a set range. The trigger filter implementation involves the use of a series of algorithms to analyze the waveform characteristics of the seismic signal. First, an amplitude check is performed to determine whether τ_p, P_d, P_v and P_a are within the set range to exclude any abnormally small or large amplitude detriggers. A “range post-trigger” parameter (R) is introduced to ensure that the signal is not a single pulse or a rapidly shifted nonseismic signal. Finally, the horizontal-to-vertical amplitude ratio is checked to prevent S wave triggers from entering the system [51]. After calculating the source parameters using the ElarmS algorithm, another module in ShakeAlert, the earthquake information to ground motion (eqInfo2GM), can calculate the shakemap based on the source parameters and provide IMs information in a map or contour format. Therefore, users can select the appropriate alert method for their target location through the application, which allows them to focus on taking action based on local risks [52].

The on-site warning system predicts IMs based on the characteristic parameters of the initial P wave observed by the seismic station. Compared with the regional warning method, the on-site warning method has fewer stations and has a relatively lower accuracy when estimating the source parameters. However, on-site warning only needs to predict future IMs values at the current location without considering the source parameters, which bypasses the uncertainty of source parameter estimation and the uncertainty of predicting the IMs field by using the GMM. This can provide a longer warning time for the area around the epicenter. Most of the current on-site warning systems use the P_d and τ_c of the initial P wave to estimate the IMs parameters. For example, the on-site algorithm in the ShakeAlert system provides warnings to locations up to 30 km away, and it provides 6 s of warning time to locations 50 km away [26]. The Italian SAVE on-site warning method has a success rate of over 80% in intensity prediction in the target area. The warning time is 8~10 s at 50 km and 15~18 s at 100 km [53].

In addition to the P_d and τ_c mentioned above, other characteristic parameters of the P wave information, such as the squared velocity integral (IV2) and cumulative absolute velocity (CAV), can also be used for on-site warnings. Wurman et al. [54] used P_v and P_d to estimate the magnitude. Odaka et al. [55] proposed that the epicenter distance can be estimated by using the waveform envelope fitting parameter, and then they constructed an empirical magnitude–amplitude relation by using the P wave amplitude. Festa et al. [56] proposed a characteristic parameter IV2 related to the energy released by the earthquake and investigated the relationship between the initial radiated energy and the magnitude of the earthquake inferred from the IV2. Alcik et al. [57] investigated the relationship between the CAV and epicentral distance and magnitude and adopted an on-site warning method based on Peak Ground Acceleration (PGA) and CAV thresholds. Additionally, Wang et al. [58] proposed a method to estimate the earthquake magnitude in real time by using a displacement squared integral (ID2) for EEWSs.

Determining how to use the initial P wave warning parameters to quickly estimate the IM parameters is an important issue in on-site warning research. The GMM that considers P wave characteristic parameters and IM parameters is constructed in the on-site early warning system so that the earthquake warning information is issued in time when the set threshold value is reached.

The P wave amplitude, frequency, IV2 and other related parameters can be used with the IM parameters to construct a GMM. For each of the P wave warning parameters, the basic GMM can be expressed as where Y represents the PGV or PGA and X is the P wave parameters (e.g., P_a, P_d, P_v). A and B are the coefficients to be obtained by the fit. For the threshold warning process, the implementation steps include (Figure 1):

• Waveform processing: when an earthquake is detected, remove the mean value and linear trend of the waveform and pick up the P waveform. Calculate the signal-to-noise ratio to eliminate data that may be contaminated by the noise for data quality control;
• P wave parameter calculation: integrate the accelerometer records once and twice to obtain the P_v and P_d records; filter them with a Butterworth high-pass filter with a cutoff frequency of 0.075 Hz to remove the low-frequency drift after the second integration; and obtain the P_d, P_v, τ_c and other parameters in the 3 s time window after the arrival of the P wave;
• Threshold setting: there is a good correlation between the seismic intensity parameter IMM and peak velocity and the early P wave peak displacement and IM parameter PGV [59]. By converting the intensity to the PGV, the threshold value of P_d is calculated by determining the empirical correlation between the P_d and PGV. Similarly, the threshold value of τ_c is determined by the correlation between τ_c and magnitude. For example, the P_d threshold and τ_c threshold are set to 0.2 cm and 0.6 s, respectively, for an earthquake with M > 6 and IMM ≥ 7 [15].
• Issue alert: judge whether the IM parameters exceed the set threshold, calculate the intensity level, determine the warning level and release the warning information.

Caruso et al. [53] employed a dataset of Italian earthquakes with magnitudes ranging from Mw 3.8–6.0 to measure the P_d and τ_c after the arrival of the P wave. They established relationships between the P_d and PGV and τ_c and Mw and derived connections between the P_d, τ_c and epicentral distance (R). Using these parameters, the GM damage potential and range in the vicinity of the station could be quickly estimated. In contrast to seismic network monitoring, the threshold warning method utilizes single or multiple stations in proximity to the ES region for monitoring. An alert is promptly issued to the epicentral area upon detecting an earthquake that exceeds the predetermined threshold. This approach circumvents the uncertainties associated with predicting IMs by using the GMM. The practical applicability of the threshold warning method has been demonstrated worldwide, such as in China [60], Japan [61] and Italy [62]. All of these studies exhibited the accuracy and efficacy of this method. Wang et al. [63] analyzed the role of two early warning parameters, τ_c and P_d, in magnitude estimation and proposed a threshold evolutionary magnitude estimation method based on τ_c and P_d. They recommended utilizing P_d within a range of 10 km to predict events in the larger magnitude range and suggested the joint use of τ_c to reduce the underestimation error for larger earthquakes.

An additional approach to constructing GMMs is multiparameter joint prediction, which predicts IMs by making full use of P wave information. By extracting multiple P wave initial characteristics, we can establish more relationships between the P wave parameters and IM parameters or magnitude, which leads to more accurate IMs predictions [64]. Meier et al. [65] proposed a new amplitude estimation algorithm named Gutenberg that improves the accuracy of on-site warning systems by utilizing the broadband frequency information of seismic signals. By implementing a novel filter bank technique that measures the absolute peak amplitude of the P wave over time, joint Bayesian estimates of the EM, ES and station distance are produced. Additionally, the correlation between the initial P wave parameters and the IM parameters is explored. Peng et al. [66] used filters of different orders to establish further relationships between initial P wave parameters and IM parameters and thereby improved the accuracy of on-site warning systems. Wang and Zhao [67] investigated the initial P wave estimation of IMs using the Wenchuan database. They selected eight characteristic parameter attributes of the P wave in a 1 to 3 s time window and established relationships with four IM parameters, and then they subsequently analyzed the accuracy of their estimation.

Artificial intelligence techniques provide another new research idea for on-site warning systems; using one or two P wave warning parameters to construct a simple GMM for fitting is no longer necessary. Multiple P wave feature parameters can be extracted from the raw waveform data by using machine learning to construct an artificial-intelligence-based GMM to establish a more complex relationship between the feature parameters and the IMs. Hsu et al. [68] extracted some P wave feature parameters from the first few seconds of an earthquake at a single station and used support vector regression methods to build a regression model to predict the PGA based on these features. Song et al. [69] used the least squares support vector machine model to construct a continuous prediction model of the PGV by selecting seven characteristic parameter inputs of the P wave. Chiang et al. [24] proposed the use of intelligent strong motion prediction to predict IMs, whereby a CNN is used to determine the relationship between the features extracted from the initial P wave and the strong motion is used to predict whether the GPA of subsequent waves exceeds a predetermined threshold. To reduce the complexity of multiparameter calculations, Hsu et al. [70] selected the two parameters of the first 3 s of two P wave signals and used an artificial neural network algorithm with the introduction of different Site Characteristic Parameters for PGA prediction.

Moreover, artificial intelligence techniques are used to predict the IMs of the target location directly from the observed P wave information. Hsu and Huang [71] performed a multiscale analysis to estimate future IMs based on the observed P wave data of the first 3 s of an earthquake detected at a single station by using the CNN approach. Jozinović et al. [72] introduced a technique for predicting the intensity of an earthquake by using a deep CNN. Using the information contained in the first 10 s P wave of the stations of neighboring earthquakes, a CNN model was used to predict the intensity of earthquakes at more distant stations, which could provide an estimate of the IM within 15–20 s after the occurrence of an earthquake. Compared with the GMM method, the CNN outperformed the GMM in terms of the residuals of the data. The CNN required no ES parameters, the receiver was located at the epicenter and its uncertainty was lower than that of the GMM. However, the weak interpretability of machine learning models can lead to false positives or false negatives, which makes it difficult to assess incorrect predictions [73]. The dataset was selected from stations in a specific geographic area for training, and if other stations are added, the model has to be retrained. Although it performed well for the complete sample, poor data quality (e.g., data loss, station failure) can occur in reality, so data from real situations are needed for model training.