An adaptive dynamic monitoring and early warning method for rock burst in coal mines

Abstract

The present application relates to a method for adaptive dynamic monitoring and early warning of rock burst in coal mines. The method comprises: collecting rock burst precursor characteristic indicators of underground working faces, obtaining a sample data set, and performing noise reduction filtering and standardized preprocessing operations to obtain a processed sample data set; using principal component analysis, projecting the processed sample data set into principal component subspace and residual subspace to obtain Hotelling's T 2 and Q statistics to obtain the T corresponding to the sample data set 2 value and Q value; use the Mann-Kendall trend test method to determine the statistical variables of the sample data set; according to the sample data set, statistical variables, and T corresponding to the sample data set 2 The value and Q value are used to analyze whether the impact risk conditions are met. When one or more of the impact risk conditions are met, it is determined that the underground working face currently has an impact risk. This improves the real-time and accuracy of the impact rock pressure monitoring and early warning.

Description

Technical Field

[0001] The present application relates to the technical field of coal-rock dynamic disaster monitoring and early warning, and in particular to a method for adaptive dynamic monitoring and early warning of coal mine rock burst. Background Art

[0002] At present, coal resources are gradually depleted. With the increasing depth of coal mining, the number of rock burst mines and the degree of rock burst danger have increased significantly, seriously restricting the safe production of coal mines. Geophysical methods such as microseismic, acoustic emission, and electromagnetic radiation are typical coal-rock dynamic disaster monitoring and early warning methods. They can directly obtain the precursor characteristic information of rock burst by monitoring the dynamic changes of relevant signals and characteristic indicators generated by the damage and fracture process of coal and rock mass, and effectively achieve the purpose of monitoring and early warning.

[0003] Most of the existing studies on rock burst early warning technology are based on critical value and trend analysis methods. Whether it is microseismic, electromagnetic radiation or acoustic emission signals, there will be an obvious upward trend in the period before rock burst occurs, accompanied by a significant increase in amplitude. When the amplitude of the real-time signal exceeds the upper limit of the critical value or produces a corresponding dynamic change trend, it is judged as a precursor feature of rock burst and a rock burst danger alarm is issued. A large number of studies on precursor features of rock burst have shown that the violent fluctuation before the impact is also one of the significant response characteristics. Even though some scholars use the mean square error, mean value and dynamic change rate of time series as indicators of abnormal fluctuations in research data, there is a lack of quantitative analysis of the correlation of multi-parameter time series features, resulting in many false alarms and missed alarms in the existing rock burst early warning technology, and the early warning accuracy is insufficient. For example, even if a shock occurs, the time series characteristics do not exceed the critical value or produce a corresponding change trend. Some have exceeded the critical value or produced a corresponding change trend, but no shock occurs.

[0004] After searching the literature of the prior art, it was found that the publication date was February 2009, the title of the paper was: Continuous Monitoring and Early Warning Software System for Electromagnetic Radiation of Rock Burst, and the journal information was: Journal of Liaoning University of Engineering and Technology (Natural Science Edition), 2009, 28(01). The self-description of the paper was: "This paper studies the continuous monitoring and early warning software system for electromagnetic radiation of rock burst, including: In view of the large amount of electromagnetic radiation and other signal data in the process of rock burst monitoring, the lack of in-depth analysis and early warning functions, and the development of a continuous monitoring and early warning software system for electromagnetic radiation of rock burst according to the characteristics of rock burst and the characteristics of precursor information. The software system can monitor and analyze electromagnetic radiation, gas and other signals, and has the functions of data storage, real-time display, historical data query, data trend analysis, statistical analysis, early warning and reporting. It adopts critical value method and dynamic trend method for early warning, and realizes real-time monitoring and early warning of rock burst." Its shortcomings are: this early warning method only involves the analysis of the critical values ​​and dynamic trends of time series characteristics, does not consider the impact of data fluctuations on the precursor information and occurrence of rock burst, and lacks research on the correlation of precursor characteristics. As a result, when using existing methods for actual rock burst early warning, it is easy to generate a large number of false alarms and missed alarms, which poses a serious hidden danger to the safe production of coal mines.

[0005] Therefore, the accuracy of current coal mine rock burst monitoring and early warning is low. Summary of the invention

[0006] Based on this, it is necessary to provide a method for adaptive dynamic monitoring and early warning of coal mine rock burst pressure in response to the above technical problems, which can improve the accuracy of coal mine rock burst pressure detection and early warning.

[0007] A method for adaptive dynamic monitoring and early warning of coal mine rock burst, the method comprising:

[0008] Collecting rock burst precursor characteristic indicators of the underground working face to obtain a sample data set, wherein the rock burst precursor characteristic indicators include acoustic emission intensity, acoustic emission pulse, electromagnetic radiation intensity and electromagnetic radiation pulse;

[0009] Performing denoising filtering and standardization preprocessing operations on the sample data set to obtain a processed sample data set;

[0010] The processed sample data set is projected into the principal component subspace and the residual subspace by using the principal component analysis method to obtain the principal component subspace and the residual subspace

[0011] According to the principal element subspace and the residual subspace Analysis of Hotelling's T 2and Q statistics to obtain the T corresponding to the sample data set 2 Value and Q value;

[0012] The Mann-Kendall trend test method is used to determine the statistical variable UF of the sample data set. k ;

[0013] According to the sample data set, the statistical variable, and T corresponding to the sample data set 2 value and Q value, analyzing whether the impact hazard condition is met, and determining that the underground working face currently has an impact hazard when one or more of the impact hazard conditions are met;

[0014] The shock hazard conditions include:

[0015] Condition 1 for the existence of shock risk: T corresponding to the sample data set 2 The value and Q value are greater than the preset threshold value T UCL 2 and / or a preset threshold value Q UCL ;

[0016] The second condition for the existence of shock risk: the original value P in the sample data set k There is a value greater than the preset upper limit of the Pauta criterion critical value P UCL ;

[0017] Condition 3 for the existence of shock risk: Statistical variable UF of monitoring data of Mann-Kendall trend test method k Greater than 0.

[0018] In one embodiment, the expression of the noise reduction filter is:

[0019]

[0020] Among them, x(k+i) is the original value of the k+i th sampling point in the sample data set, N is the number of mean calculation sampling points, N=5, and y(k) is the filtered value corresponding to the k th sampling point.

[0021] In one embodiment, the standardized expression is:

[0022]

[0023] Where X is the processed sample data set, y is the sample data set before standardization, y = [y(1),…y(k)…,y(n)], μ is the mean vector of the sample data set, μ = [μ1,…,μ i …,μ n ],μ iis the mean of the i-th sample data in the sample data set; is the variance matrix of the sample data set, is the variance of the i-th sample data in the sample data set, and n is the number of sample data in the sample data set.

[0024] In one embodiment, the principal component analysis method is used to project the processed sample data set into the principal component subspace and the residual subspace to obtain the principal component subspace and the residual subspace include:

[0025] According to the processed sample data set, a covariance matrix of the processed sample data set is obtained, and the covariance matrix S is expressed as:

[0026]

[0027] Among them, X T is the transpose of the processed sample data set;

[0028] According to the covariance matrix, based on the principal component analysis method, the processed sample data set is projected into the principal component subspace and the residual subspace to obtain the principal component subspace and the residual subspace

[0029] In one embodiment, according to the covariance matrix, based on the principal component analysis method, the processed sample data set is projected into the principal component subspace and the residual subspace to obtain the principal component subspace and the residual subspace include:

[0030] Assuming the number of pivot subspace pivots is l, the expression of the projection process is:

[0031]

[0032] Among them, t i is the i-th principal element of the original random variable, t j is the jth principal element of the original random variable, t i t j =0, is the coordinate of the new linear hyperplane on the principal subspace, is the coordinate of the new linear hyperplane on the residual subspace, is the transformation matrix on the principal subspace, is the transformation matrix on the residual subspace, The transposition of the coordinates of the i-th basis of the new linear hyperplane in the original linear space after the principal component analysis linear transformation, Transpose the coordinates of the j-th basis of the new linear hyperplane in the original linear space after the principal component analysis linear transformation;

[0033] Let x a is a sample of the processed sample data set X sampled once, and its expression is:

[0034]

[0035] in, For the sample x a The projection onto the principal subspace, For the sample x a The projection expressions on the residual subspace are:

[0036]

[0037]

[0038] Among them, the principal element subspace of the processed sample data set The residual subspace of the processed sample data set

[0039] In one embodiment, the T corresponding to the sample data 2 The calculation formula of the value and Q value is:

[0040]

[0041]

[0042] Among them, Λ is the covariance matrix of principal element t, t=[t1,…,t i …,t l ], T 2 T corresponding to the sample data 2 Value, Q is the Q value corresponding to the sample data.

[0043] In one embodiment, the determining the statistical variables of the sample data set by using the Mann-Kendall trend test method includes:

[0044] The null hypothesis H0 is that the rock burst precursor characteristic index y = [y(1), ... y(k) ..., y(n)] is a sample of n independent random variables with the same distribution. The alternative hypothesis H1 is a bilateral test. For all k, j ≤ n, and k ≠ j, the distributions of y(k) and y(j) are different.

[0045] The calculation formula of the test statistic S is:

[0046]

[0047] Among them, Z p , Z q are the observation values ​​corresponding to the pth and qth time series respectively, and p<q, sgn is the sign function;

[0048] The calculation formula of the symbolic function sgn is:

[0049]

[0050] When the Mann-Kendall test is used to test sequence mutations, a sequence S is constructed. k , calculation formula:

[0051]

[0052] in,

[0053] Statistical variables UF k The calculation formula is:

[0054]

[0055] Among them, E(S k )=k(k+1) / 4,Var(S k )=k(k-1)(2k+5) / 72.

[0056] In one embodiment, the preset threshold value T UCL 2 and the preset threshold value Q UCL The method of determining is:

[0057] Collect rock burst precursor characteristic indicators of underground working faces to obtain training sample data sets;

[0058] Performing denoising filtering and standardization preprocessing operations on the training sample data set to obtain a processed training sample data set;

[0059] Dividing the processed training sample data set into multiple training sample data subsets;

[0060] The principal component analysis method is used to project each of the training sample data subsets into the principal component subspace and the residual subspace, respectively, to obtain the principal component subspace of each of the training sample data subsets. and the residual subspace

[0061] According to the principal component subspace of each training sample data subset and the residual subspace Analysis of Hotelling's T 2and Q statistics to obtain T corresponding to each of the training sample data subsets 2 Value and Q value;

[0062] Adopting the adaptive kernel density estimation algorithm, according to the T corresponding to each of the training sample data subsets 2 The value and Q value were analyzed to determine Hotelling's T 2 and the threshold value T of the Q statistic UCL 2 and Q UCL .

[0063] In one embodiment, the preset upper limit P of the Pauta criterion critical value UCL The method of determining is:

[0064] Calculate the mean of the training sample data set:

[0065]

[0066] in, is the average value of the training sample data set, Y i is the i-th training sample data in the training sample data set, and M is the number of training sample data in the training sample data set;

[0067] The standard deviation σ of the training sample data set is:

[0068]

[0069] The premise of Pauta criterion for identifying gross errors is that the monitored time series follows a normal distribution N(μ,σ 2 ) and the amount of data in the sample data set reaches the preset amount, the upper limit of the critical value of the Pauta criterion P UCL The calculation formula is:

[0070]

[0071] Here, σ is the standard deviation.

[0072] The above-mentioned coal mine rock burst adaptive dynamic monitoring and early warning method obtains a sample data set by collecting rock burst precursor characteristic indicators of the underground working face, and the rock burst precursor characteristic indicators include acoustic emission intensity, acoustic emission pulse, electromagnetic radiation intensity and electromagnetic radiation pulse; performs noise reduction filtering and standardized preprocessing operations on the sample data set to obtain a processed sample data set; adopts the principal component analysis method to project the processed sample data set into the principal component subspace and the residual subspace to obtain the principal component subspace and the residual subspace; analyzes Hotelling's T 2and Q statistics to obtain the T corresponding to the sample data set 2 value and Q value; using the Mann-Kendall trend test method to determine the statistical variables of the sample data set; according to the sample data set, the statistical variables, and the T corresponding to the sample data set 2 The value and Q value are used to analyze whether the impact risk conditions are met. When one or more of the impact risk conditions are met, it is determined that the underground working face currently has an impact risk. As a result, the real-time and accuracy of the impact rock pressure monitoring and early warning are improved, so that reasonable anti-impact measures can be taken in a timely manner to ensure the safe production of coal mines. BRIEF DESCRIPTION OF THE DRAWINGS

[0073] Figure 1 A schematic diagram of a process flow of a method for adaptive dynamic monitoring and early warning of rock burst in a coal mine in one embodiment;

[0074] Figure 2 The effect diagram of the implementation of dynamic trend analysis method (Mann-Kendall trend test) in rock burst monitoring and early warning for the acoustic emission intensity index of this application;

[0075] Figure 3 The acoustic emission intensity index of this application adopts the critical value method (Pauta criterion) and multivariate statistical method (principal component analysis) to implement the effect diagram on rock burst monitoring and early warning;

[0076] Figure 4 The effect diagram of the implementation of the dynamic trend analysis method (Mann-Kendall trend test) for the acoustic emission pulse index of this application in the monitoring and early warning of rock burst;

[0077] Figure 5 The acoustic emission pulse index of this application adopts the critical value method (Pauta criterion) and multivariate statistical method (principal component analysis method) to implement the effect diagram on rock burst monitoring and early warning;

[0078] Figure 6 The electromagnetic radiation intensity index of this application adopts the dynamic trend analysis method (Mann-Kendall trend test) to implement the effect diagram of rock burst monitoring and early warning;

[0079] Figure 7 The electromagnetic radiation intensity index of this application adopts the critical value method (Pauta criterion) and multivariate statistical method (principal component analysis) to implement the effect diagram on rock burst monitoring and early warning;

[0080] Figure 8 The electromagnetic radiation pulse index of this application adopts the dynamic trend analysis method (Mann-Kendall trend test) to implement the effect diagram on rock burst monitoring and early warning;

[0081] Fig. 9 This is the effect diagram of the electromagnetic radiation pulse index sampling critical value method (Pauta criterion) and multivariate statistical method (principal component analysis method) of this application in rock burst monitoring and early warning;

[0082] Fig.10 This is the effect diagram of the implementation of Hotelling's statistic in rock burst monitoring and early warning using the multivariate statistical method (principal component analysis) of this application;

[0083] Fig.11 This is the effect diagram of the implementation of Hotelling's statistic, a multivariate statistical method (principal component analysis) in rock burst monitoring and early warning for this application. DETAILED DESCRIPTION

[0084] In order to make the purpose, technical solution and advantages of the present application more clearly understood, the present application is further described in detail below in conjunction with the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are only used to explain the present application and are not used to limit the present application.

[0085] In one embodiment, Figure 1 As shown, a method for adaptive dynamic monitoring and early warning of rock burst in coal mines is provided, comprising the following steps:

[0086] Step 11, collecting rock burst precursor characteristic indicators of the underground working face to obtain a sample data set, the rock burst precursor characteristic indicators include acoustic emission intensity, acoustic emission pulse, electromagnetic radiation intensity and electromagnetic radiation pulse.

[0087] Step 12: Perform denoising filtering and standardization preprocessing operations on the sample data set to obtain a processed sample data set.

[0088] In one embodiment, the expression of the noise reduction filter is:

[0089]

[0090] Among them, x(k+i) is the original value of the k+i-th sampling point in the sample data set, N is the number of mean calculation sampling points, N=5, and y(k) is the filtered value corresponding to the k-th sampling point.

[0091] In one embodiment, the normalized expression is:

[0092]

[0093] Where X is the processed sample data set, y is the sample data set before standardization, y = [y(1),…y(k)…,y(n)], μ is the mean vector of the sample data set, μ = [μ1,…,μ i …,μ n],μ i is the mean of the i-th sample data in the sample data set; is the variance matrix of the sample data set, is the variance of the i-th sample data in the sample data set, and n is the number of sample data in the sample data set.

[0094] In one embodiment, the number of sample data in the sample data set may be determined according to actual conditions, such as: n=4.

[0095] Step 13: Use principal component analysis to project the processed sample data set into the principal component subspace and the residual subspace to obtain the principal component subspace and the residual subspace

[0096] Among them, the principal component analysis method can be used to obtain the correlation changes of the precursor characteristics under normal working conditions. The correlation changes of the precursor characteristics under normal working conditions are obtained through the principal component subspace and the residual subspace reflect.

[0097] In one embodiment, the principal component analysis method is used to project the processed sample data set into the principal component subspace and the residual subspace to obtain the principal component subspace and the residual subspace include:

[0098] According to the processed sample data set, the covariance matrix of the processed sample data set is obtained. The covariance matrix S is expressed as:

[0099]

[0100] Among them, X T is the transpose of the processed sample data set;

[0101] According to the covariance matrix, based on the principal component analysis method, the processed sample data set is projected into the principal component subspace and the residual subspace to obtain the principal component subspace and the residual subspace

[0102] In one embodiment, according to the covariance matrix, based on the principal component analysis method, the processed sample data set is projected into the principal component subspace and the residual subspace to obtain the principal component subspace and the residual subspace include:

[0103] Assuming the number of pivot subspace pivots is l, the expression of the projection process is:

[0104]

[0105] Among them, t i is the i-th principal element of the original random variable, t j is the jth principal element of the original random variable, t i t j =0, is the coordinate of the new linear hyperplane on the principal subspace, is the coordinate of the new linear hyperplane on the residual subspace, is the transformation matrix on the principal subspace, is the transformation matrix on the residual subspace, The transposition of the coordinates of the i-th basis of the new linear hyperplane in the original linear space after the principal component analysis linear transformation, Transpose the coordinates of the j-th basis of the new linear hyperplane in the original linear space after the principal component analysis linear transformation;

[0106] Let x a is a sample of the processed sample data set X sampled once, and its expression is:

[0107]

[0108] in, For the sample x a The projection onto the principal subspace, For the sample x a The projection expressions on the residual subspace are:

[0109]

[0110]

[0111] Among them, the principal element subspace of the processed sample data set The residual subspace of the processed sample data set

[0112] Step 14: According to the principal element subspace and the residual subspace Analysis of Hotelling's T 2 and Q statistics to obtain the T corresponding to the sample data set 2 value and Q value.

[0113] In one embodiment, the sample data corresponds to T 2 The calculation formula of the value and Q value is:

[0114]

[0115]

[0116] Among them, Λ is the covariance matrix of principal element t, t=[t1,…,t i …,t l ], T 2 T corresponding to the sample data 2 Value, Q is the Q value corresponding to the sample data.

[0117] Step 15: Use the Mann-Kendall trend test method to determine the statistical variable UF of the sample data set k .

[0118] In one embodiment, the Mann-Kendall trend test method is used to determine the statistical variables of the sample data set, including:

[0119] The null hypothesis H0 is that the rock burst precursor characteristic index y = [y(1), ... y(k) ..., y(n)] is a sample of n independent random variables with the same distribution. The alternative hypothesis H1 is a bilateral test. For all k, j ≤ n, and k ≠ j, the distributions of y(k) and y(j) are different.

[0120] The calculation formula of the test statistic S is:

[0121]

[0122] Among them, Z p , Z q are the observation values ​​corresponding to the pth and qth time series respectively, and p<q, sgn is the sign function;

[0123] The calculation formula of the symbolic function sgn is:

[0124]

[0125] When the Mann-Kendall test is used to test sequence mutations, a sequence S is constructed. k , calculation formula:

[0126]

[0127] in,

[0128] Statistical variables UF k The calculation formula is:

[0129]

[0130] Among them, E(S k )=k(k+1) / 4,Var(S k )=k(k-1)(2k+5) / 72.

[0131] Step 16: According to the sample data set, statistical variables, and T corresponding to the sample data set 2 The Q value and the Q value are used to analyze whether the conditions for the existence of impact danger are met. When one or more of the conditions for the existence of impact danger are met, it is determined that the underground working face currently has an impact danger.

[0132] Among them, the impact hazard conditions include:

[0133] Condition 1 for the existence of shock risk: T corresponding to the sample data set 2 The value and Q value are greater than the preset threshold value T UCL 2 and / or a preset threshold value Q UCL ;

[0134] The second condition for the existence of shock risk: the original value in the sample data set is greater than the preset upper limit of the Pauta criterion critical value P UCL ;

[0135] Condition 3 for the existence of shock risk: Statistical variable UF of monitoring data of Mann-Kendall trend test method k Greater than 0.

[0136] Among them, the preset threshold value T UCL 2 and the preset threshold value Q UCL Determined by offline training.

[0137] It should be understood that the combined monitoring and early warning of rock burst based on correlation-critical value-dynamic trend analysis is achieved through real-time monitoring of precursor characteristic index data of the underground working face, that is, the original value in the sample data set is greater than the preset upper limit P of the critical value of the Pauta criterion. UCL , or the statistical variable UF of the monitoring data of the Mann-Kendall trend test method k Greater than 0, or T corresponding to the sample data set 2 The value is greater than the preset threshold value T UCL 2 , or the Q value corresponding to the sample data set is greater than the preset threshold value Q UCL , it is determined that there is an impact hazard in this area.

[0138] In one embodiment, the preset threshold value T UCL 2 and the preset threshold value Q UCL The method of determining is:

[0139] Step 21, collecting rock burst precursor characteristic indicators of the underground working face to obtain a training sample data set.

[0140] Step 22: perform denoising filtering and standardization preprocessing operations on the training sample data set to obtain a processed training sample data set.

[0141] In one embodiment, the expression of the noise reduction filter is:

[0142]

[0143] in, is the original value of the k+ith sampling point in the training sample data set, M is the number of sampling points for mean calculation, M=5, is the filtered value corresponding to the kth sampling point.

[0144] In one embodiment, the normalized expression is:

[0145]

[0146] Among them, Y is the standardized training sample data set, is the training sample data set before standardization, is the mean vector of the training sample data set, is the mean of the i-th sample data in the training sample data set; is the variance matrix of the training sample data set, is the variance of the i-th sample data in the training sample data set, and m is the number of sample data in the training sample data set.

[0147] Step 23: divide the processed training sample data set into multiple training sample data subsets.

[0148] The processed training sample data set Y is divided into multiple training sample data subsets Y g , where Y=[Y1,…,Y g ,…,Y G ], G is the number of training sample data subsets.

[0149] Step 24, using principal component analysis, project each training sample data subset into the principal component subspace and the residual subspace, and obtain the principal component subspace of each training sample data subset and the residual subspace

[0150] in, represents the principal element subspace of the g-th training sample data subset, Represents the residual subspace of the g-th training sample data subset.

[0151] In one embodiment, principal component analysis is used to project each training sample data subset into a principal component subspace and a residual subspace, respectively, to obtain a principal component subspace of each training sample data subset. and the residual subspace include:

[0152] According to the training sample data subset, the covariance matrix of the training sample data subset is obtained. The covariance matrix of the g-th training sample data subset is It is expressed as:

[0153]

[0154] Among them, Y g T is the transpose of the g-th training sample data subset, and B is the number of samples in the g-th training sample data subset;

[0155] According to the covariance matrix, based on the principal component analysis method, the training sample data subset is projected into the principal component subspace and the residual subspace to obtain the principal component subspace corresponding to the g-th training sample data subset. and the residual subspace

[0156] In one embodiment, according to the covariance matrix, based on the principal component analysis method, the g-th training sample data subset is projected into the principal component subspace and the residual subspace to obtain the principal component subspace corresponding to the g-th training sample data subset and the residual subspace include:

[0157] Assuming the number of pivot subspace pivots is l, the expression of the projection process is:

[0158]

[0159] Among them, t i is the i-th principal element of the original random variable, t j is the jth principal element of the original random variable, t i t j =0, is the coordinate of the new linear hyperplane on the principal element subspace of the g-th training sample data subset, is the coordinate of the new linear hyperplane on the residual subspace of the g-th training sample data subset, is the transformation matrix on the principal element subspace of the g-th training sample data subset, is the transformation matrix on the residual subspace of the g-th training sample data subset, is the transposition of the coordinates of the i-th basis of the new linear hyperplane formed by the principal component analysis linear transformation of the g-th training sample data subset in the original linear space, is the transposition of the coordinates of the jth basis of the new linear hyperplane formed by the principal component analysis linear transformation of the gth training sample data subset in the original linear space;

[0160] Let y g,b is the g-th training sample data subset Y g The expression of a sample sampled once is:

[0161]

[0162] in, is the sample y g,b The projection onto the principal subspace, is the sample y g,b The projection expressions on the residual subspace are:

[0163]

[0164]

[0165] Among them, the principal element subspace of the g-th training sample data subset is The residual subspace of the g-th training sample data subset B is the number of samples in the g-th training sample data subset.

[0166] Step 25: According to the principal element subspace of each training sample data subset and the residual subspace Analysis of Hotelling's T 2 and Q statistics to obtain the T corresponding to each training sample data subset 2 value and Q value.

[0167] In one embodiment, each training sample data subset corresponds to T 2 The calculation formula of the value and Q value is:

[0168]

[0169]

[0170] Among them, Λ is the covariance matrix of principal element t, t=[t1,…,t i …,t l ], T g 2 is T corresponding to the g-th training sample data subset 2 Value, Q g is the Q value corresponding to the g-th training sample data subset.

[0171] Step 26, using an adaptive kernel density estimation algorithm, according to the T corresponding to each training sample data subset 2 The value and Q value were analyzed to determine Hotelling's T 2 and the threshold value T of the Q statistic UCL 2 and Q UCL .

[0172] It should be understood that the number of sample data in the processed training sample data set is M, the M sample data in the processed training sample data set is divided into G training sample data subsets, and the corresponding Hotelling's T is calculated for each training sample data subset. 2 and Q statistics, we can obtain G Hotelling's T 2 And G Q statistics.

[0173] Among them, the analysis principle of the adaptive kernel density estimation algorithm is:

[0174] Set the initial bandwidth h0, for a set of random variables {u1,u2,...,u i ...,u G}, i∈1,2,…,G,, the local bandwidth factor τ is obtained by calculating each random variable u i , the formula is as follows:

[0175]

[0176] in, is the initial adaptive probability density estimation result, u is the current random variable, and h is the bandwidth determined during the kernel density estimation. K(·) represents the kernel function used to fit the probability distribution of the training sample data. This application selects the Gaussian kernel function, whose formula is:

[0177]

[0178] Among them, e is a natural constant and c is a variable in the Gaussian kernel function.

[0179] By the local bandwidth factor τ i Modify the initial bandwidth h0 to obtain the final adaptive probability density estimation result

[0180] Calculate each random variable u to obtain the local bandwidth factor τ i , the formula is as follows:

[0181]

[0182] Among them, α is the sensitivity factor.

[0183] By the local bandwidth factor τ i Modify the initial bandwidth h0 to obtain the final adaptive probability density estimation result Final adaptive probability density estimation result The expression is:

[0184]

[0185] Among them, at the preset confidence level Next, according to the final adaptive probability density estimation result Derivation of Hotelling's T 2 And the threshold value of Q statistic:

[0186] At a preset confidence level Next, using the final adaptive probability density estimation result Expression, analyze each Hotelling's T 2 The value of and the probability density estimation results corresponding to each Q statistic value are used to estimate the probability density using Hotelling T 2 The analytical formula for the threshold value of the Q statistic is given by each Hotelling's T 2 The corresponding threshold value is determined from the value of and each Q statistic value, Hotelling T 2 The analytical formula for the threshold value of the Q statistic is:

[0187]

[0188]

[0189] It should be understood that for the probability density function Corresponding Hotelling T 2 The analysis of the threshold value of the Q statistic sets the variable u in the adaptive probability density estimation to be the Hotelling T of the g-th training sample data subset. 2 and Q statistics, in G Hotelling's T 2 The Hotelling's T is determined by the analytical formula of its threshold value from the values ​​of and G Q statistics. 2 and the threshold value of the Q statistic.

[0190] In one embodiment, the preset upper limit of the Pauta criterion critical value P UCL The method of determining is:

[0191] Calculate the mean of the training sample data set:

[0192]

[0193] in, is the average value of the training sample data set, Y i is the i-th training sample data in the training sample data set, and M is the number of training sample data in the training sample data set;

[0194] The standard deviation σ of the training sample data set is:

[0195]

[0196] The premise of Pauta criterion for identifying gross errors is that the monitored time series follows a normal distribution N(μ,σ 2 ) and the amount of data in the sample data set reaches the preset amount, the upper limit of the critical value of the Pauta criterion P UCL The calculation formula is:

[0197]

[0198] Here, σ is the standard deviation.

[0199] The above-mentioned coal mine rock burst adaptive dynamic monitoring and early warning method obtains a sample data set by collecting rock burst precursor characteristic indicators of the underground working face, and the rock burst precursor characteristic indicators include acoustic emission intensity, acoustic emission pulse, electromagnetic radiation intensity and electromagnetic radiation pulse; performs noise reduction filtering and standardized preprocessing operations on the sample data set to obtain a processed sample data set; adopts the principal component analysis method to project the processed sample data set into the principal component subspace and the residual subspace to obtain the principal component subspace and the residual subspace; analyzes Hotelling's T 2 and Q statistics to obtain the T corresponding to the sample data set 2 value and Q value; use the Mann-Kendall trend test method to determine the statistical variables of the sample data set; according to the sample data set, statistical variables, and T corresponding to the sample data set 2 The value and Q value are used to analyze whether the impact hazard conditions are met. When one or more of the impact hazard conditions are met, it is determined that the underground working face currently has an impact hazard. As a result, the real-time and accuracy of the impact rock pressure monitoring and early warning are improved, so that reasonable anti-impact measures can be taken in a timely manner to ensure the safe production of coal mines.

[0200] To further illustrate the effect of the adaptive dynamic monitoring and early warning method for rock burst in coal mines of the present application, simulations are performed on various indicators of rock burst precursor characteristic indicators, including acoustic emission intensity, acoustic emission pulse, electromagnetic radiation intensity and electromagnetic radiation pulse, and the corresponding implementation effects are as follows:

[0201] like Figure 3 , 5, 7 and 9, we can see the distribution of the four rock burst precursor characteristic indicators and the corresponding critical value upper limits determined by the Pauta criterion. Figure 3 and Figure 7 There are over-limit situations, and Figure 5 and Fig. 9 None of them exceeded the limit and no abnormality was detected. Figure 2 , 4 , 6 and 8, the Mann-Kendall trend test method, Figure 8 Although the UF value exceeded the limit, it was not obvious and the correct abnormal moment could not be detected. Fig.10 and Fig.11 In the Hotelling's statistic test, T 2 Both Q and Q can obviously exceed the limit when the precursor characteristic indicators of impact are abnormal. Combined with the existing critical value and trend method, the high real-time and accuracy of impact rock pressure monitoring and early warning can be guaranteed.

[0202] It should be understood that although Figure 1 The steps in the flowchart are shown in sequence as indicated by the arrows, but these steps are not necessarily executed in the order indicated by the arrows. Unless otherwise specified in this document, there is no strict order restriction for the execution of these steps, and these steps can be executed in other orders. Moreover, Figure 1 At least part of the steps may include multiple sub-steps or multiple stages. These sub-steps or stages are not necessarily executed at the same time, but can be executed at different times. The execution order of these sub-steps or stages is not necessarily sequential, but can be executed in turn or alternately with other steps or at least part of the sub-steps or stages of other steps.

[0203] The technical features of the above embodiments may be combined arbitrarily. To make the description concise, not all possible combinations of the technical features in the above embodiments are described. However, as long as there is no contradiction in the combination of these technical features, they should be considered to be within the scope of this specification.

[0204] The above-mentioned embodiments only express several implementation methods of the present application, and the descriptions thereof are relatively specific and detailed, but they cannot be understood as limiting the scope of the invention patent. It should be pointed out that, for a person of ordinary skill in the art, several variations and improvements can be made without departing from the concept of the present application, and these all belong to the protection scope of the present application. Therefore, the protection scope of the patent of the present application shall be subject to the attached claims.

Claims

1. A method for adaptive dynamic monitoring and early warning of rock burst in coal mines, characterized in that: The coal mine rock burst adaptive dynamic monitoring and early warning method comprises: Collecting rock burst precursor characteristic indicators of the underground working face to obtain a sample data set, wherein the rock burst precursor characteristic indicators include acoustic emission intensity, acoustic emission pulse, electromagnetic radiation intensity and electromagnetic radiation pulse; Performing denoising filtering and standardization preprocessing operations on the sample data set to obtain a processed sample data set; The processed sample data set is projected into the principal component subspace and the residual subspace using the principal component analysis method to obtain the principal component subspace and the residual subspace According to the principal element subspace and the residual subspace Analysis of Hotelling's T 2 and Q statistics to obtain the T corresponding to the sample data set 2 Value and Q value; The Mann-Kendall trend test method is used to determine the statistical variable UF of the sample data set. k ; According to the sample data set, the statistical variable, and T corresponding to the sample data set 2 value and Q value, analyzing whether the impact hazard condition is met, and determining that the underground working face currently has an impact hazard when one or more of the impact hazard conditions are met; The shock hazard conditions include: Condition 1 for the existence of shock risk: T corresponding to the sample data set 2 The value and Q value are greater than the preset threshold value T UCL 2 and / or a preset threshold value Q UCL ; The second impact risk condition: the original value in the sample data set is greater than the preset upper limit P of the Pauta criterion critical value. UCL ; Condition 3 for the existence of shock risk: Statistical variable UF of monitoring data of Mann-Kendall trend test method k Greater than 0.

2. The method for adaptive dynamic monitoring and early warning of coal mine rock burst according to claim 1 is characterized in that: The expression of the noise reduction filter is: Among them, x(k+i) is the original value of the k+i th sampling point in the sample data set, N is the number of mean calculation sampling points, N=5, and y(k) is the filtered value corresponding to the k th sampling point.

3. The method for adaptive dynamic monitoring and early warning of coal mine rock burst according to claim 2 is characterized in that: The standardized expression is: Where X is the processed sample data set, y is the sample data set before standardization, y = [y(1),…y(k)…,y(n)], μ is the mean vector of the sample data set, μ = [μ1,…,μ i …,μ n ],μ i is the mean of the i-th sample data in the sample data set; is the variance matrix of the sample data set, is the variance of the i-th sample data in the sample data set, and n is the number of sample data in the sample data set.

4. The method for adaptive dynamic monitoring and early warning of coal mine rock burst according to claim 3 is characterized in that: The principal component analysis method is used to project the processed sample data set into the principal component subspace and the residual subspace to obtain the principal component subspace. and the residual subspace include: According to the processed sample data set, a covariance matrix of the processed sample data set is obtained, and the covariance matrix S is expressed as: Among them, X T is the transpose of the processed sample data set; According to the covariance matrix, based on the principal component analysis method, the processed sample data set is projected into the principal component subspace and the residual subspace to obtain the principal component subspace and the residual subspace 5. The method for adaptive dynamic monitoring and early warning of rock burst in coal mines according to claim 4 is characterized in that: According to the covariance matrix, based on the principal component analysis method, the processed sample data set is projected into the principal component subspace and the residual subspace to obtain the principal component subspace and the residual subspace include: Assuming the number of pivot subspace pivots is l, the expression of the projection process is: Among them, t i is the i-th principal element of the original random variable, t j is the jth principal element of the original random variable, t i t j =0, is the coordinate of the new linear hyperplane on the principal subspace, is the coordinate of the new linear hyperplane on the residual subspace, is the transformation matrix on the principal subspace, is the transformation matrix on the residual subspace, The transposition of the coordinates of the i-th basis of the new linear hyperplane in the original linear space after the principal component analysis linear transformation, Transpose the coordinates of the j-th basis of the new linear hyperplane in the original linear space after the principal component analysis linear transformation; Let x a is a sample of the processed sample data set X sampled once, and its expression is: in, For the sample x a The projection onto the principal subspace, For the sample x a The projection expressions on the residual subspace are: Among them, the principal element subspace of the processed sample data set The residual subspace of the processed sample data set 6. The method for adaptive dynamic monitoring and early warning of rock burst in coal mines according to claim 5, characterized in that: The sample data corresponds to T 2 The calculation formula of the value and Q value is: Among them, Λ is the covariance matrix of principal element t, t=[t1,…,t i …,t l ], T 2 T corresponding to the sample data 2 Value, Q is the Q value corresponding to the sample data.

7. The method for adaptive dynamic monitoring and early warning of rock burst in coal mines according to claim 1, characterized in that: The method of using the Mann-Kendall trend test method to determine the statistical variables of the sample data set includes: The null hypothesis H0 is that the rock burst precursor characteristic index y = [y(1), ... y(k) ..., y(n)] is a sample of n independent random variables with the same distribution. The alternative hypothesis H1 is a bilateral test. For all k, j ≤ n, and k ≠ j, the distributions of y(k) and y(j) are different. The calculation formula of the test statistic S is: Among them, Z p , Z q are the observation values ​​corresponding to the pth and qth time series respectively, and p<q, sgn is the sign function; The calculation formula of the symbolic function sgn is: When the Mann-Kendall test is used to test sequence mutations, a sequence S is constructed. k , calculation formula: in, Statistical variables UF k The calculation formula is: Among them, E(S k )=k(k+1) / 4,Var(S k )=k(k-1)(2k+5) / 72.

8. The method for adaptive dynamic monitoring and early warning of rock burst in coal mines according to claim 1, characterized in that: The preset threshold value T UCL 2 and the preset threshold value Q UCL The method of determining is: Collect rock burst precursor characteristic indicators of underground working faces to obtain training sample data sets; Performing denoising filtering and standardization preprocessing operations on the training sample data set to obtain a processed training sample data set; Dividing the processed training sample data set into multiple training sample data subsets; The principal component analysis method is used to project each of the training sample data subsets into the principal component subspace and the residual subspace, respectively, to obtain the principal component subspace of each of the training sample data subsets. and the residual subspace According to the principal component subspace of each training sample data subset and the residual subspace Analysis of Hotelling's T 2 and Q statistics to obtain T corresponding to each of the training sample data subsets 2 Value and Q value; Adopting the adaptive kernel density estimation algorithm, according to the T corresponding to each of the training sample data subsets 2 The value and Q value were analyzed to determine Hotelling's T 2 and the threshold value T of the Q statistic UCL 2 and Q UCL .

9. The method for adaptive dynamic monitoring and early warning of rock burst in coal mines according to claim 1, characterized in that: The preset upper limit of the critical value of the Pauta criterion P UCL The method of determining is: Calculate the mean of the training sample data set: in, is the average value of the training sample data set, Y i is the i-th training sample data in the training sample data set, and M is the number of training sample data in the training sample data set; The standard deviation σ of the training sample data set is: The premise of Pauta criterion for identifying gross errors is that the monitored time series follows a normal distribution N(μ,σ 2 ) and the amount of data in the sample data set reaches the preset amount, the upper limit of the critical value of the Pauta criterion P UCL The calculation formula is: Here, σ is the standard deviation.

Images