Theoretical model constraint and data-driven fusion-based micro-seismic quantitative prediction method

Abstract

Disclosed is a theoretical model constraint and data-driven fusion-based micro-seismic quantitative prediction method, comprising: 1, obtaining original micro-seismic data of a coal mine working face to form an initial sample data set; 2, performing data preprocessing, and dividing into a training set and a test set; 3, constructing a long short-term memory (LSTM) model; 4, sampling historical micro-seismic data for a micro-seismic event quantitative prediction sample; 5, constructing a multi-objective optimization model containing seven sub-objective functions: a micro-seismic time fractal dimension, a micro-seismic space fractal dimension, a micro-seismic energy fractal dimension, a micro-seismic prediction total energy, a micro-seismic reconstruction mining stress, a micro-seismic activity process and a micro-seismic energy-frequency power law b value; and 6, finding an optimal sampling sample to act as a prediction value. By means of the described steps, the present invention achieves theoretical model constraint and data-driven fusion-based micro-seismic event quantitative prediction. The prediction effect is good, and the prediction result can effectively guide on-site production activity process control and implement pressure relief measures and disaster early warning.

Description

A quantitative microseismic prediction method based on the fusion of theoretical model constraints and data-driven approaches Technical Field

[0001] This invention relates to a quantitative prediction method for microseismic events based on the fusion of theoretical model constraints and data-driven approaches. Background Technology

[0002] In preventing rockbursts and destructive mine tremors, disaster prediction is a prerequisite for prevention and control. How to accurately and quantitatively predict rockburst disasters, and thus achieve precise disaster prevention and control, is a pressing issue that needs to be addressed in rockburst-prone mines in my country. Existing research shows that the occurrence of rockbursts is inevitably accompanied by coal and rock mass fracturing and microseismic response, which is closely related to the size of the fracturing. In practice, this is manifested in the time, space, and energy of microseismic events obtained by microseismic monitoring systems. Traditional prediction methods use this spatiotemporal information for qualitative prediction, such as the comprehensive index method, the multi-factor coupling method, and frequency-energy evolution law prediction, to classify and categorize rockburst hazards.

[0003] As coal mining in my country delves deeper, the increasing risk of rockbursts places higher demands on rockburst prediction. Qualitative predictions are no longer sufficient to meet the needs of engineering sites, making quantitative microseismic prediction crucial for rockburst forecasting. For example, patent application number 202211011826.9 discloses a quantitative prediction method for the three elements of time, space, and intensity of microseismic events in mines prone to rockbursts. This method uses deep learning regression to quantitatively predict the occurrence time, X, Y, and Z coordinates of the epicenter, and energy magnitude of microseismic events. This solves the current problem of lacking on-site monitoring parameters and the difficulty in quantifying the mapping relationship between the time, space, and intensity of disaster events. However, in actual coal mine production, microseismic activity is constrained and influenced by various factors, such as the influence of early-stage stress and geological structure environment, as well as later-stage prevention and control measures and mining conditions. The quantitative microseismic prediction process described in this patent does not fully consider these factors, and the prediction effect needs further improvement. Summary of the Invention

[0004] To address the aforementioned problems, this invention provides a quantitative microseismic prediction method based on the fusion of theoretical model constraints and data-driven approaches. It comprehensively considers the three elements of microseismic temporal, spatial, and strong intensity, as well as their correlation with rockburst, and conducts quantitative prediction of microseismic events based on the fusion of theoretical model constraints and data-driven approaches. It employs a novel concept of sampling historical data as samples for quantitative microseismic prediction and constructs a multi-objective optimization model based on seven sub-objective functions: microseismic temporal fractal dimension, microseismic spatial fractal dimension, microseismic energy fractal dimension, total predicted energy of microseismic events, microseismic reconstruction stress, microseismic activity process, and the microseismic energy-frequency power law b-value. The model applies the three elements of temporal, spatial, and strong intensity to the sampled microseismic events, thereby achieving quantitative microseismic prediction.

[0005] To achieve the above-mentioned technical objectives and effects, the present invention is implemented through the following technical solution:

[0006] A quantitative microseismic prediction method based on the fusion of theoretical model constraints and data-driven approaches includes the following steps:

[0007] Step 1: Obtain raw microseismic data and daily coal seam mining footage during the coal mine working face production period, and statistically obtain relevant derived variables as sample features to form an initial sample dataset;

[0008] Step 2: Preprocess the initial sample dataset and divide the preprocessed dataset into training and test sets;

[0009] Step 3: Construct a Long Short-Term Memory (LSTM) network model. Use a time window sliding method based on microseismic autocorrelation to determine the lag time of a single sample and set it as the initial memory days. Set the number of LSTM layers, fully connected layers, and neurons in the model. Use Mean Squared Error (MSE) and Mean Absolute Error (MAE) as loss functions. Adjust the model structure using training set loss and test set loss to obtain the optimal model for predicting the total frequency N of microseismic events on the next day. total and total energy M total ;

[0010] Step 4: Sample historical microseismic data as samples for quantitative prediction of microseismic events, including the three elements of microseismic time, space and energy;

[0011] Step 5: Construct a multi-objective optimization model containing seven sub-objective functions: microseismic time fractal dimension, microseismic spatial fractal dimension, microseismic energy fractal dimension, total energy predicted by microseismic activity, microseismic reconstruction mining stress, microseismic activity procedure, and microseismic energy-frequency power law b-value, to achieve constraints on the sampled microseismic events in the three dimensions of time, space, and intensity.

[0012] Step 6: Find the optimal sample as the predicted value.

[0013] Preferably, in step 1, the raw microseismic data includes the occurrence time and energy of microseismic events, and the relevant derived variables include the daily total frequency, daily total energy, daily average energy, daily maximum energy, daily minimum energy, daily total energy magnitude, daily average energy magnitude, daily maximum energy magnitude, and daily minimum energy magnitude of microseismic events.

[0014] The daily total frequency, daily total energy, daily average energy, daily maximum energy, and daily minimum energy of microseismic events compared to the previous day;

[0015] The daily total frequency of microseismic events, daily total energy magnitude, daily average energy magnitude, daily maximum energy magnitude, and daily minimum energy magnitude are the averages of the previous three days.

[0016] Preferably, in step 2, the frequency and energy data are transformed to eliminate outliers according to the following formula: M E =lg(E S +10), M N =log2(N) Frequency +1)

[0017] In the formula, ME represents the energy magnitude, ES represents the corresponding energy-related data, MN represents the total frequency level, and NFrequency represents the daily total frequency data.

[0018] The processed data is then normalized for each column using the maximum-minimum value normalization method. The calculation method is as follows:

[0019] In the formula, Here, xmax is the normalized value, xmin is the maximum value of the sample data, and x is the original value of the sample.

[0020] 80% of the preprocessed dataset was divided into the training set and 20% into the test set.

[0021] Preferably, the time window sliding method based on microseismic autocorrelation includes:

[0022] The correlation between microseismic data and its hysteresis data is calculated using the autocorrelation function (ACF). The calculation method is as follows:

[0023] In the formula, C ΔI representing microseismic data x at time i i Microseismic data x at time i+ΔI i+ΔI The covariance, ΔI is the lag in the time dimension, N is the total number of data to be calculated for autocorrelation, and μ is the mean of all data.

[0024] Autocorrelation (ACF) of two sets of microseismic data with a time lag ΔI ΔI The calculation method is as follows:

[0025] In the formula, ACF ΔI Let C0 be the autocorrelation of two sets of microseismic data with a time difference of ΔI, and C0 be the autocovariance when ΔI = 0, i.e., the initial term of the time window x. i The autocovariance;

[0026] ACF ΔI The value is between -1 and 1, when ACF ΔI Below the upper confidence interval C I Previously, the microseismic data sequences were correlated, when ACF ΔI Reaching CI At the upper limit, ΔI represents the maximum correlation period of the microseismic data. max When the time window lag time is greater than ΔI max At that time, it was assumed that the microseismic data had no correlation, and ΔI was used. max As the initial lag time length for a single sample in the LSTM model;

[0027] C I The calculation method is as follows:

[0028] The methods for calculating the mean squared error (MSE) and the mean absolute error (MAE) are as follows:

[0029] In the formula, N is the number of samples, Y i For the true value, These are predicted values.

[0030] Preferably, in step 4, the reliability function derived based on the Weibull probability distribution is used to verify the microseismic prediction simulation experiment:

[0031] The probability density functions of the three elements of time, space, and energy

[0032] In the formula, x is the value of the three elements, x≥0; β is the proportional parameter, β>0; α is the shape parameter, α>0;

[0033] Cumulative distribution function

[0034] The reliability function R(Δt) represents the probability that the system will operate normally within a time interval Δt. It is the complement of the cumulative distribution function, specifically expressed as: R(Δt)=1-F(Δt;β,α)

[0035] Right now

[0036] The absolute values ​​of the pairwise differences between the three key variables of a microseismic event are used as independent variables to measure the failure rate of the entire system. The reliability function between microseismic events is as follows:

[0037] Δlg E=|lg E i -lg E j |(i≠j)

[0038] In the formula, ΔT, ΔS, and Δlg E represent the time difference, horizontal distance difference, and energy level difference between microseismic events, respectively; R(ΔT), R(ΔS), and R(Δlg E) represent the failure rates of microseismic events in terms of time difference, horizontal distance difference, and energy level difference, respectively.

[0039] Construct the reliability objective function: max f(i,j)=ω1·R(ΔT)+ω2·R(ΔS)+ω3·R(Δlg E)

[0040] In the formula, max f(i,j) is the reliability between the i-th microseismic event in the real microseismic data and the j-th microseismic data in the historical microseismic event dataset. The larger the value, the greater the reliability between the two microseismic events. ω1, ω2, and ω3 are the weights of the reliability in terms of time difference, distance difference in planar location, and energy level difference, respectively.

[0041] Preferably, the models for the microseismic time fractal dimension, microseismic spatial fractal dimension, microseismic energy fractal dimension, total predicted energy of microseismic events, microseismic reconstruction mining stress, microseismic activity procedures, and the microseismic energy-frequency power law b-value are as follows:

[0042] In the formula, These represent the capacity dimensions of the microseismic data from the past seven days in terms of temporal fractals, spatial fractals, and microseismic fractals. Add sampling to the microseismic data of the past seven days to predict the capacity dimension of microseismic events in terms of temporal fractal, spatial fractal, and microseismic fractal capacity dimensions; min f(D T ),min f(D S ),min f(D E The values ​​represent the error magnitudes in three dimensions between the recent seven-day microseismic data and the sampled microseismic data; min f(Energy) represents the error magnitude between the predicted total energy of the next day's microseismic events Mtotal and the total energy of the sampled Ntotal microseismic data. i Let i be the energy magnitude of the i-th microseismic event in the sampled microseismic dataset. The total energy of the sampled microseismic event dataset; min f(N) MS This refers to the constraint on the microseismic time element by the early warning value of the frequency of microseismic activity processes quantified based on the correlation analysis of production activities. This is the early warning value for the frequency of microseismic activity processes for the i-th microseismic event in the sampled microseismic dataset. This is the sum of the early warning values ​​for the frequency of microseismic activity processes in the sampled microseismic event set. is the mean of the early warning values ​​for the frequency of microseismic activity processes in the sampled microseismic event set; min f(AI) is the constraint of the impact risk assessment value quantified based on the mining-induced stress reconstruction theory on the microseismic spatial elements, AI i Let i be the risk assessment value of the i-th microseismic event in the sampled microseismic dataset. The sum of the risk assessment values ​​for the sampled microseismic event set. denoted as the mean of the risk assessment values ​​of the sampled microseismic dataset; b0 is the b-value of the microseismic energy-frequency power law fitting of the microseismic data of the past seven days; b1 is the b-value of the microseismic energy-frequency power law fitting of the microseismic data of the past seven days with the addition of the sampled predicted microseismic events; min f(b) is the quantified error between the b-values ​​of the microseismic data of the past seven days and the microseismic energy-frequency power law fitting of the microseismic data with the addition of the sampled microseismic data.

[0043] Preferably, the formula for calculating the microseismic energy-frequency power law b-value is: lg N(≥lg E)=a-blg E

[0044] In the formula, lg E represents the microseismic energy level; N (≥lg E) represents the number of microseismic events with an energy level greater than or equal to lg E; a and b are constants, and their values ​​are calculated using the least squares method.

[0045] In the formula, m is the total number of energy level divisions; lgEi is the i-th energy level; Ni is the actual number of microseisms in the i-th energy level; and a and b are the estimated parameters obtained by least squares fitting.

[0046] Preferably, in step 6, the optimal sample is obtained as the predicted value using the multi-objective particle swarm optimization algorithm (MOPSO).

[0047] Preferably, the multi-objective particle swarm optimization (MOPSO) algorithm includes the following steps:

[0048] 1) Generate an initial particle swarm, each particle containing a position and velocity vector, as well as the corresponding objective function value. At the same time, initialize an external archive to store non-dominated solutions.

[0049] 2) For each iteration, perform the following steps:

[0050] ① Update particle velocity based on inertia weight, individual optimal solution pBest, and global optimal solution gBest: v id (t+1)=ω·v id (t)+c1·rand1·(pBest id -x id (t))+c2·rand2·(gBest id -x id (t))

[0051] In the formula, v id (t) is the velocity of the particle at iteration t, x id (t) is the position of the particle, pBest id It is the individual optimal position of the particle, gBest id ω is the global optimal position, c1 and c2 are learning factors, and rand1 and rand2 are random factors.

[0052] ② Adjust the particle position according to the updated velocity: x id (t+1)=x id (t)+v id (t+1)

[0053] ③ Perform non-dominated sorting on the individuals in the population to find the current non-dominated solution set;

[0054] ④ Compare the current non-dominated solution set with the solutions in the external archive, and update the external archive;

[0055] ⑤ Select the leader particle gBest from the external archive based on crowding distance and probability;

[0056] 3) The algorithm terminates when the maximum number of iterations or other termination conditions are reached; otherwise, step 2) is repeated.

[0057] 4) Output the non-dominated solution set in the external archive as the optimal solution set for the multi-objective function.

[0058] The beneficial effects of this invention are:

[0059] (I) This invention establishes a prediction model for the total energy and frequency of microseismic events considering the regulation of production processes; it constructs a multi-objective optimized quantitative prediction model for microseismic events containing seven sub-objective functions: time fractal dimension, spatial fractal dimension, energy fractal dimension, predicted total energy, reconstructed mining stress, activity process, and energy-frequency power law b-value. Through the above steps, this invention achieves quantitative prediction of microseismic events based on the fusion of theoretical model constraints and data-driven approaches. The prediction results are excellent and can effectively guide the regulation of on-site production processes, the implementation of pressure relief measures, and disaster early warning.

[0060] (II) By calculating the time, space, and energy fractal dimensions of microseismic events, it was verified that microseismic events exhibit self-similarity in time, space, and energy. Microseismic events similar to real events can be found in historical data in the three dimensions of time, space, and energy. Therefore, this invention proposes a new concept of sampling historical data as samples for quantitative prediction of microseismic events. A theoretical constraint model based on production activity correlation analysis, mining stress reconstruction, and energy-frequency power law distribution was constructed, realizing the constraint of sampled microseismic events in the three dimensions of time, space, and intensity.

[0061] (III) This invention introduces empirical knowledge, numerical simulation and other theoretical driving forces to drive information such as early stress and geological structure environment and later prevention and mining conditions. It combines microseismic data to drive the reconstruction of information such as mining stress and seismic wave CT detection stress. It uses mathematical models such as Weibull and normal probability density distribution to reconstruct the mining stress index. The average value of the corresponding index value of the predicted sample is taken. The larger the value, the closer it is to the real microseismic spatial location. Engineering examples have proved that the prediction effect of this invention is good. Attached Figure Description

[0062] Figure 1 is a flowchart of a quantitative microseismic prediction method based on the fusion of theoretical model constraints and data-driven approaches according to the present invention.

[0063] Figure 2 is a typical ACF diagram;

[0064] Figure 3 is a schematic diagram of the optimal model structure applicable to the prediction of total energy of microseismic events the following day in an embodiment of the present invention;

[0065] Figure 4 is a schematic diagram of the optimal model structure applicable to the prediction of the total frequency of microseismic events on the next day in an embodiment of the present invention;

[0066] Figure 5 is a schematic diagram of prediction simulation using the derived reliability function in an embodiment of the present invention;

[0067] Figure 6 is a schematic diagram comparing the actual microseismic data and the predicted microseismic data in terms of time variables, energy variables, and planar position variables in an embodiment of the present invention.

[0068] Figure 7 is a schematic diagram of the prediction simulation effect when the weights of the three variables ω1:ω2:ω3 = 2 / 7:4 / 7:1 / 7 in an embodiment of the present invention;

[0069] Figure 8 is a schematic diagram showing the results of the entire optimal solution set of three constraint samplings in an embodiment of the present invention, with three constraints: energy, risk assessment value, and early warning value of frequency of microseismic activity processes.

[0070] Figure 9 is a schematic diagram comparing the location of the three-sample prediction dataset and the actual microseismic data in an embodiment of the present invention. Detailed Implementation

[0071] The technical solution of the present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments, so that those skilled in the art can better understand the present invention and implement it. However, the embodiments are not intended to limit the present invention.

[0072] As shown in Figure 1, a quantitative microseismic prediction method based on the fusion of theoretical model constraints and data-driven approaches includes the following steps:

[0073] Step 1: Obtain raw microseismic data and daily coal seam mining footage during the coal mine working face production period, and statistically obtain relevant derived variables as sample features to form an initial sample dataset, which is used in Step 3 to construct a Long Short-Term Memory (LSTM) deep learning model to predict the total frequency and total energy of microseismic events on the next day.

[0074] Preferably, in step 1, the raw microseismic data includes the occurrence time and energy of microseismic events. Step 1 not only uses the daily total frequency and energy of microseismic events as input parameters for the model, but also incorporates daily coal seam mining footage data and related derived variables, including further processing of the frequency and energy data and mean data considering geological factors. The related derived variables are obtained statistically from the raw microseismic data, including:

[0075] Daily total frequency, daily total energy, daily average energy, daily maximum energy, daily minimum energy, daily total energy magnitude, daily average energy magnitude, daily maximum energy magnitude, and daily minimum energy magnitude of microseismic events;

[0076] The daily total frequency, daily total energy, daily average energy, daily maximum energy, and daily minimum energy of microseismic events compared to the previous day;

[0077] The daily total frequency of microseismic events, daily total energy magnitude, daily average energy magnitude, daily maximum energy magnitude, and daily minimum energy magnitude are the averages of the previous three days.

[0078] Step 2: Preprocess the initial sample dataset and divide the preprocessed dataset into training and test sets.

[0079] When predicting total frequency and total energy, the minimum and maximum values ​​often differ significantly, and the statistical frequency and energy values ​​may contain zero values, affecting the accuracy of the prediction. Therefore, the frequency and energy data are transformed to eliminate outliers from the original data.

[0080] Preferably, in step 2, the frequency and energy data are transformed to eliminate outliers according to the following formula: M E =lg(E S +10), M N =log2(N) Frequency +1)

[0081] In the formula, ME represents the energy magnitude, ES represents the corresponding energy-related data (such as the daily total energy, daily average energy, daily maximum energy, and daily minimum energy of microseismic events), MN represents the total frequency level, and NFrequency represents the daily total frequency data.

[0082] The processed data is then normalized for each column using the maximum-minimum value normalization method. The calculation method is as follows:

[0083] In the formula, Here, xmax is the normalized value, xmin is the maximum value of the sample data, and x is the original value of the sample.

[0084] The preprocessed dataset can be divided into 80% as the training set and 20% as the test set.

[0085] Step 3: Construct a Long Short-Term Memory (LSTM) network model. The time window sliding method based on microseismic autocorrelation is used to determine the lag time of a single sample and set it as the initial memory days. The number of LSTM layers, fully connected layers, and neurons in the model are then set. The LSTM layers, acting as long short-term memory layers, control the information flow through input, forget, and output gates, allowing for flexible control of information transmission and forgetting when processing time-series data. The fully connected layers convert the output features of the LSTM layers into the final prediction results. Mean Squared Error (MSE) and Mean Absolute Error (MAE) are used as loss functions. The model structure is adjusted using training set loss and test set loss to obtain the optimal model for predicting the total frequency N of microseismic events for the next day. total and total energy M total .

[0086] Preferably, the time window sliding method based on microseismic autocorrelation includes:

[0087] The correlation between microseismic data and its hysteresis data is calculated using the autocorrelation function (ACF). The calculation method is as follows:

[0088] In the formula, C ΔI representing microseismic data x at time i i Microseismic data x at time i+ΔI i+ΔI The covariance, ΔI is the lag in the time dimension, N is the total number of data to be calculated for autocorrelation, and μ is the mean of all data.

[0089] Autocorrelation (ACF) of two sets of microseismic data with a time lag ΔI ΔI The calculation method is as follows:

[0090] In the formula, ACF ΔI Let C0 be the autocorrelation of two sets of microseismic data with a time difference of ΔI, and C0 be the autocovariance when ΔI = 0, i.e., the initial term of the time window x. i The autocovariance;

[0091] ACF ΔI The value of is between -1 and 1, and it decreases as ΔI increases. Figure 2 shows a typical ACF plot. When ACF ΔI Below the upper confidence interval CI Previously, the microseismic data sequences were correlated, when ACF ΔI Reaching C I At the upper limit, ΔI represents the maximum correlation period of the microseismic data. max When the time window lag time is greater than ΔI max At that time, it was assumed that the microseismic data had no correlation, and ΔI was used. max As the initial lag time length for a single sample in the LSTM model;

[0092] C I The calculation method is as follows:

[0093] The methods for calculating the mean squared error (MSE) and the mean absolute error (MAE) are as follows:

[0094] In the formula, N is the number of samples, Y i For the true value, These are predicted values.

[0095] Step 4: Sample historical microseismic data as a quantitative prediction sample for microseismic events, including three elements: time, space, and energy. The time element is the time of occurrence excluding the date, such as 12:34:52. The space element is the relative coordinates with respect to the current working face mining footage center position, not the absolute coordinates. The energy element is the energy level, i.e., the logarithm of the absolute value of the microseismic energy.

[0096] Using historical microseismic data as samples for quantitative prediction of microseismic events is a novel concept proposed in this invention. Preferably, in step 4, a reliability function derived based on the Weibull probability distribution is used to verify the microseismic prediction simulation experiment.

[0097] The probability density functions of the three elements of time, space, and energy

[0098] In the formula, x is the value of the three elements, x≥0; β is the proportional parameter, β>0; α is the shape parameter, α>0;

[0099] Cumulative distribution function

[0100] The reliability function R(Δt) represents the probability that the system will operate normally within a time interval Δt. It is the complement of the cumulative distribution function, specifically expressed as: R(Δt)=1-F(Δt;β,α)

[0101] Right now

[0102] The absolute values ​​of the pairwise differences between the three key variables of a microseismic event are used as independent variables to measure the failure rate of the entire system. The reliability function between microseismic events is as follows:

[0103] Δlg E=|lg E i -lg E j |(i≠j)

[0104] In the formula, ΔT, ΔS, and Δlg E represent the time difference, horizontal distance difference, and energy level difference between microseismic events, respectively; R(ΔT), R(ΔS), and R(Δlg E) represent the failure rates of microseismic events in terms of time difference, horizontal distance difference, and energy level difference, respectively.

[0105] Construct the reliability objective function: max f(i,j)=ω1·R(ΔT)+ω2·R(ΔS)+ω3·R(Δlg E)

[0106] In the formula, max f(i,j) is the reliability between the i-th microseismic event in the real microseismic data and the j-th microseismic data in the historical microseismic event dataset. The larger the value, the greater the reliability between the two microseismic events. ω1, ω2, and ω3 are the weights of the reliability in terms of time difference, distance difference in planar location, and energy level difference, respectively.

[0107] Step 5: Construct a multi-objective optimization model containing seven sub-objective functions: microseismic time fractal dimension, microseismic spatial fractal dimension, microseismic energy fractal dimension, total predicted energy of microseismic events, microseismic reconstruction mining stress, microseismic activity procedures, and microseismic energy-frequency power law b-value. This model aims to constrain the sampled microseismic events in the three dimensions of time, space, and intensity.

[0108] Preferably, the models for the microseismic time fractal dimension, microseismic spatial fractal dimension, microseismic energy fractal dimension, total predicted energy of microseismic events, microseismic reconstruction mining stress, microseismic activity procedures, and the microseismic energy-frequency power law b-value are as follows:

[0109] In the formula, These represent the capacity dimensions of the microseismic data from the past seven days in terms of temporal fractals, spatial fractals, and microseismic fractals. Add sampling to the microseismic data of the past seven days to predict the capacity dimension of microseismic events in the time fractal, spatial fractal, and microseismic fractal dimensions; min f(D T ),min f(D S ),min f(D E The values ​​represent the error magnitudes in three dimensions between the recent seven-day microseismic data and the sampled microseismic data; min f(Energy) represents the error magnitude between the predicted total energy of the next day's microseismic events Mtotal and the total energy of the sampled Ntotal microseismic data. iLet i be the energy magnitude of the i-th microseismic event in the sampled microseismic dataset. The total energy of the sampled microseismic event dataset; min f(N) MS (N) represents the constraint on the microseismic time element by the early warning value of the frequency of microseismic activity processes quantified based on the correlation analysis of production activities. MSi This is the early warning value for the frequency of microseismic activity processes for the i-th microseismic event in the sampled microseismic dataset. This is the sum of the early warning values ​​for the frequency of microseismic activity processes in the sampled microseismic event set. is the mean of the early warning values ​​for the frequency of microseismic activity processes in the sampled microseismic event set; min f(AI) is the constraint of the impact risk assessment value quantified based on the mining-induced stress reconstruction theory on the microseismic spatial elements, AI i Let i be the risk assessment value of the i-th microseismic event in the sampled microseismic dataset. The sum of the risk assessment values ​​for the sampled microseismic event set. denoted as the mean of the risk assessment values ​​of the sampled microseismic dataset; b0 is the b-value of the microseismic energy-frequency power law fitting of the microseismic data of the past seven days; b1 is the b-value of the microseismic energy-frequency power law fitting of the microseismic data of the past seven days with the addition of the sampled predicted microseismic events; min f(b) is the quantified error between the b-values ​​of the microseismic data of the past seven days and the microseismic energy-frequency power law fitting of the microseismic data with the addition of the sampled microseismic data.

[0110] The formulas for calculating the time fractal dimension, spatial fractal dimension, and energy fractal dimension of microseisms are as follows:

[0111] In the formula, D0 is the fractal capacity dimension, ε represents the box measurement scale in the fractal dimension calculation method - the box counting method, and N(ε) represents the minimum number of boxes of length ε required to cover all microseismic events on the time, space, and energy scales respectively. For detailed calculations, please refer to the invention patent (CN201910210964.1, A comprehensive intelligent method for predicting shock hazard levels based on microseismic fractals).

[0112] The calculation of mining-induced stress in microseismic reconstruction incorporates information on pre-existing stress and geological structural environment, as well as post-existing prevention and mining generalization, driven by empirical knowledge and numerical simulation. It also integrates information on damage reconstruction mining-induced stress and seismic wave velocity distribution obtained through CT detection, driven by microseismic monitoring data. For detailed calculations, please refer to the invention patent (CN202111337335.9, A Method for Assessing Rockburst Risk Based on the Fusion of Theory and Data).

[0113] The calculation of microseismic activity sequence employs a time window sliding method to divide the microseismic sequence into equal time window segments. Simultaneously, an accumulation method is used to calculate the microseismic frequency per hour within each time window, and the hourly warning value is obtained based on the standard deviation coefficient. For detailed calculations, please refer to the invention patent (CN202210872839.9, An Impact Warning Method Based on Microseismic Activity Sequence Analysis).

[0114] Preferably, the formula for calculating the microseismic energy-frequency power law b-value is: lg N(≥lg E)=a-blg E

[0115] In the formula, lg E represents the microseismic energy level; N (≥lg E) represents the number of microseismic events with an energy level greater than or equal to lg E; a and b are constants, and their values ​​are calculated using the least squares method.

[0116] In the formula, m is the total number of energy level categories; lgEi is the i-th energy level; Ni is the actual number of microseisms in the i-th energy level; and a and b are estimated parameters obtained by least squares fitting. The value of a reflects the magnitude of microseismic activity: the larger the value, the stronger the microseismic activity. The value of b reflects the magnitude of microseismic intensity: the larger the value, the greater the proportion of low-energy microseismic events.

[0117] Step 6: Find the optimal sample as the predicted value.

[0118] Preferably, in step 6, the optimal sample is obtained as the predicted value using the Multi-Objective Particle Swarm Optimization (MOPSO) algorithm. The MOPSO algorithm includes the following steps:

[0119] 1) Generate an initial particle swarm, each particle containing a position and velocity vector, as well as the corresponding objective function value. At the same time, initialize an external archive to store non-dominated solutions.

[0120] 2) For each iteration, perform the following steps:

[0121] ① Update particle velocity based on inertia weight, individual optimal solution pBest, and global optimal solution gBest: v id (t+1)=ωv id (t)+c1·rand1·(pBest id -x id (t))+c2·rand2·(gBest id -x id (t))

[0122] In the formula, v id (t) is the velocity of the particle at iteration t, x id (t) is the position of the particle, pBest id It is the individual optimal position of the particle, gBestid ω is the global optimal position, c1 and c2 are learning factors, and rand1 and rand2 are random factors.

[0123] ② Adjust the particle position according to the updated velocity: x id (t+1)=x id (t)+v id (t+1)

[0124] ③ Perform non-dominated sorting on the individuals in the population to find the current non-dominated solution set;

[0125] ④ Compare the current non-dominated solution set with the solutions in the external archive, and update the external archive;

[0126] ⑤ Select the leader particle gBest from the external archive based on crowding distance and probability;

[0127] 3) The algorithm terminates when the maximum number of iterations or other termination conditions are reached; otherwise, step 2) is repeated.

[0128] 4) Output the non-dominated solution set in the external archive as the optimal solution set for the multi-objective function.

[0129] The present invention will be further illustrated below with reference to the embodiments. In this example, the microseismic monitoring system data of the 7302 working face of Zhaolou Coal Mine from February 2, 2021 to October 31, 2022 is selected as the original data for LSTM model prediction, and a multi-objective optimization model is constructed to realize quantitative prediction of microseismic events.

[0130] The specific steps for implementing this invention are as follows:

[0131] (1) Based on the original microseismic data, the daily total frequency, daily total energy, daily coal seam mining footage, daily average energy, daily maximum energy, daily minimum energy, daily total energy magnitude, daily average energy magnitude, daily maximum energy magnitude, daily minimum energy magnitude, frequency and energy series data increments compared to the previous day, and the three-day average of frequency and energy magnitude series data were statistically obtained to form an initial sample set.

[0132] (2) Preprocess the initial sample dataset to remove outliers from the original data, eliminate zero values ​​in frequency and energy values, normalize the data after outlier processing according to the maximum and minimum value normalization method, and divide 80% of the preprocessed dataset into the training set and 20% into the test set.

[0133] (3) Determine the maximum correlation period ΔI of the microseismic data based on the time window sliding method of microseismic autocorrelation. max=16 days. The initial memory days for a single LSTM sample are set to 16 days. An LSTM model is constructed, and MSE and MAE are used as model evaluation metrics. The model is trained to predict the total daily frequency and the total daily energy, respectively. The structure is adjusted to obtain the optimal model structure for predicting the total energy of the next day as [Memory days: 16, LSTM: 2, Dense: 6, Units: 256], and the optimal model structure for predicting the total daily frequency as [Memory days: 20, LSTM: 6, Dense: 4, Units: 1024]. The model structure is shown in Figures 3 and 4, where the angle brackets indicate the vector dimension of the parameters in each layer.

[0134] (4) Samples are taken from historical microseismic data as samples for quantitative prediction of microseismic events, including three elements: time, space and energy of microseismic events. The time element is the time of occurrence without the date, such as 12:34:52. The space element is the relative coordinates of the current working face mining footage center position, rather than the absolute coordinates. The energy element is the energy level, that is, the logarithm of the absolute value of the microseismic energy.

[0135] The microseismic event that occurred at the mine on July 15, 2022 was selected for prediction simulation to verify the feasibility of the proposed concept.

[0136] 1000 microseismic data points prior to July 15 were selected, and the reliability function was statistically analyzed and calculated, as shown in Figure 5. The reliability function was then used for prediction simulation.

[0137] First, simulation predictions are performed using a reliability function for a single variable, including time, space, and energy. As shown in Figure 6, the predictions for time, energy, and planar location variables show that the actual microseismic data events are almost identical to the predicted data. That is, under single-variable constraints, future microseismic events can be completely simulated and predicted from events sampled from historical microseismic data.

[0138] Secondly, the prediction simulation results were compared using multivariate constraints: Figure 7 shows the effect of the prediction simulation when the weights of the three variables ω1:ω2:ω3 = 2 / 7:4 / 7:1 / 7. Spatially, the approximate location area of ​​the microseismic event can be predicted, and the time and energy level prediction errors between the actual microseismic event and the predicted microseismic event data are very small.

[0139] The prediction results show that the reliability function can effectively constrain variables for prediction. In univariate constraint prediction, microseismic events that closely match the actual microseismic events can be predicted. That is, microseismic events similar to the actual microseismic events can be found in historical microseismic data from the three dimensions of time, space, and energy. When using multivariate constraint prediction, the weights between the three variables can be adjusted according to the requirements of the variables to achieve a more accurate prediction effect. This verifies the concept proposed in this invention of sampling historical microseismic data as samples for quantitative microseismic prediction.

[0140] (5) Construct a multi-objective optimization model containing seven sub-objective functions: microseismic time fractal dimension, microseismic spatial fractal dimension, microseismic energy fractal dimension, total energy of microseismic prediction, microseismic reconstruction mining stress, microseismic activity process, and microseismic energy-frequency power law b value, to achieve constraints on the sampled microseismic event samples in the three dimensions of time, space, and intensity.

[0141] (6) The MOPSO multi-objective particle swarm optimization algorithm was used to sample and predict the occurrence of microseismic events on July 15, 2022. The total frequency and total energy of microseismic events predicted by the LSTM model on that day were used as the number of samples and the total energy in the historical microseismic events in the algorithm. Seven sub-objectives were used as constraints: three capacity dimensions in time, space and energy capacity, total energy, early warning value of microseismic activity process frequency, risk assessment value and b value of GR power law fitting. The sampling dataset was selected from April 16 to July 14, 2022. The MOPSO algorithm was iterated and optimized 100 times. The algorithm calculated multiple objective functions to find the set of all optimal solutions. Figure 8 shows the set of all optimal solutions for the three constraint sampling conditions of energy, risk assessment value and early warning value of microseismic activity process frequency. Figure 9 is a comparison of the location of the three sampling prediction datasets and the actual microseismic data.

[0142] (1) The simulation verification results of microseismic prediction show that the microseismic events predicted by univariate constraint are almost identical to the real microseismic events, which confirms that microseismic events similar to real events can be found from historical data in the three dimensions of time, space and energy. Multivariate constraint can achieve quantitative prediction of microseismic events by adjusting the weights of the reliability functions of the three variables.

[0143] (2) The verification results of the quantitative prediction of microseismic events show that the total energy prediction loss MAE is 0.0827 and MSE is 0.0144 based on the LSTM deep learning model, and the total frequency prediction loss MAE is 0.1672 and MSE is 0.0451, which shows good prediction effect. On this basis, the MOPSO multi-objective particle swarm algorithm is further used to quantitatively predict the three elements of time, space and energy of microseismic events. Among them, the time error is mainly controlled in the range of 0 to 4h, the energy error mainly falls in the range of 0 to 1000J, and the spatial error accounts for an average of 45% in the same location area.

[0144] The above are merely preferred embodiments of the present invention and do not limit the patent scope of the present invention. Any equivalent structural or procedural transformations made based on the content of the present invention's specification and drawings, or direct or indirect applications in other related technical fields, are similarly included within the patent protection scope of the present invention.

Claims

1. A quantitative microseismic prediction method based on the fusion of theoretical model constraints and data-driven approaches, characterized in that, Includes the following steps: Step 1: Obtain raw microseismic data and daily coal seam mining footage during the coal mine working face production period, and statistically obtain relevant derived variables as sample features to form an initial sample dataset; Step 2: Preprocess the initial sample dataset and divide the preprocessed dataset into training and test sets; Step 3: Construct a Long Short-Term Memory (LSTM) network model. Use a time window sliding method based on microseismic autocorrelation to determine the lag time of a single sample and set it as the initial memory days. Set the number of LSTM layers, fully connected layers, and neurons in the model. Use Mean Squared Error (MSE) and Mean Absolute Error (MAE) as loss functions. Adjust the model structure using training set loss and test set loss to obtain the optimal model for predicting the total frequency N of microseismic events on the next day. total and total energy M total ; Step 4: Sample historical microseismic data as samples for quantitative prediction of microseismic events, including the three elements of microseismic time, space and energy; Step 5: Construct a multi-objective optimization model containing seven sub-objective functions: microseismic time fractal dimension, microseismic spatial fractal dimension, microseismic energy fractal dimension, total energy predicted by microseismic activity, microseismic reconstruction mining stress, microseismic activity procedure, and microseismic energy-frequency power law b-value, to achieve constraints on the sampled microseismic events in the three dimensions of time, space, and intensity. Step 6: Find the optimal sample as the predicted value.

2. The method for quantitative prediction of microseismic events based on the fusion of theoretical model constraints and data-driven approaches according to claim 1, characterized in that: In step 1, the raw microseismic data includes the occurrence time and energy of microseismic events, and the relevant derived variables include the daily total frequency, daily total energy, daily average energy, daily maximum energy, daily minimum energy, daily total energy magnitude, daily average energy magnitude, daily maximum energy magnitude, and daily minimum energy magnitude of microseismic events. The daily total frequency, daily total energy, daily average energy, daily maximum energy, and daily minimum energy of microseismic events compared to the previous day; The daily total frequency of microseismic events, daily total energy magnitude, daily average energy magnitude, daily maximum energy magnitude, and daily minimum energy magnitude are the averages of the previous three days.

3. The microseismic quantitative prediction method based on the fusion of theoretical model constraints and data-driven approaches according to claim 2, characterized in that: In step 2, the frequency and energy data are transformed using the following formula to eliminate outliers: M E =lg(E S +10),M N =log2(N Frequency +1) In the formula, ME represents the energy magnitude, ES represents the corresponding energy-related data, MN represents the total frequency level, and NFrequency represents the daily total frequency data. The processed data is then normalized for each column using the maximum-minimum value normalization method. The calculation method is as follows: In the formula, Here, xmax is the normalized value, xmin is the maximum value of the sample data, and x is the original value of the sample. 80% of the preprocessed dataset was divided into the training set and 20% into the test set.

4. The microseismic quantitative prediction method based on the fusion of theoretical model constraints and data-driven approaches according to claim 3, characterized in that: The time window sliding method based on microseismic autocorrelation includes: The correlation between microseismic data and its hysteresis data is calculated using the autocorrelation function (ACF). The calculation method is as follows: In the formula, C ΔI representing microseismic data x at time i i Microseismic data x at time i+ΔI i+ΔI The covariance is given by ΔI, where ΔI is the time lag, N is the total number of data points for which autocorrelation needs to be calculated, and μ is the mean of all data. The autocorrelation coefficient (ACF) between the two sets of microseismic data with a time lag of ΔI is also given. ΔI The calculation method is as follows: In the formula, ACF ΔI Let C0 be the autocorrelation of two sets of microseismic data with a time difference of ΔI, and C0 be the autocovariance when ΔI = 0, i.e., the initial term of the time window x. i The autocovariance; ACF ΔI The value is between -1 and 1, when ACF ΔI Below the upper confidence interval C I Previously, the microseismic data sequences were correlated, when ACF ΔI Reaching C I At the upper limit, ΔI represents the maximum correlation period of the microseismic data. max When the time window lag time is greater than ΔI max At that time, it was assumed that the microseismic data had no correlation, and ΔI was used. max As the initial lag time length for a single sample in the LSTM model; C I The calculation method is as follows: The methods for calculating the mean squared error (MSE) and the mean absolute error (MAE) are as follows: In the formula, N is the number of samples, Y i For the true value, These are predicted values.

5. The microseismic quantitative prediction method based on the fusion of theoretical model constraints and data-driven approaches according to claim 4, characterized in that: In step 4, the reliability function derived based on the Weibull probability distribution is used to verify the microseismic prediction simulation experiment: The probability density functions of the three elements of time, space, and energy In the formula, x is the value of the three elements, x≥0; β is the proportional parameter, β>0; α is the shape parameter, α>0; Cumulative distribution function The reliability function R(Δt) represents the probability that the system will operate normally within a time interval Δt. It is the complement of the cumulative distribution function, specifically expressed as: R(Δt) = 1 - F(Δt; β, α) Right now The absolute values ​​of the pairwise differences between the three key variables of a microseismic event are used as independent variables to measure the failure rate of the entire system. The reliability function between microseismic events is as follows: ΔT=|T i -T j |(i≠j) ΔS=|S i -S j |(i≠j) Δlg E=|lg E i -lg E j |(i≠j) In the formula, ΔT, ΔS, and Δlg E represent the time difference, horizontal distance difference, and energy level difference between microseismic events, respectively; R(ΔT), R(ΔS), and R(Δlg E) represent the failure rates of microseismic events in terms of time difference, horizontal distance difference, and energy level difference, respectively. Construct the reliability objective function: max f(i,j)=ω1·R(ΔT)+ω2·R(ΔS)+ω3·R(Δlg E) In the formula, maxf(i,j) is the reliability between the i-th microseismic event in the real microseismic data and the j-th microseismic data in the historical microseismic event dataset. The larger the value, the greater the reliability between the two microseismic events. ω1, ω2, and ω3 are the weights of the reliability in terms of time difference, distance difference in planar location, and energy level difference, respectively.

6. The microseismic quantitative prediction method based on the fusion of theoretical model constraints and data-driven approaches according to claim 5, characterized in that: The models for the time fractal dimension, spatial fractal dimension, energy fractal dimension, total predicted energy of microseismic events, mining-induced stress of microseismic reconstruction, microseismic activity procedures, and the energy-frequency power law b-value of microseismic events are as follows: In the formula, These represent the capacity dimensions of the microseismic data from the past seven days in terms of temporal fractals, spatial fractals, and microseismic fractals. Add sampling to the microseismic data of the past seven days to predict the capacity dimension of microseismic events in terms of temporal fractal, spatial fractal, and microseismic fractal capacity dimensions; min f(D T ),min f(D S ),min f(D E The values ​​represent the error magnitudes in three dimensions between the recent seven-day microseismic data and the sampled microseismic data; min f(Energy) represents the error magnitude between the predicted total energy of the next day's microseismic events Mtotal and the total energy of the sampled Ntotal microseismic data. i Let i be the energy magnitude of the i-th microseismic event in the sampled microseismic dataset. The total energy of the sampled microseismic event dataset; min f(N) MS This refers to the constraint on the microseismic time element by the early warning value of the frequency of microseismic activity processes quantified based on the correlation analysis of production activities. This is the early warning value for the frequency of microseismic activity processes for the i-th microseismic event in the sampled microseismic dataset. This is the sum of the early warning values ​​for the frequency of microseismic activity processes in the sampled microseismic event set. is the mean of the early warning values ​​for the frequency of microseismic activity processes in the sampled microseismic event set; min f(AI) is the constraint of the impact risk assessment value quantified based on the mining-induced stress reconstruction theory on the microseismic spatial elements, AI i Let i be the risk assessment value of the i-th microseismic event in the sampled microseismic dataset. The sum of the risk assessment values ​​for the sampled microseismic event set. denoted as the mean of the risk assessment values ​​of the sampled microseismic dataset; b0 is the b-value of the microseismic energy-frequency power law fitting of the microseismic data of the past seven days; b1 is the b-value of the microseismic energy-frequency power law fitting of the microseismic data of the past seven days with the addition of the sampled predicted microseismic events; min f(b) is the quantified error between the b-values ​​of the microseismic data of the past seven days and the microseismic energy-frequency power law fitting of the microseismic data with the addition of the sampled microseismic data.

7. The method for quantitative prediction of microseismic events based on the fusion of theoretical model constraints and data-driven approaches according to claim 6, characterized in that: The formula for calculating the b-value of the microseismic energy-frequency power law is as follows: lg N(≥lg E)=ab lg E In the formula, lg E represents the microseismic energy level; N (≥lg E) represents the number of microseismic events with an energy level greater than or equal to lg E; a and b are constants, and their values ​​are calculated using the least squares method. In the formula, m is the total number of energy level divisions; lgEi represents the i-th energy level; Ni represents the actual microseismic number of the i-th energy level; and a and b are the estimated parameters obtained by least squares fitting.

8. The method for quantitative prediction of microseismic events based on the fusion of theoretical model constraints and data-driven approaches according to claim 7, characterized in that: In step 6, the optimal sample is obtained as the predicted value using the multi-objective particle swarm optimization algorithm (MOPSO).

9. The microseismic quantitative prediction method based on the fusion of theoretical model constraints and data-driven approaches according to claim 8, characterized in that: The multi-objective particle swarm optimization (MOPSO) algorithm includes the following steps: 1) Generate an initial particle swarm, each particle containing a position and velocity vector, as well as the corresponding objective function value. At the same time, initialize an external archive to store non-dominated solutions. 2) For each iteration, perform the following steps: ① Update particle velocities based on inertia weights, individual optimal solutions pBest, and global optimal solutions gBest: v id (t+1)=ω·v id (t)+c1·rand1·(pBest id -x id (t))+c2·rand2·(gBest id -x id (t)) In the formula, v id (t) is the velocity of the particle at iteration t, x id (t) is the position of the particle, pBest id It is the individual optimal position of the particle, gBest id ω is the global optimal position, c1 and c2 are learning factors, and rand1 and rand2 are random factors. ② Adjust the particle position according to the updated velocity: x id (t+1)=x id (t)+v id (t+1) ③ Perform non-dominated sorting on the individuals in the population to find the current non-dominated solution set; ④ Compare the current non-dominated solution set with the solutions in the external archive, and update the external archive; ⑤ Select the leader particle gBest from the external archive based on crowding distance and probability; 3) The algorithm terminates when the maximum number of iterations or other termination conditions are reached; otherwise, step 2) is repeated. 4) Output the non-dominated solution set in the external archive as the optimal solution set for the multi-objective function.

Images