Coal mine goaf surface deformation time series prediction and trend analysis method based on transformer model

Abstract

The application discloses a coal mine goaf surface deformation time sequence prediction and trend analysis method based on a Transformer model, relates to the technical field of mine safety monitoring and geological disaster early warning, and comprises the following steps: obtaining surface deformation time sequence, mining engineering dynamics and static geological data and preprocessing; constructing and training a dynamic stability adaptive prediction model, wherein the model comprises a point-by-point time sequence coding subnetwork, a graph attention knowledge coding subnetwork, a point-by-point cross-attention fusion module, a dynamic stability gate unit, a double-branch evolution module and an output module, and outputs a deformation prediction value and an instability risk index; and the model is trained by using a two-stage training strategy. Through the double-branch dynamic fusion architecture, the nonlinear disturbance branch weight is automatically enhanced when the system is unstable, the prediction stability is improved, and multi-source data fusion and instability early warning are realized.

Description

Technical Field

[0001] This invention relates to the field of mine safety monitoring and geological disaster early warning technology, specifically to a method for time-series prediction and trend analysis of surface deformation in coal mine goaf based on the Transformer model. Background Technology

[0002] After coal mining creates goaf areas, the overlying strata shift and break down, causing surface subsidence, collapse, and other deformations that seriously threaten buildings, roads, the ecological environment, and the safety of people in the mining area. Therefore, accurate monitoring and prediction of surface deformation in goaf areas is crucial.

[0003] In existing technologies, SBAS-InSAR technology can effectively monitor large-scale, long-term surface deformation, but its results are easily affected by factors such as spatiotemporal decoherence and atmospheric delay, and it cannot directly predict future deformation trends. Traditional prediction models, such as the probability integral method, are based on physical and mechanical mechanisms, but they have poor adaptability to complex geological conditions and multi-face mining, and parameter inversion is difficult. Conventional deep learning methods, such as LSTM or standard Transformer models, mainly have the following technical problems: poor long-term prediction stability, especially the error is amplified sharply in the unstable stage when deformation changes from uniform to accelerated; lack of utilization of geomechanical laws, and pure data-driven models are prone to learning spurious correlations; and difficulty in effectively integrating multi-source heterogeneous information such as InSAR deformation data, mining progress data, and geological structure data.

[0004] Therefore, there is an urgent need for a method for predicting and analyzing the temporal distribution and trends of surface deformation in coal mine goaf areas that can integrate multi-source heterogeneous data, embed geomechanical stability constraints, and maintain high-precision prediction during the deformation instability stage. Summary of the Invention

[0005] This invention aims to address the following problems in existing technologies: First, how to provide stable short- to medium-term predictions even when the surface deformation system tends towards instability; second, how to embed the concept of stability from geomechanics into a deep learning model so that its predicted trajectory conforms to physical laws; and third, how to construct an engineering-feasible multi-source heterogeneous data fusion framework. To this end, this invention proposes a method that integrates InSAR monitoring data and mining engineering data, utilizing an improved Transformer deep learning model to predict and analyze the trend of surface deformation in coal mine goaf areas.

[0006] To achieve the above objectives, the following technical solution is adopted:

[0007] This invention provides a method for time-series prediction and trend analysis of surface deformation in coal mine goaf based on the Transformer model, comprising the following steps:

[0008] The surface deformation time series data, mining engineering dynamic data and static geological data of each monitoring point in the coal mine goaf are acquired and preprocessed to obtain the surface deformation time series, dynamic mining characteristics and static geological characteristics of each monitoring point.

[0009] A dynamic stability adaptive prediction model is constructed and trained. The model includes: a point-by-point temporal encoding sub-network for extracting temporal feature vectors from the surface deformation time series; a graph attention knowledge encoding sub-network for fusing the static geological features and the dynamic mining features to generate a knowledge embedding vector; a point-by-point cross-attention fusion module for cross-attention fusion of the temporal feature vector and the knowledge embedding vector to obtain a fused feature vector; a dynamic stability gating unit for pooling the fused feature vector in the time dimension, inputting it into a fully connected network, and outputting stable weights and nonlinear weights; a two-branch evolution module, including a stable evolution branch and a nonlinear perturbation branch, for performing stable evolution prediction and perturbation evolution prediction respectively, and combining the stable weights and nonlinear weights to obtain the evolved fused feature vector; and an output module for mapping the evolved fused feature vector to deformation prediction values ​​for multiple future time steps and outputting an index characterizing the instability risk.

[0010] Furthermore, the surface deformation time series is obtained by processing with SBAS-InSAR technology; the dynamic mining characteristics are obtained by converting the daily working face advance distance into the same time resolution as the surface deformation time series using sliding window statistics; the dynamic mining characteristics include at least one of the cumulative advance distance within the window, the maximum daily advance speed within the window, or the variance of the advance distance within the window; the static geological characteristics include the distance from the monitoring point to the nearest fault, lithology coding, and the overlying strata strength index calculated based on lithology and thickness.

[0011] Furthermore, the point-by-point temporal coding sub-network adopts the PatchTST architecture: the land deformation time series of each monitoring point is subjected to first-order difference detrending processing; the detrended land deformation time series is divided into multiple time slices, and each time slice is input into a multi-layer Transformer encoder after linear embedding and location encoding, and outputs a temporal feature vector.

[0012] Furthermore, the graph attention knowledge encoding sub-network specifically includes: using all monitoring points as nodes, concatenating the static geological features and dynamic mining features of each monitoring point as the initial feature vector of that node; constructing a graph structure using a k-nearest neighbor graph strategy, where each node is only connected to a preset number of nodes with the closest spatial distance; calculating edge weights based on the distance between nodes and whether they are located in the same fault influence zone, wherein edge weights are increased between nodes located in the same fault influence zone; inputting the graph structure and the initial feature vectors of each node into a multi-layer graph attention network, and outputting the knowledge embedding vector of each monitoring point.

[0013] Furthermore, for each monitoring point, the point-to-point cross-attention fusion module uses its knowledge embedding vector as the query vector and its temporal feature vector as the key vector and value vector, performs cross-attention calculation to obtain a conditional feature vector, and then uses a gating mechanism to perform weighted fusion of the conditional feature vector and the temporal feature vector to obtain a fused feature vector.

[0014] Furthermore, the point-by-point cross-attention fusion module independently performs the following specific operations for each monitoring point: using the knowledge embedding vector of each monitoring point as the query vector; obtaining the key vector and value vector from the temporal feature vector of each monitoring point through linear transformation; calculating the similarity between the query vector and the key vector, and obtaining the cross-attention weight after normalization; weighted summing the cross-attention weight and the value vector to obtain the conditional feature vector; and weighted fusion of the conditional feature vector and the temporal feature vector of each time step through a gating mechanism to obtain the fused feature vector of each time step; wherein, the gating mechanism adopts a learnable gating scalar, which is obtained by concatenating the temporal feature vector of the current time step and the conditional feature vector, followed by linear transformation and a Sigmoid activation function.

[0015] Furthermore, the dynamic stability gating unit specifically includes: performing average pooling on the fused feature vector of each monitoring point in the time dimension to obtain a pooling vector; inputting the pooling vector into a multilayer fully connected network; and processing the output of the last fully connected network through the Softmax function to obtain stable weights and nonlinear weights, wherein the sum of the stable weights and the nonlinear weights is 1.

[0016] Furthermore, the dual-branch evolution module independently performs the following operations for each monitoring point: The fused feature vector at the current moment is input into the stable evolution branch, and a linear transformation is performed using a linear Koopman operator to obtain a linear evolution prediction value, which serves as the stable evolution prediction result; the fused feature vector at the current moment is subtracted from the linear evolution prediction value obtained by transforming the fused feature vector at the previous moment using a linear Koopman operator to obtain the linear prediction residual at the current moment; the linear prediction residuals from the current moment and a preset number of previous moments are taken to form a residual sequence, which is then input into the nonlinear perturbation branch and nonlinearly transformed using a Transformer decoder to obtain the perturbation evolution prediction result; based on the stable weights and nonlinear weights output by the dynamic stability gating unit, the stable evolution prediction result and the perturbation evolution prediction result are weighted and summed to obtain the fused feature vector at the next moment, i.e., the evolved fused feature vector.

[0017] Furthermore, the index characterizing the risk of instability is the nonlinear weight output by the dynamic stability gating unit; when the nonlinear weight exceeds the first warning threshold within a consecutive preset number of time steps, a yellow warning is issued, and when it exceeds the second warning threshold, a red warning is issued; and the Monte Carlo Dropout method is used to estimate the prediction uncertainty, and the mean and standard deviation of the prediction value are calculated through multiple forward propagations to output the confidence interval.

[0018] Furthermore, a two-stage training strategy is adopted to train the dynamic stability adaptive prediction model, including: a pre-training stage: identifying deformation stationary periods from the training data, and using only the data from the stationary periods to pre-train the pointwise temporal coding sub-network and the linear Koopman operator in the stable evolution branch, with the mean square error between the predicted deformation value and the actual deformation value as the loss function; a joint training stage: using all training data, training the complete dynamic stability adaptive prediction model with a total loss function, the total loss function including prediction loss, instability loss, and spectral radius constraint loss; wherein, the prediction loss is used to calculate the error between the predicted deformation value and the actual deformation value, the instability loss is used to supervise the alignment of the nonlinear weights output by the dynamic stability gating unit with the actual deformation acceleration label, and the spectral radius constraint loss is used to constrain the spectral radius of the linear Koopman operator to not exceed a preset upper limit.

[0019] Compared with the prior art, the present invention has the following beneficial effects:

[0020] First, this invention uses a dual-branch dynamic fusion architecture. When the system is stable, the model mainly uses linear Koopman prediction, and when the system tends to become unstable, it automatically increases the weight of the nonlinear disturbance branch, which alleviates the contradiction of using a global stable model to predict instability to a certain extent.

[0021] Second, this invention designs a method for converting high-frequency mining data to InSAR temporal resolution (sliding window statistics), which fuses static geological and dynamic mining features through a graph attention network, and achieves an organic combination of knowledge embedding and temporal features through point-by-point cross-attention and gating fusion, providing an engineering-feasible multi-source heterogeneous data fusion scheme.

[0022] Third, the nonlinear weights output by the gating network of this invention are related to the deformation acceleration, which can be used as a reference for the real-time instability risk index, providing auxiliary information for mine safety management and early warning decision-making.

[0023] It should be understood that the description in the Summary of the Invention is not intended to limit the key or essential features of the embodiments of the present invention, nor is it intended to restrict the scope of the invention. Other features of the invention will become readily apparent from the following description. Attached Figure Description

[0024] The above and other features, advantages, and aspects of the various embodiments of the present invention will become more apparent from the accompanying drawings and the following detailed description. The drawings are provided for a better understanding of the invention and are not intended to limit the scope of the invention. In the drawings, the same or similar reference numerals denote the same or similar elements, wherein:

[0025] Figure 1 This is a flowchart of the overall method for predicting and analyzing the surface deformation of coal mine goaf in an embodiment of the present invention.

[0026] Figure 2 This is a block diagram of the overall architecture of the dynamic stability adaptive prediction model according to an embodiment of the present invention;

[0027] Figure 3 This is a schematic diagram of the two-stage training strategy process according to an embodiment of the present invention. Detailed Implementation

[0028] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0029] Furthermore, the term "and / or" in this article is merely a description of the relationship between related objects, indicating that three relationships can exist. For example, A and / or B can represent: A existing alone, A and B existing simultaneously, or B existing alone. Additionally, the character " / " in this article generally indicates that the preceding and following related objects have an "or" relationship.

[0030] Figure 1 This is a flowchart illustrating the overall method for predicting and analyzing surface deformation trends in coal mine goaf areas according to an embodiment of the present invention. For details, see [link to relevant documentation]. Figure 1 A method for time-series prediction and trend analysis of surface deformation in coal mine goaf based on the Transformer model includes the following steps:

[0031] S1: Acquire and preprocess the time series data of surface deformation, dynamic data of mining engineering, and static geological data of each monitoring point in the coal mine goaf area to obtain the time series of surface deformation, dynamic mining characteristics, and static geological characteristics of each monitoring point.

[0032] Step S1 is used to acquire and preprocess data.

[0033] (1) Surface deformation time series data and preprocessing: Sentinel-1A up-orbit images covering the mining area were acquired (the revisit period is 12 days in single-satellite mode; if a dual-satellite network is used, the revisit period can be 6 days. This embodiment uses 12 days for single-satellite mode as an example). The images were processed using SBAS-InSAR technology to obtain the surface deformation time series (line-of-sight deformation) of high coherence points. Taking a two-year monitoring period as an example, approximately 100 to 120 time points can be obtained.

[0034] (2) Dynamic data and preprocessing of mining operations: Collect daily working face advance distance and thickness data. In order to align with the time resolution of the surface deformation time series data, a sliding window statistical method is used to convert the daily data into features with the same time resolution as the surface deformation time series data, including: cumulative advance distance within the window, maximum daily advance speed within the window, and variance of advance distance within the window.

[0035] (3) Static geological data and digital preprocessing: extract fault location and lithological zoning information, calculate the distance to the nearest fault, lithological code, and overlying stratum strength index calculated by empirical formula based on lithology and thickness for each monitoring point.

[0036] S2: Construct and train a dynamic stability adaptive prediction model, such as Figure 2 The diagram shown is an overall architecture block diagram of the dynamic stability adaptive prediction model according to an embodiment of the present invention. The dynamic stability adaptive prediction model includes the following sub-networks and modules:

[0037] Point-by-point temporal coding subnetwork 210 is used to extract temporal feature vectors from the land deformation time series;

[0038] (a) Constructing a point-by-point temporal coding subnetwork 210:

[0039] The point-by-point temporal coding subnetwork 210 adopts the PatchTST architecture: First-order difference detrending processing is performed on the land deformation time series of each monitoring point; the detrended land deformation time series is divided into multiple time slices, and each time slice, after linear embedding and location encoding, is input into a multi-layer Transformer encoder, outputting a temporal feature vector. Specifically:

[0040] For each monitoring point, its historical deformation sequence is subjected to first-order differencing to remove the linear trend, resulting in a detrended sequence. Where L is the historical time step. In this embodiment, L=96 (corresponding to approximately 1152 days, slightly longer than two years of data). When the actual sequence length is less than L, a sliding window method is used to construct training samples of length 96 from the original sequence.

[0041] The PatchTST architecture is adopted: the sequence of length L is divided into multiple time slices, with a time slice length of 16 time steps and a step size of 8 time steps, resulting in the number of time slices N = (96-16) / 8+1 = 11. Each time slice is mapped to a 128-dimensional space through linear embedding and a learnable positional encoding is added. The embedded sequence is then fed into a multi-layer Transformer encoder with 2 layers, 4 attention heads, and a feedforward network of 256 dimensions, outputting a feature vector. To align with subsequent knowledge features, linear interpolation is used to extend the features of N time slices back to L time steps, resulting in a temporal feature vector. .

[0042] Graph attention knowledge encoding subnetwork 220 is used to fuse the static geological features and the dynamic mining features to generate a knowledge embedding vector;

[0043] (II) Constructing a graph attention knowledge encoding subnetwork 220:

[0044] The graph attention knowledge encoding subnetwork 220 uses all monitoring points as nodes, concatenating the static geological features and dynamic mining features of each monitoring point as the initial feature vector for that node; it constructs a graph structure using a k-nearest neighbor graph strategy, where each node is connected only to a predetermined number of nodes with the closest spatial distance; it calculates edge weights based on the distance between nodes and whether they are located in the same fault influence zone, with additional edge weights added between nodes located in the same fault influence zone; the graph structure and the initial feature vectors of each node are input into a multi-layer graph attention network, outputting a knowledge embedding vector for each monitoring point. Specifically:

[0045] Let all monitoring points be nodes, and the number of nodes be M. Each monitoring point... The feature vector is composed of static geological features and dynamic mining features: static geological features include the distance to the nearest fault, lithological code, and overlying strata strength index, with a dimension of 3; dynamic mining features include the cumulative advance distance and average advance speed within a 12-day time window, with a dimension of 2. Therefore, the node feature dimension is 5.

[0046] Constructing a graph structure A k-nearest neighbor (kNN) strategy is used, where each node is connected only to its k nearest neighbors (k=10) to control graph density. The edge weight calculation formula is as follows:

[0047]

[0048] in, represents the distance between nodes i and j; nodes i and j represent different monitoring points. This is the distance decay factor (Gaussian kernel bandwidth), used to control the rate at which the edge weights between nodes decrease with distance; This is the fault influence coefficient, used to control the additional enhancement of edge weights between nodes located within the same fault influence zone; in this embodiment, =500 meters, =0.5; This is an indicator function; it takes the value 1 when both nodes are located within the same fault influence zone (defined as within 500 meters on either side of the fault), and 0 otherwise.

[0049] A graph attention network is used to extract knowledge embeddings: the network has two layers, with the first layer having an output dimension of 64 and the second layer having an output dimension of 128, and four attention heads. The output is a knowledge embedding vector. M represents the total number of highly coherent points (number of nodes) in the monitoring area.

[0050] The point-to-point cross-attention fusion module 230 is used to perform cross-attention fusion of the temporal feature vector and the knowledge embedding vector to obtain a fused feature vector.

[0051] (III) Constructing a point-by-point cross-attention fusion module 230:

[0052] For each monitoring point, the point-to-point cross-attention fusion module 230 uses its knowledge embedding vector as the query vector and its temporal feature vector as the key and value vectors. It performs cross-attention calculation to obtain a conditional feature vector, and then uses a gating mechanism to weightedly fuse the conditional feature vector and the temporal feature vector to obtain a fused feature vector. Further, the point-to-point cross-attention fusion module 230 independently performs the following specific operations for each monitoring point: using the knowledge embedding vector of each monitoring point as the query vector; obtaining the key and value vectors from the temporal feature vectors of each monitoring point through linear transformation; calculating the similarity between the query vector and the key vector, and obtaining the cross-attention weights after normalization; weightedly summing the cross-attention weights and the value vectors to obtain the conditional feature vector; and weightedly fusing the conditional feature vector with the temporal feature vector of each time step through a gating mechanism to obtain the fused feature vector of each time step. The gating mechanism uses a learnable gating scalar, which is obtained by concatenating the temporal feature vector and the conditional feature vector of the current time step, followed by a linear transformation and a sigmoid activation function. Specifically:

[0053] For each monitoring point The value of i ranges from 1 to M; perform the following operations independently:

[0054] Extract the vector of monitoring point i from the knowledge embedding as the query vector. The shape is (1, 128). Query vector Extracted from the knowledge embedding vector of the i-th monitoring point.

[0055] From the temporal characteristics of this point The key vector is obtained through linear transformation. Sum value vector All of them have a shape of (L, 128). It is a time series feature matrix with shape (L, 128), representing the time series features of the i-th monitoring point after interpolation; L: historical time steps (L=96 in this embodiment), used to construct a fixed-length residual sequence; : Key vector matrix, transformed from get; Value vector matrix, transformed from... get.

[0056] Calculate the cross-attention weights:

[0057]

[0058] in, : Cross-attention weight matrix, with shape (1,L), obtained by softmax normalization, representing the attention weight of the query vector at each time step.

[0059] The conditional feature vector is obtained by weighted summation of the attention weights and the value vector:

[0060]

[0061] in, The conditional feature vector, with a shape of (1,128), is obtained by weighted summation of attention weights and value vectors, and incorporates the overall attention information of knowledge embedding to temporal features.

[0062] This vector incorporates the overall attention information of knowledge embedding to temporal features.

[0063] To preserve timing details, Gated fusion with the original temporal features at each time step:

[0064]

[0065] in It is a learnable gated scalar, normalized by Sigmoid. The final result is the fused feature vector for each monitoring point. .in, : Time step index, with a value range from 1 to L; : No. The fused feature vector of each monitoring point at time step t has a shape of (1, 128); : A learnable gating scalar, normalized by the Sigmoid function and taking values ​​between (0,1), used to control the fusion ratio of conditional feature vectors and temporal features; ⊙: Element-wise multiplication (Hadamard product). : The learnable weight matrix of the gating mechanism maps the concatenated feature vectors to scalars; [·;·] : Vector concatenation operation, which... and Concatenate along the feature dimension; Sigmoid(·): activation function that maps the output of the gated network to the (0,1) interval; : No. The complete fusion feature vector matrix of each monitoring point has a shape of (L, 128);

[0066] The set of fused feature vectors from all monitoring points is denoted as .

[0067] in, : The fused feature vector set of all monitoring points, with shape (M, L, 128). M: The total number of highly coherent points (number of nodes) in the monitoring area.

[0068] The dynamic stability gating unit 240 is used to pool the fused feature vector in the time dimension, input it into the fully connected network, and output stable weights and nonlinear weights.

[0069] (iv) Constructing a dynamic stability gating unit 240:

[0070] The dynamic stability gating unit 240 specifically includes: performing average pooling on the fused feature vector of each monitoring point along the time dimension to obtain a pooled vector; inputting the pooled vector into a multi-layer fully connected network; and processing the output of the last fully connected network through a Softmax function to obtain stable weights and nonlinear weights. Details are as follows:

[0071] For each monitoring point, its fused feature vector is averaged and pooled over the time dimension to obtain a pooled vector. .Will The input is a multi-layer fully connected network with a structure of 128-dimensional → 64-dimensional → 2-dimensional, using ReLU as the activation function, and the last layer outputs stable weights via the Softmax function. and nonlinear weights ,and Each monitoring point receives an independent weight pair.

[0072] The dual-branch evolution module 250 includes a stable evolution branch and a nonlinear perturbation branch, which respectively perform stable evolution prediction and perturbation evolution prediction, and combine the stable weights and nonlinear weights to obtain the evolved fusion feature vector;

[0073] (v) Constructing a dual-branch evolution module 250:

[0074] The dual-branch evolution module 250 independently performs the following operations for each monitoring point: It inputs the fused feature vector at the current moment into the stable evolution branch, performs a linear transformation using the linear Koopman operator to obtain a linear evolution prediction value, which serves as the stable evolution prediction result; it subtracts the linear evolution prediction value obtained by transforming the fused feature vector at the previous moment using the linear Koopman operator from the fused feature vector at the current moment to obtain the linear prediction residual at the current moment; it takes the linear prediction residuals from the current moment and a preset number of previous moments to form a residual sequence, inputs the residual sequence into the nonlinear perturbation branch, performs a nonlinear transformation using the Transformer decoder to obtain the perturbation evolution prediction result; based on the stable weights and nonlinear weights output by the dynamic stability gating unit, it performs a weighted sum of the stable evolution prediction result and the perturbation evolution prediction result to obtain the fused feature vector at the next moment, i.e., the evolved fused feature vector. Specifically:

[0075] For each monitoring point Perform the following evolutionary operations independently (subscripts omitted below) ):

[0076] Let the fused feature vector at the current time t be... That is, fusing feature vectors The shape is (128,) and the dimension is 128;

[0077] (1) First branch (stable evolution): using the linear Koopman operator The operator adopts an orthogonal-diagonal-orthogonal parameterization form:

[0078]

[0079] Where K: linear Koopman operator, a 128×128 matrix used to describe linear evolution in the latent space; U: orthogonal matrix, a 128×128 matrix, the left orthogonality factor of the Koopman operator; V: orthogonal matrix, a 128×128 matrix, the right orthogonality factor of the Koopman operator. A diagonal matrix whose diagonal elements are composed of... The calculated shape is 128×128; : Diagonal matrix The The diagonal elements have values ​​ranging from (0, ..., ...) ); =0.95, This is the activation function, used to map the input to the (0,1) interval; These are learnable parameters. : Upper limit of spectral radius, used to constrain the maximum singular value of the Koopman operator; the spectral radius constraint is achieved by limiting the maximum element of the diagonal matrix to be less than 1, that is, the spectral radius is the maximum singular value of the matrix, and here the constraint is less than 1.

[0080] The linear evolution prediction is:

[0081]

[0082] Since the input sequence has already undergone detrended processing, the Koopman operator is applied to the detrended residual components, and its convergence properties do not contradict the long-term monotonic accumulation of deformation.

[0083] (2) Second branch (nonlinear perturbation): Calculates the linear prediction residual sequence for historical time points. Residual is defined as:

[0084]

[0085] Take the residuals from the current time step and the previous L-1 time steps to form a residual sequence of length L. Input this residual sequence into a lightweight Transformer decoder with 2 layers, 4 attention heads, and a 256-dimensional feedforward network, outputting a perturbation evolution prediction. .

[0086] Fusion Evolution: Fusion Feature Vector at the Next Time Step

[0087]

[0088] Stable evolution prediction, i.e., the eigenvector after linear transformation by the Koopman operator, has a shape of (128,). : The fused feature vector of the previous time step (t-1), with shape (128,); Current moment The linear prediction residual, i.e. the difference between the fused feature vector and the prediction value of the Koopman operator for the previous time step, has the shape (128,). : Perturbation evolution prediction, i.e. the feature vector output by the Transformer decoder, with a shape of (128,); Stable weights, output by the gating unit, represent the contribution ratio of linear evolution; Nonlinear weights, output by the gating unit, represent the contribution ratio of the nonlinear disturbance evolution, satisfying the following conditions: + =1; The fused feature vector at the next time step (t+1), i.e. the evolved fused feature vector, is obtained by weighted summation of stable evolution prediction and perturbation evolution prediction, and has a shape of (128,).

[0089] The output module 260 is used to map the evolved fused feature vector into deformation prediction values ​​for multiple future time steps and output an index characterizing the risk of instability.

[0090] (vi) Constructing output module 260:

[0091] Output module 260 uses a linear output layer to process the evolved fused feature vector. Deformation predictions are mapped to the next H time steps using a linear layer. In this embodiment, the prediction step size H=6, corresponding to 72 days.

[0092] (vii) Model Training

[0093] 1. Data preparation: In this embodiment of the invention, historical surface deformation time series data, mining engineering dynamic data and static geological data of each monitoring point in the coal mine goaf are collected and preprocessed. Then, the data are divided into training set, validation set and test set in chronological order, with proportions of 70%, 15% and 15%, respectively, for model training, hyperparameter tuning and performance evaluation.

[0094] 2. This model employs a two-stage training strategy:

[0095] like Figure 3 The diagram shown illustrates the two-stage training strategy of this invention. The two-stage training strategy includes: a pre-training stage: identifying deformation stationary periods from the training data, and pre-training the pointwise temporal coding sub-network and the linear Koopman operator in the stable evolution branch using only the data from these stationary periods. Specifically, the absolute value of the second-order difference of the deformation is identified from the training data as being less than 0.01 (unit: ...). The time step is 12 days, corresponding to approximately 0.0007. A period of more than three consecutive time steps (approximately 36 days) is considered a stationary period. Using only stationary period data, a pointwise temporal coding subnetwork and a linear Koopman operator are trained. The loss function is the mean squared error between the predicted and actual deformation values. The training run consists of 100 epochs with a learning rate of [missing information]. .

[0096] Joint training phase: Using all training data, a complete dynamic stability adaptive prediction model is trained with a total loss function, which includes prediction loss, instability loss, and spectral radius constraint loss, as detailed below:

[0097]

[0098] in: To predict the loss, Huber loss (δ=1.0) is used to calculate the error between the predicted deformation value and the actual deformation value;

[0099] For instability loss, a binary cross-entropy loss is used to supervise the alignment of the nonlinear weights of the gating unit with the actual deformation acceleration labels. It should be noted that the deformation acceleration phase is continuous in time; the acceleration label at the current moment can approximately reflect the evolution pattern at the next moment. Therefore, using the current acceleration label to supervise the nonlinear weights of the gating unit is reasonable. When the acceleration is greater than 0.02 (unit: , corresponding to When the label is positive, it should be used.

[0100] Here, is the spectral radius constraint term, where The spectral radius (i.e., the maximum modulus of the matrix eigenvalues) of the linear Koopman operator K is used to constrain the spectral radius of the linear Koopman operator to not exceed 0.95, thus ensuring the stability of the evolving system. , , , ..., These are the eigenvalues ​​of matrix K.

[0101] This invention approximates the calculation using the power iteration method or eigenvalue decomposition. And a spectral radius constraint term is introduced into the loss function to ensure < 0.95, thus ensuring the stability of the evolutionary system.

[0102] Using the AdamW optimizer, the learning rate... The batch size is 32, and the training rounds are 200. An early stopping strategy on the validation set is adopted (the training stops if the validation set loss does not decrease for 20 consecutive rounds).

[0103] Finally, the historical data of the monitoring points to be predicted are input into the trained dynamic stability adaptive prediction model to obtain the deformation prediction results. Simultaneously, the nonlinear weights output by the dynamic stability gating unit are... As a real-time instability risk index, it enables deformation trend analysis and early warning.

[0104] The thresholds for the early warning rules need to be calibrated based on historical deformation data and safety levels of the actual mining area. A yellow warning is issued when the nonlinear weight exceeds the first warning threshold within a preset number of time steps, and a red warning is issued when it exceeds the second warning threshold. Furthermore, the Monte Carlo Dropout method is used to estimate the prediction uncertainty. As an example, the following rule is used in simulation verification: when... A yellow alert is issued when the value exceeds 0.6 for three consecutive time steps; a red alert is issued when the value exceeds 0.8 for three consecutive time steps. In practical applications, the value should be recalibrated based on the geological conditions and safety requirements of the mining area.

[0105] The Monte Carlo Dropout method is used to estimate prediction uncertainty. The mean and standard deviation of the predicted values ​​are calculated through multiple forward propagations, and the confidence interval is output. The Dropout layer is applied to the output of the feedforward network and the output of the attention layer of the Transformer encoder. The Dropout rate is set to 0.1, and the forward propagation is performed 30 times. The mean and standard deviation of the predicted values ​​are calculated, and the output is... Confidence interval.

[0106] Example data description:

[0107] This embodiment uses simulation data for proof-of-concept, aiming to demonstrate the complete process and feasibility of the technical solution of the present invention. The simulation parameters are based on the typical deformation law of coal mines described in published literature: 100 monitoring points are generated, each simulating a three-stage deformation sequence of uniform settlement → accelerated settlement → stabilization, lasting 2 years with 110 time points (corresponding to the 12-day revisit cycle of a single Sentinel-1A satellite). Mining engineering data simulates daily advance records, and the cumulative advance distance and average advance speed are statistically analyzed over a 12-day window. Static geological data is randomly generated to include fault distance, lithological codes, and overlying strata strength indices.

[0108] It should be noted that SBAS-InSAR measured data contains complex noise such as atmospheric delay, decoherence, and phase unwrapping errors, and simulation data cannot fully simulate these interferences. Therefore, the quantitative results in this embodiment are only used to illustrate the relative performance of the model under ideal conditions and do not represent its absolute performance in actual mining areas. In actual deployment, the model should be retrained and its parameters calibrated using real InSAR monitoring data.

[0109] The model parameters are set as shown in Table 1 below:

[0110] Table 1

[0111]

[0112] Proof of concept results:

[0113] The model of this invention is compared with the standard Transformer model based on simulation data:

[0114] During the uniform settling stage, the root mean square error of the model of this invention and the standard Transformer are both in the range of 0.8~1.2 mm, and their performance is similar.

[0115] During the accelerated settlement phase, the maximum prediction error of the model in this invention is reduced by about 25% compared to the standard Transformer (reference results under simulation conditions), indicating that it has a certain adaptability to deformation mode transitions.

[0116] The nonlinear weights output by the gating network show a significant upward trend about 2 to 3 time steps (24 to 36 days) before the deformation acceleration begins to increase, which can be used as an auxiliary reference for early warning.

[0117] It should be noted that this model requires the InSAR input sequence to have high spatiotemporal continuity. In severely decoherent regions, phase unwrapping and interpolation preprocessing are necessary. When the number of monitoring points exceeds 2000, it is recommended to use a subgraph sampling strategy for the graph attention network to control the computational load. The warning threshold should be calibrated based on historical deformation data and engineering safety requirements of the actual mining area. The values ​​of 0.6 and 0.8 in this embodiment are merely examples and do not constitute a limitation on actual applications. The first threshold in the training parameters (0.01) ), second threshold (3 time steps), loss weight , This can be adjusted based on the actual data distribution. Generally, the recommended range is: first threshold 0.005~0.05, second threshold 3~10. , .

[0118] In summary, the method for time-series prediction and trend analysis of surface deformation in coal mine goaf based on the Transformer model provided in this embodiment of the invention effectively solves the problems of poor long-term prediction stability, lack of physical law constraints, and difficulty in fusing multi-source data in the existing technology through a dual-branch dynamic fusion architecture, a multi-source heterogeneous data fusion framework, and a two-stage training strategy. It has high engineering application value.

[0119] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

[0120] It should also be noted that, in the embodiments of this application, relational terms such as "first" and "second" are used only to distinguish one entity or operation from another entity or operation, and do not necessarily require or imply any such actual relationship or order between these entities or operations. Furthermore, the terms "comprising," "including," or any other variations thereof are intended to cover non-exclusive inclusion, such that a process, method, article, or apparatus that comprises a list of elements includes not only those elements but also other elements not expressly listed, or elements inherent to such a process, method, article, or apparatus. Without further limitations, an element defined by the phrase "comprising one..." does not exclude the presence of other identical elements in the process, method, article, or apparatus that includes said element.

[0121] The above description of the disclosed embodiments enables those skilled in the art to make or use this application. Various modifications to these embodiments will be readily apparent to those skilled in the art, and the general principles defined in the embodiments of this application may be implemented in other embodiments without departing from the spirit or scope of this application. Therefore, this application is not to be limited to the embodiments shown in this application, but is to be accorded the widest scope consistent with the principles and novel features disclosed in the embodiments of this application.

Claims

1. A method for time-series prediction and trend analysis of surface deformation in coal mine goaf based on the Transformer model, characterized in that, Includes the following steps: The surface deformation time series data, mining engineering dynamic data and static geological data of each monitoring point in the coal mine goaf are acquired and preprocessed to obtain the surface deformation time series, dynamic mining characteristics and static geological characteristics of each monitoring point. Construct and train a dynamic stability adaptive prediction model, the model comprising: The point-to-point temporal coding subnetwork is used to perform first-order difference detrending processing on the land deformation time series of each monitoring point. The detrended land deformation time series is divided into multiple time slices. Each time slice is input into a multi-layer Transformer encoder after linear embedding and location coding, and outputs a temporal feature vector. The graph attention knowledge encoding subnetwork is used to construct a graph structure with all monitoring points as nodes. It concatenates the static geological features and dynamic mining features of each monitoring point to form the initial feature vector of that node. A k-nearest neighbor graph strategy is employed to construct the graph structure, where each node is connected only to a predetermined number of nodes with the closest spatial distance. Edge weights are calculated based on the distance between nodes and whether they are located within the same fault influence zone, with additional edge weights added between nodes located within the same fault influence zone. The graph structure and the initial feature vectors of each node are input into a multi-layer graph attention network, which outputs a knowledge embedding vector for each monitoring point. The point-by-point cross-attention fusion module is used to perform cross-attention fusion of the temporal feature vector and the knowledge embedding vector to obtain a fused feature vector. The dynamic stability gating unit is used to pool the fused feature vector in the time dimension, input it into the fully connected network, and output stable weights and non-linear weights. The dual-branch evolution module includes a stable evolution branch and a nonlinear perturbation branch. For each monitoring point, it independently performs the following operations: The fused feature vector at the current moment is input into the stable evolution branch and linearly transformed using a linear Koopman operator to obtain a linear evolution prediction value, which serves as the stable evolution prediction result. The current fused feature vector is subtracted from the linear evolution prediction value obtained by transforming the fused feature vector at the previous moment using a linear Koopman operator, resulting in the current linear prediction residual. The linear prediction residuals from the current moment and a predetermined number of previous moments are used to construct a residual sequence. This residual sequence is input into the nonlinear perturbation branch and nonlinearly transformed using a Transformer decoder to obtain the perturbation evolution prediction result. Based on the stable weights and nonlinear weights output by the dynamic stability gating unit, the stable evolution prediction result and the perturbation evolution prediction result are weighted and summed to obtain the fused feature vector at the next moment, i.e., the evolved fused feature vector. The output module is used to map the evolved fused feature vector into deformation prediction values ​​for multiple future time steps and output an index characterizing the risk of instability.

2. The method according to claim 1, characterized in that, The surface deformation time series was obtained by processing with SBAS-InSAR technology; the dynamic mining characteristics were obtained by converting the daily working face advance distance into the same time resolution as the surface deformation time series using sliding window statistics; the dynamic mining characteristics include at least one of the cumulative advance distance within the window, the maximum daily advance speed within the window, or the variance of the advance distance within the window; the static geological characteristics include the distance from the monitoring point to the nearest fault, lithology coding, and the overlying strata strength index calculated based on lithology and thickness.

3. The method according to claim 1, characterized in that, For each monitoring point, the point-to-point cross-attention fusion module uses its knowledge embedding vector as the query vector and its temporal feature vector as the key vector and value vector, performs cross-attention calculation to obtain a conditional feature vector, and then uses a gating mechanism to perform weighted fusion of the conditional feature vector and the temporal feature vector to obtain a fused feature vector.

4. The method according to claim 3, characterized in that, The point-by-point cross-attention fusion module independently performs the following specific operations for each monitoring point: Use the knowledge embedding vector of each monitoring point as the query vector; The key vector and value vector are obtained from the temporal feature vector of each monitoring point through linear transformation; The similarity between the query vector and the key vector is calculated, and the cross-attention weights are obtained after normalization. The conditional feature vector is obtained by weighted summation of the cross-attention weights and the value vector; The conditional feature vector and the temporal feature vector of each time step are weighted and fused through a gating mechanism to obtain the fused feature vector of each time step; The gating mechanism employs a learnable gating scalar, which is obtained by concatenating the temporal feature vector of the current time step with the conditional feature vector, followed by a linear transformation and a Sigmoid activation function.

5. The method according to claim 1, characterized in that, The dynamic stability gating unit specifically includes: The fused feature vector of each monitoring point is averaged and pooled over time to obtain a pooled vector. The pooling vector is input into a multilayer fully connected network; The output of the last fully connected network is processed by the Softmax function to obtain stable weights and nonlinear weights, and the sum of the stable weights and the nonlinear weights is 1.

6. The method according to claim 1, characterized in that, The index characterizing the risk of instability is the nonlinear weight output by the dynamic stability gating unit; when the nonlinear weight exceeds the first warning threshold within a preset number of time steps, a yellow warning is issued, and when it exceeds the second warning threshold, a red warning is issued; furthermore, the Monte Carlo Dropout method is used to estimate the prediction uncertainty, and the mean and standard deviation of the prediction value are calculated through multiple forward propagations to output the confidence interval.

7. The method according to claim 1, characterized in that, A two-stage training strategy is employed to train a dynamically stable adaptive prediction model, including: Pre-training phase: Identify deformation stationary periods from the training data, and pre-train the pointwise temporal coding sub-network and the linear Koopman operator in the stable evolution branch using only the data from the stationary periods. The loss function is the mean square error between the predicted deformation value and the actual deformation value. Joint training phase: Using all training data, a complete dynamic stability adaptive prediction model is trained with a total loss function, which includes prediction loss, instability loss, and spectral radius constraint loss. The prediction loss is used to calculate the error between the predicted deformation value and the actual deformation value. The instability loss is used to supervise the alignment of the nonlinear weights output by the dynamic stability gating unit with the actual deformation acceleration label. The spectral radius constraint loss is used to constrain the spectral radius of the linear Koopman operator to not exceed a preset upper limit.

Images