Research on the extraction of pre-seismic anomalies in borehole strain data of the Maduo earthquake based on the SVMD-Informer model

Earthquake is a major natural disaster triggered by the accumulation and release of crustal stress, and the accurate extraction of pre-seismic anomaly signals is crucial to improve the earthquake prediction capability. In this study, an anomaly detection method for borehole strain data based on the combination of Segmented Variational Modal Decomposition (SVMD) and Informer network is proposed, and a pre-seismic anomaly extraction study is carried out for the 2021 Maduo M_s 7.4 earthquake in Qinghai. The SVMD method effectively solves the memory limitation problem of traditional Variational Modal Decomposition (VMD) when dealing with large-scale data through the sliding-window mechanism, and at the same time maintains the correlation between the data. The Informer network significantly reduces the computational complexity of the long-series prediction and realizes the high-precision one-time long-time series prediction by utilizing its ProbSparse self-attention mechanism and self-attention distillation. By analyzing the borehole strain data from the Menyuan station, this study identifies the accelerated anomaly accumulation phenomenon in the two stages before the Maduo earthquake: in the first stage, the number of anomalous days shows an accelerated growth starting from about 3 months before the earthquake (13 February 2021); in the second stage, the anomalous accumulation tendency is further intensified since the second month before the earthquake (the end of March 2021), and the accumulation curve shows a typical S-shape growth characteristic. The results are highly consistent with the time windows of the index of microwave radiation anomaly (IMRA), outward long-wave radiation (OLR) and geoelectric field anomalies, and with the subsurface-to-atmosphere multilayer anomalies (e.g., Benionff strain, CO concentration, electron concentration anomalies), which indicate that the borehole strain anomalies are closely related to the gestation process of the Maduo earthquake. This study provides a new method for the extraction of pre-seismic anomalies based on machine learning, and provides an important basis for understanding earthquake precursors.

1 Introduction

Earthquakes are caused by the release of accumulated stresses in the Earth's crust beyond the strength of the rocks during plate interactions and collisions. Once the stress in the earth's crust exceeds a certain critical threshold, the crust will rupture, releasing seismic waves reflected through the ground, which can cause great damage (Kanamori and Brodsky, 2001; Fan et al., 2021). They can damage infrastructure such as ground, transportation, and buildings, and may lead to secondary disasters such as volcanic eruptions (Seropian et al., 2021), tsunamis (Koshimura and Shuto, 2015), and landslides (Fan et al., 2019). Meanwhile, seismic hazards not only threaten human lives (Potter et al., 2015), but also have far-reaching impacts on socioeconomic development and quality of life (Peptan et al., 2023). Studies have shown that more than 80 % of earthquakes of magnitude 5.0 and above in mainland China in 2022 occur within the annual seismic hazard zone (Yu et al., 2024a). Therefore, it is crucial to study the anomalies that precede earthquakes.

Researchers around the world have now explored a wide range of phenomena before and after earthquakes, covering different structural levels of the Earth, including the subsurface, surface and atmospheric realms. Ma et al. (2022) identified geodetic anomalies prior to the Maduo earthquake through observed GPS values and computationally analyzed b-values, and verified that the pre-seismic anomalies of GPS data, b values, and stress and strain accumulations were associated with the Maduo earthquake. Liu et al. (2023) used outward longwave radiation (OLR) data to find that thermal infrared anomalies synchronized with the tidal stress cycle preceded the 2023 M 7.8 earthquake in Turkey, possibly reflecting the thermal response of tectonic stresses as they accumulate to a critical state. Guo et al. (2015) found significant anomalous ionospheric disturbances prior to the 11 April 2012 Sumatra M_s 8.6 and Mexico M_s 6.7 earthquakes using global total ionospheric electron content (TEC) data. In addition, other scholars have also studied fields such as geomagnetism (Li et al., 2019), microwave bright temperature (MBT) (Qi et al., 2021), CO (Cui et al., 2024), and electron density (Han et al., 2023). The above studies provide abundant data support and theoretical basis for the exploration of earthquake precursors, and are of great significance to our understanding of earthquake mechanisms and potential impacts.

Since the implementation of the U.S. Earth Lens Program, a large amount of open data has been accessed and utilized, and the Plate Boundary Observations (PBO) program has contributed to the development of strain techniques. Borehole strain observations have received unprecedented attention because of their high resolution and sensitivity (Lou and Tian, 2022; Roeloffs, 2010; Barbour and Agnew, 2012). Borehole strain observation is superior to GPS and laser strain gauges in capturing short- and medium-term strain changes and pre-seismic anomalies (Qiu and Shi, 2004), and is an important means to study crustal deformation and stress field changes. As the main observation equipment of China's digital seismic observation network, China's self-developed YRY-4 four-component borehole strain gauge is usually installed at the bottom of bedrock at 40 m, and has the capability of minute-level strain sampling, which can continuously record high-resolution stress and strain changes (Qiu et al., 2013). The high-resolution recordings provided by borehole strain gauges allow us to capture minute strain changes, thus providing accurate data to gain insight into crustal deformation processes (Lou and Tian, 2022). In addition, the borehole strain gauge not only provides four-component data, but also records ancillary observations such as solid tides, air temperature, and air pressure (Chi, 2009; Tang et al., 2023b). In addition to the application in earthquake precursor research, the borehole strain observation data also play an irreplaceable role in the research fields of slow earthquakes, co-seismic stress triggering, tremor, earth free oscillation, and seismic wave propagation, and many unique results have been achieved (Qiu, 2014).

Numerous scholars have accumulated rich research experience and results in extracting and identifying pre-seismic anomalous signals using borehole strain observation data. Chi (2013) inferred that these anomalies were closely related to the strain precursors in the preparation of two earthquakes by analyzing the tidal anomalies recorded by borehole strain gauges during the preparation of Wenchuan and Lushan strong earthquakes. Zhu et al. (2020) studied the anomalous characteristics of the borehole strain data before the Wenchuan earthquake by using principal component analysis. By analyzing the first eigenvalue and the first eigenvector of the borehole strain data, the characteristics of pre-earthquake crustal strain changes are revealed. Yu et al. (2020) constructed a complex network using multi-channel singular spectrum analysis (MSSA) using borehole strain data from the southwestern terminus of the Sichuan-Yunnan Longmenshan rift zone, and the results showed that the network provides a powerful tool for earthquake precursor monitoring. Wu (2012) revealed the intrinsic evolutionary characteristics of the seismic process through quantitative simulation of cluster sub-statistical equations and established a link between them and the dynamic change of earthquake precursors. Liu et al. (2014) used the S-transform method to analyze the observation data of four-component borehole strain gauges at Guzantai before and after the Lushan M_s 7.0 earthquake, and found that two clusters of high-energy clusters in the time-frequency domain may be related to the seismic activity in Lushan. In addition, Chi et al. (2019) extracted the Wenchuan earthquake borehole strain data anomalies by decomposing the borehole strain signal into multiple modes using VMD and using a new state-space model to determine the number of decomposed modes, and then calculating eigenvalues used to detect anomalies associated with the earthquake by using the anomaly detection method of principal component analysis (PCA). These studies show that borehole strain observations have significant advantages in precursor anomaly extraction and seismic correlation analysis.

Borehole strain observation can capture the subtle phenomena during seismic activities in a timely manner, and its observation data can reflect the stress-strain changes in rocks, thus providing a potential method for extracting strain anomalies before earthquakes. However, due to the high precision and wide bandwidth characteristics of borehole strain observation, it is highly susceptible to interference from external environmental factors. To address this problem, researchers have conducted in-depth studies on how to remove external interference from borehole strain data. Qiu et al. (2011) extracted pre-earthquake borehole strain anomalies by using high-pass filtering to remove the interference from seasonal variations and long-period signals. Subsequently, Zhang et al. (2019) performed first-order difference processing on the borehole strain data to significantly attenuate the low-frequency effects such as solar and lunar gravitational tidal forces in the measured data, and meanwhile, wavelet transform based on the first-order difference data enhanced the short-periodic component of the strain signal. In a recent study, Zhu et al. (2020) used harmonic analysis to eliminate the effects of solid tidal and seasonal trends on borehole strain data. Yu et al. (2021) used state-space modeling to remove the strain response due to seasonal trends, barometric pressure, solid tides, and water level variations, thus effectively isolating non-tectonic disturbances. Li et al. (2025) successfully extracted the pre-seismic anomalies of the M 7.0 earthquake in Jiuzhaigou by removing the effects of seasonal trends and tides on the borehole observations based on the Variational Modal Decomposition (VMD) and combining with the Graph WaveNet model to process the multi-station data. In addition, Zhu et al. (2024) firstly removed the long-term background trend and cyclic trend of the borehole strain observation data before the Lushan earthquake by time-series decomposition, secondly, decomposed the data by multi-channel singular spectrum analysis to eliminate the strain response due to water level and air pressure, and finally extracted the pre-earthquake negentropic anomaly in the borehole strain data.

Variational Modal Decomposition (VMD), as an adaptive signal decomposition method, is able to effectively extract the features of nonlinear and nonsmooth signals in the frequency domain, which is widely used in the analysis of complex waveforms such as steps, jumps, and burrs, and outperforms the traditional Empirical Modal Decomposition (EMD) and its derivatives in seismic signal processing (Rao et al., 2024; Li et al., 2018; Xue et al., 2019). However, with the increasing size and complexity of observed data, VMD has limitations such as high memory overhead in large-scale data processing. For this reason, this paper adopts the Segmented Variational Modal Decomposition (SVMD) method, which effectively solves the memory limitation problem through the sliding window mechanism while ensuring the correlation between the data. The studies of Chi et al. (2023) and Li et al. (2024) have shown that SVMD outperforms the traditional VMD method in terms of computational efficiency and memory utilization when dealing with large-scale data.

In recent years, machine learning techniques have become promising tools for studying earthquake precursor data. Studies have been conducted using Random Forest (RF) models (Tsuchiya et al. (2024), decision trees (Sikder and Munakata, 2009), Gated Recurrent Units (GRU) (Chi et al., 2023), Long Short-Term Memory (LSTM) (Zhang and Wang, 2023, 2024), neural networks (Kail et al., 2022; Bilal et al., 2022), and other methods that show significant potential in the prediction of seismic events, for example, estimating earthquake location and magnitude, while reducing the false alarm rate. Despite the remarkable progress, challenges remain in terms of data quality, scarcity, and heterogeneity. With the continuous development and optimization of machine learning models and the continued availability of high-quality data, the accuracy and reliability of earthquake prediction is expected to be further improved.

In this paper, we propose an innovative method based on SVMD-Informer network for extracting pre-seismic anomaly information from borehole strain data, and take the Maduo earthquake as an example to analyze the borehole strain data of Menyuan Terrace. The SVMD-Informer network used in this paper is not only good at capturing the features of the borehole strain data, but also demonstrates high prediction accuracy. The prediction intervals are constructed by statistical methods, thus realizing efficient anomaly detection. The research in this paper can provide more accurate pre-earthquake anomaly information for the earthquake early warning system and improve the accuracy and timeliness of earthquake prediction.

In order to analyze the correlation between the borehole strain data from Menyuan station and the Maduo earthquake, the process shown in Fig. 1 was used in this study. As shown in Fig. 1, first, the four-component borehole data were converted to strain data. Subsequently, the borehole strain data were decomposed using the Segmented Variational Modal Decomposition (SVMD) method, and the decomposition results were fused to remove the interference of influencing factors such as the annual trend and the solid tidal response. The fused data were used as inputs for training, validation and prediction using the Informer model. In the prediction stage, prediction intervals for the upper and lower bounds are established based on the model output, and anomalous signals are identified by comparing the prediction intervals with the original data. Finally, the pre-seismic anomalous cumulative values in the borehole strain data from Menyuan station were statistically analyzed, and the significant anomalous signals were successfully extracted.

Figure 1 Framework for borehole strain data processing and pre-seismic anomaly extraction.

2 Observational data and earthquakes

2.1 Borehole strain data

In recent years, several studies have verified the reliability of the sampled data from high-component borehole strain gauges, demonstrating the ability of four-component borehole strain gauges to record strained seismic waves from earthquakes of different magnitudes (Tang et al., 2023a; Qiu et al., 2015). As an ideal tool for crustal motion observation, the observation technology of the four-component borehole strain gauge is gradually maturing, which provides important support for the identification of earthquake precursors (Qiu et al., 2009). The YRY-4 four-component borehole strain gauge is one of the core instruments for the observation of crustal deformation in China, and it has the significant advantages of high sensitivity, wide observation band, good data consistency, and superior long-term stability. Borehole Strain Gauge is to place the sensor in the borehole for observation, relative to the earth its observation is an extremely small part of the crust deformation, can be approximated as a point of deformation observation results. The four-component borehole strain observation represents a relative observation that is capable of detecting changes in the target observation, but does not provide a complete measurement of the target observation. Its characterization is determined by the underlying principles of its model. The four probes on the borehole strain gauge are evenly spaced at 45° intervals to ensure multi-directional monitoring of crustal deformation. The measured value of any one element is recorded as S₁, rotated by 45° in turn, and the remaining three elements are recorded as S₂, S₃, and S₄. The amount of change in the observed values of the four elements satisfies the self-consistent Eq. (1):

$$\begin{array}{}\text{(1)}& S_{\mathrm{1}}+S_{\mathrm{3}}=k\left(S_{\mathrm{2}}+S_{\mathrm{4}}\right)\end{array}$$

This equation can be used to estimate the reliability of the data. k is the self-consistency coefficient, and in the ideal case k = 1. We consider the data to be reliable when k ≥ 0.95. Under plane strain conditions at or near the Earth's surface, only three independent variables are considered. Therefore, we can derive the various strains from the Menyuan station record using Eq. (2) as follows:

$$\begin{array}{}\text{(2)}& \left\{\begin{array}{l}{S}_{\mathrm{13}}=S_{\mathrm{1}}-S_{\mathrm{3}}\\ S_{\mathrm{24}}=S_{\mathrm{2}}-S_{\mathrm{4}}\\ S_{\mathrm{a}}=\left(S_{\mathrm{1}}+S_{\mathrm{2}}+S_{\mathrm{3}}+S_{\mathrm{4}}\right)/\mathrm{2}\end{array}\right.\end{array}$$

All three substitutions are important. Among them, S₁₃ and S₂₄ represent shear strains independent of each other, and S_a is the surface strain. Compared with the shear strains S₁₃ and S₂₄, the surface strain S_a is more representative of the four components measured by the YRY-4 borehole strain gauge, and therefore the data characterization of the surface strain S_a is taken as the object of study in this paper.

2.2 Maduo earthquake

According to the China Earthquake Networks Center (CENC), the M_s 7.4 magnitude earthquake occurred at 02:04:03 local time (LT) on 22 May 2021, in Maduo County, Qinghai Province, with an epicenter located at 98.34° E, 34.59° N, about 70–80 km from the East Kunlun Fracture Zone, at a depth of 17 km (Zhu et al., 2021). This great earthquake is located in the Bayankala block, which is one of the main blocks of the Qinghai Tibet Plateau orogeny and one of the regions with the most frequent seismic activity in China (Lü et al., 2022). The sudden-onset earthquake was the largest earthquake to occur in China since the 2008 Wenchuan M_s 7.9 earthquake. As of 08:00 LT 9 d after the earthquake, the China Earthquake Administration recorded a total of 2979 aftershocks, including one M_s 5.1, 13 M_s 4.0–4.9, and 63 M_s 3.0–3.9, which caused severe damage to buildings and roads in the area (Wang et al., 2021).

Dobrowolski's estimate of the radius of influence of precursors for earthquakes of different magnitudes is shown in Eq. (3) (Dobrovolsky et al., 1979):

$$\begin{array}{}\text{(3)}& \mathit{\rho }={\mathrm{10}}^{\mathrm{0.43}M}\mathrm{km},\end{array}$$

where M denotes the magnitude and ρ denotes the radius of influence of an M magnitude earthquake. The radius of influence of the Maduo earthquake is about 1377.2 km. In order to analyze the seismic effects near the epicenter, we selected the Menyuan station, which is closer to the epicenter of the Maduo earthquake, with a distance of 422.06 km from the epicenter. By analyzing the data collected at Menyuan station, we extracted 34 anomalous signals. This result indicates that the Menyuan station has the potential and ability to monitor earthquake-related anomalies. Figure 2 visually depicts the geographic location of the Menyuan station relative to the epicenter of the Maduo earthquake.

Figure 2 Location of the Menyuan observatory relative to the epicenter of the Maduo earthquake. The blue triangle represents the location of the borehole strain observatory. The purple star represents the location of the epicenter of the Maduo earthquake, while the red line inscribes the East Kunlun Fault Zone. This map was generated by GMT software, v. 6.0.0rc5 (https://gmt-china.org/, last access: 4 August 2025).

3 Method

3.1 Segmented variational modal decomposition (SVMD)

Since the introduction of VMD by Dragomiretskiy and Zosso (2014), it has achieved better results in dealing with data nonlinear, nonsmooth, and nonsmoothed signals. VMD is a fully implicit, adaptive, and completely nonrecursive approach to modal variational and signal processing. VMD is developed based on the classical Wiener filter denoising and Fourier transform, which aims to decompose the time series data into a series of intrinsic modal functions (IMFs) with finite bandwidth. The decomposition process essentially involves solving the variational problem, and the variational model can be expressed as follows:

$$\begin{array}{}\text{(4)}& \begin{array}{l}\stackrel{min}{\left\{u_k\right\},\left\{{\mathit{\omega }}_k\right\}}{\sum }_{k=\mathrm{1}}^K{∥{\partial }_t\left\{\left[\left(\mathit{\delta }\left(t\right)+\frac{j}{\mathit{\pi }t}\right)\times u_k\left(t\right)\right]e^{-j{\mathit{\omega }}_{k}t}\right\}∥}_{\mathrm{2}}^{\mathrm{2}}\\ \text{s.t.}{\sum }_{k=\mathrm{1}}^K{u}_k\left(t\right)=f\left(t\right)\end{array}\end{array}$$

In Eq. (4) {u_k} = $u_{\mathrm{1}}{u}_{\mathrm{2}},\mathrm{\dots },u_K$ and {ω_k} = ${\mathit{\omega }}_{\mathrm{1}}{\mathit{\omega }}_{\mathrm{2}},\mathrm{\dots },{\mathit{\omega }}_K$ are the corresponding center frequencies of the k modal functions and signal decompositions, respectively. Similarly, ${\sum }_k:={\sum }_{k=\mathrm{1}}^K$ is a shorthand notation for all modes and their center frequencies, respectively.

The quadratic penalty term α and the Lagrange multiplier λ(t) are introduced in the VMD in order to solve the constraint problem, which benefits from both the good convergence properties of the finite power quadratic penalty and the strict enforcement of the constraint by the Lagrange multiplier. Thus, the augmented Lagrangian is introduced:

$$\begin{array}{}\text{(5)}& \begin{array}{l}L\left(\left\{u_k\right\},\left\{{\mathit{\omega }}_k\right\},\mathit{\lambda }\right)\\ =\mathit{\alpha }{\sum }_{k=\mathrm{1}}^K{∥{\partial }_t\left[\left(\mathit{\delta }\left(t\right)+\frac{j}{\mathit{\pi }t}\right)\times u_k\left(t\right)\right]e^{-j{\mathit{\omega }}_{k}t}∥}_{\mathrm{2}}^{\mathrm{2}}\\ +{∥f\left(t\right)-{\sum }_{k=\mathrm{1}}^K{u}_k\left(t\right)∥}_{\mathrm{2}}^{\mathrm{2}}+〈\mathit{\lambda }\left(t\right),f\left(t\right)-{\sum }_{k=\mathrm{1}}^K{u}_k\left(t\right)〉\end{array}\end{array}$$

In the equation L denotes the Lagrange promotion operator, α denotes the data fidelity constraint function, and λ denotes the Lagrange multiplier. Meanwhile the alternating direction multiplier method (ADMM) is used in the VMD to find u_k, ω_k, λ and to solve the updated iterative optimization of Eq. (5).

Although VMD has the advantages of good performance and greater robustness to sampling and noise relative to existing mode decomposition models, when faced with a large number of data computations, problems such as slow data processing and computer memory limitations may arise because the VMD method requires a global search and solves the variational segmentation problem. Therefore, in this paper, we adopt a method that applies Variational Modal Decomposition (VMD) to data segments by combining sliding windows (SVMD). This method not only effectively solves the above problems in VMD. And through the sliding window mechanism, the correlation between data points can be effectively maintained (Chi et al., 2023).

The principle of SVMD method is shown in Fig. 3. We choose the sliding window mechanism with a size of 7 d and a sliding step of 1 d to realize the data segmentation. First we set the initial window as all the data from the first day to the seventh day and perform VMD decomposition of the data within that window. Starting from the second sliding window, only the results of the VMD decomposition of the current window are retained and superimposed with the decomposition results of the previous window, and in this logical order, the data are processed sequentially, and finally the complete dataset processed by SVMD is obtained.

Figure 3 Schematic diagram of the Segmented Variational Modal Decomposition (SVMD).

3.2 Informer network

In recent years, deep learning techniques have been widely used in the field of seismic hazard recognition and prediction, injecting new vigor into the development of earthquake science. Although traditional sequence models, such as long short-term memory networks (LSTM), gated recurrent units (GRU), and recurrent neural networks (RNN), have made significant contributions to time series prediction, they still have limitations in capturing long-range dependencies when handling long sequence time series prediction (LSTF) tasks. In contrast, the Transformer model, with its self-attentive mechanism, is able to reduce the maximum path length of network signals to a theoretical minimum while avoiding the complexity of the loop structure, thus showing great potential for LSTF problems (Vaswani et al. (2017). However, there are some serious issues with the traditional Transformer model that hinder its application to LSTF, including quadratic time complexity, high memory usage, inherent limitations of the encoder-decoder architecture, and the drawback of needing to rely on the previous prediction result for the current prediction. To address these issues, this paper employs the Informer model, which significantly improves the performance of long-series time series prediction while maintaining efficient computation.

Since the Informer model was proposed by Zhou et al. (2021), it has been widely used in the fields of river runoff time series prediction (Tepetidis et al., 2024), long-term prediction of indoor air quality (Long et al., 2023), wind power series prediction (Wei et al., 2023), financial time series prediction (Zhang et al., 2024), automatic driving trajectory prediction (Chen et al., 2023), etc. These applications fully show that Informer model has significant prediction performance advantages. The efficiency and accuracy of Informer model provide innovative ideas for earthquake research, especially in the analysis and prediction of earthquake precursor data, which shows great potential.

As shown in Fig. 4, the Informer model consists of three main components: encoder layer, decoder layer and prediction layer. Compared to the traditional self-attention methods, the encoder mainly deals with longer sequence inputs by using sparse self-attention. Self-attention distilling refers to distilling extracted from attention operations, which not only capture long-range dependencies in long input sequences but also reduce the model's dimension and network parameters. Multi-Head ProbSparse self-attention denotes the self-attention extracted block, which is used to stitch together the feature maps. Finally, the decoder processes the inputs from the long time series and zeroes out the target elements, to get the predicted results.

Figure 4 Framework diagram of Informer structure including encoder (Encoder) and decoder (Decoder).

The attention mechanism calculates the relevance between positions in a sequence by generating query (q), key (k), and value (v) vectors for the input sequence. Its core lies in obtaining attention scores through scaled dot products, thereby allocating attention weights. Based on this, Informer uses formula (7) instead of formula (6) in traditional Transformer for the calculation of attention. This enables the ProbSparse self-attention mechanism to significantly reduce the cost of self-attention computation while maintaining its effectiveness.

$$\begin{array}{}\text{(6)}& {\displaystyle }\text{Attention}\left(\mathbf{Q},\mathbf{K},\mathbf{V}\right)=\text{softmax}\left({\displaystyle \frac{{\mathbf{QK}}^T}{\sqrt{d}}}\right)\mathbf{V}\text{(7)}& {\displaystyle }\text{Attention}\left(\mathbf{Q},\mathbf{K},\mathbf{V}\right)=\text{softmax}\left({\displaystyle \frac{\stackrel{\mathrm{‾}}{Q}{\mathbf{K}}^T}{\sqrt{d}}}\right)\mathbf{V}\end{array}$$

$\stackrel{\mathrm{‾}}{\mathbf{Q}}$, K, V in Eq. (7) are the input matrices of the attention mechanism, where d denotes the dimensionality of the query and key vectors (dimensionality of the inputs) and $\stackrel{\mathrm{‾}}{Q}$ in Eq. (7) contains only Top-u queries under the sparsity metric $\stackrel{\mathrm{‾}}{M}({\mathbf{q}}_i,\mathbf{K})$. where for $\stackrel{\mathrm{‾}}{M}({\mathbf{q}}_i,\mathbf{K})$ an empirical method is proposed to effectively measure query sparsity, which is computed by a method similar to the Kullback–Leibler scattering. The calculation formula is as follows:

$$\begin{array}{}\text{(8)}& \stackrel{\mathrm{‾}}{M}\left({\mathbf{q}}_i,\mathbf{K}\right)=\underset{j}{max}\left\{{\displaystyle \frac{{\mathbf{q}}_i{\mathbf{k}}_j^T}{\sqrt{d}}}\right\}-{\displaystyle \frac{\mathrm{1}}{L_K}}\sum _{j=\mathrm{1}}^{L_K}{\displaystyle \frac{{\mathbf{q}}_i{\mathbf{k}}_j^T}{\sqrt{d}}}\end{array}$$

Here for every (query) q_i∈R^d and k_i∈R^d in the set of key K, we have the bound as $\ln L_K\le \stackrel{\mathrm{‾}}{M}\left({\mathbf{q}}_i\mathbf{K}\right)\le {max}_j\mathit{\left\{}{\mathbf{q}}_i{\mathbf{k}}_j^T/\sqrt{d}\mathit{\right\}}-\frac{\mathrm{1}}{L_K}{\sum }_{j=\mathrm{1}}^{L_K}\mathit{\left\{}{\mathbf{q}}_i{\mathbf{k}}_j^T/\sqrt{d}\mathit{\right\}}+\ln L_K$. The formula also holds when q_i∈K.

As shown in Fig. 5, the encoder contains three parts such as attention block(Attention Block), convolutional layer (Conv1d) and maximum pooling layer (MaxPool) for encoding the input data. After the main stack (Attention Block 1), Attention Block 2 and Attention Block 3 halve the input successively, thus increasing the reliability of the distillation operation, and the whole process continues by gradually reducing to $\mathrm{1}/\mathrm{4}$ of the original length. At the end of the encoder, all feature maps are concatenated and the output of the encoder is directed to the decoder. In addition, a self-attention mechanism is employed to evaluate the main features and generate a unified self-attention feature map for subsequent levels. Essentially, this process reduces the complex self-attention mechanism in the Transformer model to a simpler, smaller form that is suitable for integration into the Informer model.

Figure 5 Informer's encoder structure.

In Fig. 4, the decoder consists of two identical multi-head attention layers. The main difference is that predictions are generated through a process called generative inference, which greatly speeds up the long-term prediction process. The decoder provides the following vectors:

$$\begin{array}{}\text{(9)}& {\mathit{X}}_{\mathrm{de}}^t=\text{Concat}\left({\mathit{X}}_{\mathrm{token}}^t,{\mathit{X}}_{\mathrm{0}}^t\right)\in R^{\left(L_{\mathrm{token}}+L_y\right)\times d_{\mathrm{model}}}\end{array}$$

where ${\mathit{X}}_{\mathrm{token}}^t\in R^{L_{\mathrm{token}}\times d_{\mathrm{model}}}$ is the start marker and ${\mathit{X}}_{\mathrm{token}}^t\in R^{L_y\times d_{\mathrm{model}}}$ is a placeholder for the target sequence (set to 0). The computation of ProbSparse self-attention is then tuned by setting the inner product to −∞ via the masked multi-head attention mechanism, thus avoiding the involvement of each position in subsequent positions. Finally, a fully connected layer produces the final prediction.

For Informer model, mean square error (MSE) is used as a loss function for back propagation during training. A 0.2 dropout is applied during training to enhance the model generation. The Adam optimizer is used to update the weights with a learning rate of 0.000001 and the number of layers of the model is 128. The number of training rounds for this model is set to 6. When the model does not decline in 3 rounds of loss, the training of the model is stopped. This configuration allows the model to decay and effectively prevents overfitting by reducing the parameter magnitude, resulting in a better trained model.

Statistically based methods are widely used in constructing prediction intervals, and they are capable of detecting anomalies and analyzing anomalies in seismic precursors. Zhang et al. (2023b) constructed confidence intervals for the predicted values by calculating the residual values of raw infrared long-wave radiation (OLR) data. When the data exceeded the upper and lower limits of the confidence interval, they were determined as cold or hot anomalies, respectively. Chi et al. (2023) used the lower and upper bound estimation (LUBE) method to directly construct prediction intervals to extract anomalies in seismic precursors. Li et al. (2025) utilized a nonparametric method of constructing prediction intervals using data extremes to construct upper and lower bounds of prediction intervals directly. In this paper, the upper and lower limits of the predicted sequence output from the decoder are determined using a normal distribution method. The upper and lower bounds of the prediction intervals in the network are determined using the following Eq. (10):

$$\begin{array}{}\text{(10)}& \begin{array}{l}\text{Lower}=\text{Predicition}-Z\times \text{mad},\\ \text{Upper}=\text{Predicition}+Z\times \text{mad},\end{array}\end{array}$$

where Predicition is the predicted value; Z is the Z score of the normal distribution, which is approximately 1.44 at the 85 % confidence level; and mad denotes the median absolute deviation.

4 Data processing

The borehole strain data used in this study were provided by the Beijing Seismological Bureau. For short-term and proseismic anomaly extraction from borehole strain data. The borehole strain data of Menyuan station from 1 January 2020 to 31 May 2021 were selected for extracting the anomalous signals before the Maduo earthquake, taking into account the time of the earthquake that occurred in the Maduo earthquake. The borehole strain data can effectively extract anomalies in the short and medium term. In this paper, the relatively stable data from the first 8 months of 2020 were used as the training set and validation set, while the prediction set was selected from the data from the 7 months prior to the earthquake. The data from Menyuan station were validated by the self-consistent equation, and the four-component borehole strain data from Menyuan station were converted into two shear strain components S₁₃ and S₂₄ and one surface strain component S_a by using the strain conversion equation as shown in Fig. 6.

Subsequently, the S_a data from the Menyuan station were decomposed using the Segmented Variational Modal Decomposition (SVMD). Five modes are selected in the decomposition process, the decomposition parameters are set to 2000 bandwidths, and the convergence accuracy is set to 10⁻⁷. The decomposition results are shown in Fig. 6, where k1 represents the annual trend component, k2 represents the tidal component, and the remaining three components are summed up to obtain the final SVMD results. The decomposition results are compared with the relevant influencing factors, and the decomposition reveals that the method effectively eliminates the influence of seasonal trends and solid tides on the borehole strain data, and significantly improves the extraction of anomalous signals.

Figure 6 Data for S_a, S₁₃, S₂₄ at the Menyuan station and the results of the SVMD decomposition of S_a. k1 denotes the trend term, k2 denotes the earth tides, and k3, k4, and k5 denote the strains associated with the earthquakes. The red dashed line indicates the date of the earthquake.

5 Results and discussion

In this study, we use the SVMD-Informer network to extract the pre-seismic anomaly signals of the Maduo earthquake from the borehole strain data of the Menyuan station. The analysis focuses on identifying pre-seismic anomalies based on the obtained results. The anomalies were recognized when the raw data exceeded the corresponding upper or lower bounds. The prediction results of the SVMD-Informer network are shown in Fig. 7. From the figure, it can be seen that the real values are highly compatible with the prediction intervals, especially in the peak and valley regions of the real data, and the prediction interval can better capture the change trend of the data. This result indicates that we have a high reliability in predicting the borehole strain data using the SVMD-Informer network.

Figure 7 Detail of the prediction results of the Informer model for the borehole strain data at the Menyuan station. The red line indicates the actual data and the gray area indicates the prediction interval.

To identify anomalies in the prediction of borehole strain data, we used the following criteria: (a) detecting 15 points outside the interval within a 30 min window; (b) identifying interval bandwidths where the difference between the center of the predicted interval and the actual value is more than 1.5 times the interval bandwidth (more than three such points occurring in the same 30 min period). Days that meet these conditions are considered anomalous (Chi et al., 2023). This criterion can effectively distinguish random fluctuations from significant anomalies, thus significantly improving the accuracy of anomaly identification and providing a more precise basis for earthquake precursor analysis.

De Santis et al. (2017) first proposed an anomaly accumulation analysis method based on the S-shaped fitting function by studying the Swarm satellite data of the 25 April 2015 Nepal earthquake, and found that the S-shaped fitting function performs more superiorly in describing the anomaly accumulation process compared to the traditional linear fitting method. Similarly, Yu et al. (2024b) analyzed the geomagnetic data by using a geomagnetic model converter (ATGM) based on the self-attention mechanism and successfully extracted magnetic anomalies (MG-anomalies) in the Wenchuan, Lushan, and Kangding earthquakes, and found that the number of anomalies accumulating conformed to the S-type growth law. Li et al. (2024) used Graph WaveNet to extract pre-earthquake anomalies from borehole strain data of several stations near the earthquake source, from which they found that the anomalies accumulated at several stations before the Lushan earthquake showed an S-type accelerated upward trend, and further explained their correlation with the occurrence of the Lushan earthquake. Fan et al. (2024) used the Sigmoid function to fit the anomaly accumulation results of electron concentration, geomagnetic data, and their fusion parameters in the Maduo earthquake, and the results showed that the S-shaped curves can well reflect the process of pre-seismic acceleration and post-acceleration deceleration recovery of the earthquake. These results provide important theoretical support and reference basis for the anomaly accumulation analysis in this paper.

As shown in Fig. 8, the cumulative results of anomalous days at the Menyuan station from before to the time of the Maduo earthquake show a two-phase continuum of change. The first phase shows that the number of anomalous days from 13 February 2021 to the middle of the 2 months prior to the earthquake accelerated from the beginning of the anomaly on 13 February to the middle of March, and then levelled off. The second phase shows that the number of anomalous days has been in an accelerated increasing trend from the end of March until the earthquake. The timing of the first anomaly at Menyuan station coincides with the timing of the first-order anomaly of the index of microwave radiation anomaly (IMRA) found by Jing et al. (2022), which both appeared in mid-February. The timing of the second anomaly coincides well with the outward long-wave radiation (OLR) anomalies reported by Zhang et al. (2023a, b), and the observation of the thermal anomalies reveals that the Menyuan station is located within the range of the thermal anomalies of 21 and 22 March, which suggests that there was interplate thermal motion on the surface at this time. Moreover, the correlation coefficients of the geoelectric field began to show a significant decrease 2 months before the Menyuan earthquake (Xin and Zhang, 2021). Meanwhile Fan et al. (2024) calculated the Benionff strain S 90 d before the earthquake, and the results showed that the time of acceleration in the first stage was basically the same as the time when our second stage began to accelerate, inferring that there might be frequent lithospheric activities at this time. Therefore, the anomalies we found in the borehole strain data have some reliability and research value.

Figure 8 Cumulative results of anomalous days of borehole data at Menyuan station (MY). The red dashed line indicates the date of the earthquake, and the blue and red curves indicate the results of the S-fit function for the first and second phases, respectively.

In addition, the two-stage accelerated growth of the borehole strain accumulation results may reveal two preparatory mechanisms prior to the mainshock. This is consistent with the theory of fault synergistic process of Ma et al. (2014). They found that the occurrence of earthquakes is closely related to the three-stage synergistic evolution of faults through an indoor experimental study of plane-walk-slip faults. In the first stage, the initial stress nonlinear divergence leads to localized weakening and the formation of discrete strain release zones. The second stage is characterized by an increase in stress and a widening of the strain release zone. In the final stage, the expansion of the strain-release region and the rapid increase in the strain level in the strain-accumulation region. The anomalous cumulative acceleration about 3 months before the Maduo earthquake corresponds to the first and second stages in the theory, which is manifested by the deviation of the stress curve from linearity and the beginning of the formation and slow expansion of the discrete release zone. The acceleration 2 months before the earthquake reflects the characteristic changes of the third stage, which is characterized by the accelerated expansion of the release zone and the sharp increase of the strain in the accumulation zone. Therefore, we believe that the anomalies observed before the Maduo earthquake are related to the process of earthquake incubation.

Figure 9 Summary of our study of borehole strain anomalies (light blue) versus subsurface-to-atmosphere anomalies (from light green to gray) for the Maduo M_s 7.4 earthquake. The red dashed line indicates the date of the earthquake.

As shown in Fig. 9, related researchers have studied a variety of anomalies that occurred prior to the Maduo earthquake. From the observations at the surface, the magnetic field showed two accelerated phases, with the first phase showing an accelerated growth trend from 38 to 24 d before the main earthquake, and the second phase showing a significant increase in the accumulation of magnetic field anomalies from the beginning of the 24th day before the main earthquake to the time of the earthquake (Fan et al. (2024). In the transition from the surface to the atmosphere, the CO concentration at the 600 hPa altitude fluctuated from 17 February, increased again after a small decrease on 5 March, and showed the first anomaly peak on 21 March. The maximum anomaly of near-surface CO appeared in April, with a sudden increase on 14 April, and the anomaly increased to the maximum level before and after the earthquake in late April (Shi et al., 2024). It is noteworthy that the timing of the CO anomalies is highly consistent with the timing of our borehole strain anomalies. The temperature of brightness blackbody (TBB) started to increase significantly in the direction of the fracture zone north of the epicenter about 1.5 months before the earthquake, reached its maximum intensity on 17 May, and the anomaly gradually weakened during the earthquake (Yang et al., 2024). In the atmosphere, the electron density and electron temperature showed anomalous activity about 40 and 20 d before the earthquake (Tian et al., 2023; Fan et al., 2024; Du and Zhang, 2022). The synthesis in Fig. 9 contains anomalies covering a wide range of phenomena from our borehole strains in the subsurface to electron concentrations and temperatures in the atmosphere. It can be found that the borehole strain anomalies show a more comprehensive temporal coverage in the short-term anomalies and proseismic anomalies. In this regard, our study concludes that anomalies existed three to two months before the Maduo earthquake. Therefore, we hypothesize that the borehole strain anomalies may be related to the Maduo earthquake.

Figure 10 (a) Three-hourly regional variations of barometric pressure, temperature and precipitation in the Menyuan area during the period from 1 January 2020 to 31 May 2021. (b) Differential processing results of three-hourly regional variations of barometric pressure, air temperature and precipitation in Menyuan area during the period from 1 January 2020 to 31 May 2021. The red dashed line indicates the date of the earthquake, the two green dashed lines indicate the first phase anomaly, and the two royal blue dashed lines indicate the second phase anomaly.

Although borehole strain monitoring techniques are capable of accurately capturing crustal microstrain induced by plate tectonic movements or seismic activity, they are susceptible to factors such as temperature, air pressure, and rainfall. For this purpose we analyzed regional three-hourly variations of barometric pressure, air temperature and rainfall in the Menyuan area (35.97 to 39.97° N, 99.4 to 103.4° E) from NASA's Giovanni-4 platform (https://giovanni.gsfc.nasa.gov/giovanni, last access: 9 May 2025), with a time frame of 1 January 2020 to 31 December 2021. In Fig. 10a, the barometric pressure and air temperature fluctuate inversely within a certain range, while the rainfall gradually decreases year by year after peaking in summer, reflecting a significant annual cyclicity. In the two abnormal phases marked by dotted lines, meteorological factors showed obvious periodic fluctuations, but there was no obvious correlation with the borehole strain anomalies extracted in this paper. In addition, in order to minimize the influence of external factors on the borehole strain data, we performed differential processing on the three-hourly regional averages of barometric pressure, air temperature, and rainfall in the Menyuan area. Differential processing was utilized to remove the effect of cyclic variations, thereby highlighting anomalies in the data. The results are shown in Fig. 10b. In the two abnormal periods indicated by the dotted lines, the 3 h regional data for air pressure, temperature, and precipitation in the Menyuan area did not show any abnormal changes corresponding to this paper. Therefore, we have excluded the influence of pressure, temperature, and precipitation on the anomalies observed in the pre-earthquake borehole data from Maduo. Thus, we have reason to believe that the anomalies we extracted before the Maduo earthquake are related to the earthquake genesis process.

6 Conclusion

In this study, a SVMD-Informer network-based anomaly detection method for borehole strain data is proposed, and the 2021 Maduo M_s 7.4 earthquake is used as an example for pre-seismic anomaly extraction. The method addresses the issues of slow computation speed and memory limitations in traditional VMD by adopting SVMD, and significantly improves the accuracy and stability of long-sequence time series prediction by integrating the Informer network. By analyzing the borehole strain data from the Menyuan station, we successfully identified two distinct phases of anomalous cumulative acceleration preceding the Maduo earthquake, occurring approximately three and two months before the event, respectively. The cumulative anomaly curves exhibited a characteristic S-shaped growth pattern. This finding is consistent with the fault synergy process theory proposed by Ma et al. (2014), further supporting the correlation between borehole strain anomalies and the Maduo earthquake. With the continued progress of machine learning technology and the ongoing accumulation of seismic observation data, this method is expected to provide higher-precision technical support for earthquake prediction and help reduce seismic disaster risk.