Pipe leakage aperture identification method based on apso-hmm

Abstract

The present application relates to the technical field of pipeline leakage, in particular to a pipeline leakage aperture identification method based on APSO-HMM. The method comprises: constructing multiple training sets, each corresponding to a pipeline leakage aperture; inputting each training set into an HMM model for training iteration, and using APSO algorithm to optimize the initial value of the observation value under each state, i.e. the initial value of B, in the training process; bringing the optimal initial value obtained by the APSO algorithm into the HMM model, and using the Baum-Welch algorithm to iteratively calculate the HMM model to obtain the trained APSO-HMM model; extracting the feature vector T of the infrasound wave original signal leakage signal x(t) to be identified, inputting the feature vector T into each trained APSO-HMM model respectively, and outputting the pipeline leakage aperture corresponding to the APSO-HMM model with the maximum log-likelihood probability as the identification result. The present application can avoid improper selection of the initial value of the model, thereby improving the identification effect of the pipeline leakage aperture.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of pipeline leakage, in particular to a pipeline leakage aperture identification method based on APSO-HMM. BACKGROUND

[0002] Most of the urban pipe networks are distributed in a network shape, involve a wide range and have a large number of pipe network nodes. For the small leakage signals of the urban buried pipeline, the collected signals are mixed with a large amount of environmental noise, medium noise and other non-related signals. The small pipeline leakage, especially the persistent pinhole leakage, is difficult to detect. Therefore, how to effectively detect and locate the initial small leakage of the pipeline has good economic value and social significance.

[0003] The hidden Markov model (HMM) is a time series model, which is a non-memory non-stationary random process. The hidden Markov model has a very strong ability to process the statistical information and time sequence information contained in the sound signal, can well represent the time-varying signal, and can establish the hidden Markov model under different conditions according to different signal characteristic quantities of the pipeline leakage, so as to classify and identify the pipeline leakage and has a good recognition rate, especially for the small persistent leakage.

[0004] The Baum-Welch algorithm is an important method for solving the parameters of the hidden Markov model and is used when the hidden sequence cannot be obtained. However, the existing HMM model has a poor recognition effect on the pipeline. SUMMARY

[0005] The technical problem to be solved by the present application is to overcome the defects of the prior art and provide a pipeline leakage aperture identification method based on APSO-HMM, which can avoid improper selection of the initial value of the model, so as to improve the recognition effect of the pipeline leakage aperture.

[0006] In order to solve the above technical problems, the technical scheme of the present application is as follows: a pipeline leakage aperture identification method based on APSO-HMM, comprising:

[0007] A plurality of training sets are constructed, and each training set corresponds to a pipeline leakage aperture. The construction process of each training set comprises: obtaining a plurality of infrasound original leakage signals x(t) corresponding to the pipeline leakage aperture; extracting a feature vector T of each infrasound original leakage signal x(t), and a set formed by the feature vectors T of the plurality of infrasound original leakage signals x(t) is used as a training set of the HMM model for the pipeline leakage aperture;

[0008] Each training set is input into an HMM model for training iteration, and the APSO algorithm is used to optimize the initial value of the observation value under each state, i.e. the initial value of B, in the training process; the optimal initial value obtained by the APSO algorithm is brought into the HMM model, and the Baum-Welch algorithm is used to iteratively calculate the HMM model, so as to obtain the trained APSO-HMM model;

[0009] The feature vector T of the to-be-identified infrasonic wave original signal leakage signal x(t) is extracted, the feature vector T is input into each trained APSO-HMM model, and the pipeline leakage aperture corresponding to the APSO-HMM model with the maximum log-likelihood probability is output as the identification result;

[0010] The step of extracting the feature vector T of the infrasonic wave original leakage signal x(t) is as follows: the infrasonic wave original leakage signal x(t) is subjected to Gaussian low-pass filtering denoising processing to obtain a denoised signal x'(t); the denoised signal x'(t) is decomposed by using the EEMD algorithm, and the effective IMF component is selected from the IMFs obtained by decomposition based on the correlation coefficient; and the energy entropy of the effective IMF component is taken as an element to form the feature vector T of the infrasonic wave original leakage signal x(t).

[0011] Further, the denoised signal x'(t) is decomposed by using the EEMD algorithm; specifically:

[0012] Step 11, Gaussian white noise is introduced into the denoised signal x'(t), and the calculation formula is as follows:

[0013] e i (t) is the signal obtained after adding the i th noise, n i (t)

[0014] e i (t) is the signal obtained after adding the i th noise, n i (t) is the i th white noise signal;

[0015] Step 12, according to the calculation method in EMD, the signal after adding white noise is decomposed by EMD to obtain a group of IMFs:

[0016]

[0017] e ij (t) is the j th IMF component obtained by decomposing the signal e i (t) added with the i th noise, r ij (t) is the j th residual component obtained by decomposing the signal added with the i th noise;

[0018] Step 13, return to step 11 and repeat N times;

[0019] Step 14, the mean of each IMF component is calculated to obtain the final decomposition result, and the calculation formula is:

[0020]

[0021] In the formula, m j (t) is the jth IMF component obtained by EEMD decomposition, and r(t) is the residual component obtained by EEMD decomposition.

[0022] Further, the effective IMF component is selected from the IMF components obtained by decomposition based on the correlation coefficient; specifically,

[0023] The correlation coefficient between each IMF component and the denoised signal x'(t) is calculated using the correlation coefficient formula;

[0024]

[0025] In the formula, cov(x, y) is the covariance of x and y; D(x) and D(y) are the variances of x and y, respectively; and p xy The value range is [0, 1], and x and y represent the IMF component to be calculated and the denoised signal x'(t), respectively;

[0026] Each IMF component is sorted in descending order of the correlation coefficient, and the first J IMF components are taken as the effective IMF components, wherein the sum of the energy values of the first J-1 IMF components is less than a preset percentage of the total energy value of the denoised signal x'(t), and the energy values of the first J IMF components are not less than a preset percentage of the total energy value of the denoised signal x'(t).

[0027] Further, the calculation formula of the energy entropy is:

[0028]

[0029] In the formula, m j (t) is the jth effective IMF component, E j is the energy entropy of the jth effective IMF component.

[0030] Further, the feature vector T is normalized; wherein,

[0031] In the training process, after the normalization processing of each feature vector T in the training set, the HMM model is inputted;

[0032] In the recognition process, after the normalization processing of the feature vector T of the to-be-recognized infrasound wave original signal leakage signal x(t), each trained APSO-HMM model is inputted;

[0033] T = [E1, E2,..., E J ];

[0034]

[0035] T' = [E1 / E, E2 / E,..., E J / E]

[0036] wherein E1, E2,..., E J represent the energy entropy of the first to the Jth effective IMF component, and T' represents the normalized vector.

[0037] Further, the APSO algorithm is used to optimize the initial value of the observation value in each state, i.e., the initial value of B; specifically,

[0038] Step 21, first define the number of groups as M, and each particle searches and optimizes in a D-dimensional space at a certain speed, and the optimal position of the nth particle is defined as pbest n ; the group optimal position is defined as gbest, in each iteration process, each particle and group find its optimal position, and the speed and position are updated according to the formula, thereby completing an iteration process;

[0039] The speed updating formula is:

[0040]

[0041] The position updating formula is:

[0042]

[0043] wherein: is the speed and position of the nth particle in the d-dimensional space in the tth iteration process; is the speed and position of the nth particle in the d-dimensional space in the t-1th iteration process; is the optimal position of the nth particle in the t-1th iteration process; is the optimal position of the group in the t-1th iteration process; w is the inertia weight, c1 and c2 are learning factors, and r1 and r2 are random numbers in [0, 1]; the following variables are defined:

[0044]

[0045] wherein: is the current global optimal value, is the last global optimal value, and α t is the population evolution speed: the ratio of the current global optimal value to the last global optimal value, and α tThe smaller the value of w, the faster the population evolution speed, w should be increased to improve the global search performance of the population, and vice versa; is the average fitness of the current population, β t is the population convergence speed factor at the tth iteration, representing the population convergence degree: the ratio of the global optimum to the average fitness of the current population, the greater the value, the more likely the population falls into local optimum, w should be increased to make the particles diverge, so that the population is easy to jump out of the local area, and vice versa; γ t is the population premature convergence factor at the tth iteration, representing the value of γ t , the smaller the value, the more likely the population tends to premature convergence;

[0046] Step 22, when w shows linear or nonlinear decline, the execution efficiency of the algorithm is higher, the mapping relationship between w and γ is designed as

[0047]

[0048] In the formula: w t is the inertia weight in the tth iteration;

[0049] When the number of iterations is less than the preset value T rec , that is, t≤T rec , w always keeps a larger value to make the particles search in a larger range, and the expansion method is:

[0050]

[0051] When the preset iteration value is reached, the inertia weight factor should be reduced adaptively, that is, t>T rec ,

[0052]

[0053] In the formula: w min and w max are the upper and lower bounds of the inertia weight factor, q1, q2 are empirical constants, T max is the maximum number of iterations, τ is the threshold value;

[0054] Step 23, in order to facilitate analysis and expression, the problem space is simplified to one dimension, define pbest i = p i , gbest = g i , then the speed and position update formula is

[0055]

[0056] Through the iterative relationship, we get

[0057]

[0058] wherein: is the velocity and position of the i-th particle in d-dimension at the t-th iteration; is the velocity and position of the i-th particle in d-dimension at the t+1-th iteration; is the velocity and position of the i-th particle in d-dimension at the t+2-th iteration;

[0059] Combining the above, the second-order difference equation for the particle velocity is obtained as:

[0060]

[0061] Stability analysis of the equation gives the following adaptive factor:

[0062]

[0063] wherein, is the value of the learning factor c1, c2 at the t-th iteration;

[0064] Step 24: The optimal solution B is obtained after the improvement of the weight factor in step 22 and the learning factor in step 23.

[0065] Further, the step of obtaining the log-likelihood probability by the model is:

[0066] Step 41: Local state initialization:

[0067] δ t (i) = π i b i (O1), i = 1, 2, 3,... N

[0068] wherein, δ t (i) is the maximum probability of all transition paths of hidden state i at t time.

[0069]

[0070] wherein, is the previous node hidden state of the maximum probability of all transition paths of hidden state i at t time;

[0071] Step 42: Perform dynamic programming recursion for the local state at time t (t = 1, 2, 3,... T);

[0072]

[0073] Step 43: Solve the maximum probability of δ T (i) at T time, and solve the maximum δ T (i) at T time

[0074] P * = max 1≤j≤N delta T (i)

[0075]

[0076] Step 44, by Start backtracking, that is

[0077]

[0078] The most probable hidden state sequence is obtained

[0079]

[0080] The recursive expression of delta is

[0081]

[0082] The recursive expression of estimated is

[0083]

[0084] Finally, the maximum probability value P* of the last time T is calculated, that is, the log-likelihood probability.

[0085] After the above technical solution, the adaptive particle swarm optimization algorithm is used to optimize and solve the parameter B, which can well solve the problem of fast convergence during search, obtain the optimal solution, avoid the problem of premature convergence, and thus improve the recognition effect of the pipeline leakage aperture. BRIEF DESCRIPTION OF DRAWINGS

[0086] Figure 1 is an experimental field diagram of the present application;

[0087] Figure 2 is a pipeline leakage simulation test schematic diagram of the present application;

[0088] Figure 3 is a flow chart of the pipeline leakage aperture recognition method based on APSO-HMM of the present application;

[0089] Figure 4 is a flow chart of extracting the characteristic vector T from the infrasound wave original leakage signal x(t) of the present application;

[0090] Figure 5 is a time domain waveform diagram of an infrasound wave original leakage signal x(t);

[0091] Figure 6 is Figure 4The time-domain waveform diagram of the infrasound wave original leakage signal x(t) after denoising in the middle;

[0092] Figure 7 The IMF component correlation coefficient diagram obtained by performing EEMD analysis on the denoised signal;

[0093] Figure 8 The initial value diagram of the initial value of the observation value in each state, that is, the initial value of B, is optimized by using an adaptive particle swarm optimization (APSO) algorithm;

[0094] Figure 9 The differences in model fitness convergence value and convergence speed between the APSO-HMM model and the conventional HMM model under a 1mm aperture are compared;

[0095] Figure 10 The differences in model fitness convergence value and convergence speed between the APSO-HMM model and the conventional HMM model under a 2mm aperture are compared;

[0096] Figure 11 The differences in model fitness convergence value and convergence speed between the APSO-HMM model and the conventional HMM model under a 3mm aperture are compared. DETAILED DESCRIPTION

[0097] In order to make the content of the present application more easily understood, the present application will be further described in detail below according to specific embodiments and in conjunction with the accompanying drawings.

[0098] As shown in Figure 3 , a pipeline leakage aperture identification method based on APSO-HMM includes:

[0099] Firstly, a plurality of training sets are constructed, and each training set corresponds to a pipeline leakage aperture; wherein the construction process of each training set includes: obtaining a plurality of infrasound wave original leakage signals x(t) corresponding to the pipeline leakage aperture; extracting a feature vector T of each infrasound wave original leakage signal x(t); and taking a set formed by the feature vectors T of the plurality of infrasound wave original leakage signals x(t) as a training set of an HMM model for the pipeline leakage aperture;

[0100] Secondly, each training set is input into an HMM model for training iteration, and the initial value of the observation value in each state, that is, the initial value of B, is optimized by using an APSO algorithm in the training process; the optimal initial value obtained by the APSO algorithm is brought into the HMM model, and the HMM model is iteratively calculated by using a Baum-Welch algorithm to obtain a trained APSO-HMM model;

[0101] The third step is to extract the feature vector T of the original infrasound signal leakage x(t) to be identified. The feature vector T is input into each trained APSO-HMM model, and the pipeline leakage aperture corresponding to the APSO-HMM model with the largest log-likelihood probability is output as the identification result.

[0102] Among them, such as Figure 4 As shown in FIG, the steps of extracting the characteristic vector T of the original infrasound leakage signal x(t) are as follows: performing Gaussian low-pass filtering denoising processing on the original infrasound leakage signal x(t) to obtain the denoised signal x'(t); since the EEMD algorithm can completely decompose the non-stationary signal into independent IMF signals of different scales, the decomposed signal still has good signal integrity, therefore, the EEMD algorithm is used to decompose the denoised signal x'(t); based on the correlation coefficient, the effective IMF component is selected from the decomposed IMF components; and the energy entropy of the effective IMF component is used as an element to form the characteristic vector T of the original infrasound leakage signal x(t).

[0103] It's important to note that the leak aperture is continuous and can be set to a number of point values ​​as needed. During the training set construction process, the point value to which the actual leak aperture value is closest is considered the leak aperture. After training, during the identification process, the APSO-HMM model with the largest log-likelihood probability corresponds to the point value of the pipeline leak aperture that is closest to the actual leak aperture value and is therefore considered the identified leak aperture.

[0104] Specifically, the Baum-Welch algorithm is highly dependent on the initial value λ = (π, A, B). Until now, the user has provided empirical estimates of each parameter at the initialization of the hidden Markov model. Through continuous iteration, the parameters gradually converge to more reasonable, optimal values. The inventors of this patent have discovered that parameters π and A have little influence on the selection of initial values ​​and can be selected randomly or uniformly, as long as the random constraints are met. However, the selection of B directly affects the accuracy of model training. In this embodiment, to avoid improper initial value selection for model parameters, an adaptive particle swarm optimization algorithm is used to optimize parameter B. This improves existing pipeline leak identification methods, overcoming the shortcomings of traditional hidden Markov search methods, which rapidly converge to the optimal solution and reduce the accuracy of leak signals due to premature convergence, thereby improving the accuracy of urban pipeline leak identification, preventing problems from arising, and safeguarding urban public safety. This improves the current situation of inaccurate urban pipeline leak identification, reduces potential energy waste and harm, minimizes urban environmental damage and pollution, conserves energy, and improves economic efficiency. Identify pipeline leaks in a timely manner and take measures to control them, thereby eliminating potential accident hazards, reducing or controlling the occurrence of pipeline accidents, protecting urban public safety, and reducing the loss of life and property of urban residents.

[0105] In one embodiment, the infrasound wave original leakage signal x(t) is subjected to Gaussian low-pass filter denoising processing, wherein the expression formula of the transfer function of the Gaussian low-pass filter is:

[0106]

[0107] In the formula, G(x, y) is the transfer function of the Gaussian low-pass filter; δ is the standard deviation of the normal distribution, and the size of δ will affect the final processing effect. The larger δ is, the more obvious the processing effect on the noise will be, but the image will also be more blurred. This value is preferably selected by comparison in actual application; (x, y) is the coordinate of the pixel point, and is the point coordinate, which can be considered as an integer in image processing.

[0108] In one embodiment, EEMD is improved on the basis of the EMD method, and is a method of assisting data decomposition by adding noise signals. The method introduces Gaussian white noise into the original signal, and the noise will be eliminated after multiple average calculations because its mean value is zero. This realizes the automatic distribution of the signal on a suitable time scale, and the final integrated mean value is the calculation result, effectively solving the mode aliasing problem. The EEMD method shows better feature extraction capability in processing non-stationary and nonlinear complex time series data. The superposition effect of the decomposed data, that is, the superposition of the original data and the Gaussian white noise, is calculated by the ensemble average of the obtained IMF components, so that the added Gaussian white noise can be offset. When the iteration number is large enough, it can be regarded as the decomposition result of the original signal. The EEMD algorithm is used to decompose the denoised signal x'(t); specifically:

[0109] Step 11, Gaussian white noise is introduced into the denoised signal x'(t), and the calculation formula is:

[0110] e i (t)=x'(t)+n i (t)(2)

[0111] In the formula, e i (t) is the signal obtained after adding the i-th noise, n i (t) is the i-th white noise signal;

[0112] Step 12, according to the calculation method in EMD, the signal after adding white noise is subjected to EMD decomposition to obtain a group of IMF components:

[0113]

[0114] In the formula, m ij (t) is the j-th IMF component obtained by EMD decomposing the signal e i (t) added with the i-th noise, and r ij(t) is the jth residual component obtained by decomposing the signal with the ith added noise;

[0115] Step 13, return to step 11, repeat N times;

[0116] Step 14, in order to offset the effect of the added Gaussian white noise in each calculation, the mean of each IMF component is calculated to obtain the final decomposition result, and the calculation formula is:

[0117]

[0118] In the formula, m j (t) is the jth IMF component obtained by EEMD decomposition, and r(t) is the residual component obtained by EEMD decomposition.

[0119] After the infrasound wave signal is decomposed by the EEMD algorithm, a plurality of different scale independent IMFs are obtained, each IMF has different physical meaning in its own scale, in order to avoid misleading the recognition result caused by the redundant signal in the signal and reduce the redundancy of the signal, the decomposed IMF signal needs to be reasonably selected. In an embodiment, based on the correlation coefficient, the effective IMF component is selected from the decomposed IMF component, and the redundant component is removed; specifically,

[0120] The correlation coefficient between each IMF component and the denoised signal x'(t) is calculated by using the correlation coefficient formula;

[0121]

[0122] In the formula, cov(x, y) is the covariance of x and y; D(x) and D(y) are the variances of x and y, respectively; and p xy The value range is [0, 1], the correlation coefficient is close to 0, indicating low correlation; the correlation coefficient is close to 1, indicating high correlation; x and y represent the IMF component to be calculated and the denoised signal x'(t), respectively;

[0123] Each IMF component is sorted according to the correlation coefficient from large to small, and the first J IMF components are taken as the effective IMF components, wherein the sum of the energy values of the first J-1 IMF components is less than a preset percentage of the total energy value of the denoised signal x'(t), and the energy values of the first J IMF components are not less than the preset percentage of the total energy value of the denoised signal x'(t). The preset percentage can be but is not limited to 90%.

[0124] The calculation formula of the energy entropy is:

[0125]

[0126] In the formula, m j (t) is the jth effective IMF component, and Ej is the energy entropy of the jth effective IMF component.

[0127] In one embodiment, considering that the energy is sometimes large, in order to facilitate subsequent analysis and processing, T is adjusted by normalizing the features. Therefore, the method further includes normalizing the feature vector T; wherein,

[0128] During the training process, each feature vector T in the training set is normalized and then input into the HMM model training;

[0129] During the recognition process, the feature vector T of the original infrasound signal leakage signal x(t) to be identified is normalized and then input into each trained APSO-HMM model;

[0130] T=[E1,E2,...,E J ];

[0131]

[0132] T'=[E1 / E,E2 / E,...,E J / E] (8)

[0133] Where E1, E2, ..., E J represents the energy entropy of the 1st to Jth effective IMF components, and T' represents the normalized vector.

[0134] In one embodiment, the APSO algorithm is used to optimize the initial value of the observation value under each state, that is, the initial value of B; specifically,

[0135] Step 21: First, define the number of particles in the population as M. Each particle searches for the optimal position in the D-dimensional space at a certain speed. The optimal position of the nth particle is defined as pbest. n The optimal position of the group is defined as gbest. In each iteration, each particle and group finds its optimal position and updates its speed and position according to the formula, thus completing an iteration process.

[0136] The speed update formula is:

[0137]

[0138] The position update formula is:

[0139]

[0140] Where: is the velocity and position of the nth particle in the tth iteration in d dimensions; The velocity of the d-dimension and the position of the nth particle in the t-1th iteration; The optimal position of the nth particle in the t-1th iteration; The optimal position of the population in the t-1th iteration; w is the inertia weight, c1, c2 are learning factors, and r1, r2 ∈ [0, 1] are random numbers.

[0141] In practical applications, premature convergence and improper selection of fixed parameters w and c1, c2 can cause the optimization process to fall into local optimum, affecting the performance of the PSO algorithm. Therefore, an adaptive parameter mechanism will be used to overcome this problem. Define the following variables:

[0142]

[0143] In the formula: is the current global optimal value, is the last global optimal value, and α t is the population evolution speed: the ratio of the current global optimal value to the last global optimal value, α t The smaller the value of α, the faster the population evolution speed, and w should be increased to improve the global search performance of the population, and vice versa; is the average fitness of the current population, β t is the population aggregation speed factor at the tth iteration, representing the degree of population aggregation: the ratio of the global optimal value to the average fitness of the current population. The larger the value, the more likely the population is to fall into local optimum, and w should be increased to make the particles diverge, so that the population can easily jump out of the local area, and vice versa; γ t is the premature convergence factor of the population at the tth iteration, representing γ t The smaller the value, the more likely the population is to tend to premature convergence;

[0144] When w shows linear or nonlinear decline, the execution efficiency of the algorithm is higher, and the mapping relationship between w and γ is designed as

[0145]

[0146] In the formula: w t is the inertia weight at the tth iteration;

[0147] When the number of iterations is less than the preset value T rec , that is, t ≤ T rec , w always maintains a larger value to make the particles search and expand in a larger range, and the expansion method is:

[0148]

[0149] When the preset iteration value is reached, the inertia weight factor should be adaptively reduced, that is, t > Trec Time,

[0150]

[0151] where w min and w max are the upper and lower bounds of the inertia weight factor, q1, q2 are empirical constants, T max is the maximum number of iterations, and τ is a threshold value, τ = 0.002.

[0152] Step 23, in the PSO algorithm, the changes of self and social learning factors are not synchronized, c1 increases and c2 decreases at the beginning of iteration, and the population has strong "self-cognition" ability. In the last stage of iteration, c1 decreases and c2 increases. This embodiment will adopt the adaptive learning factor improvement method. In order to facilitate analysis and expression, the problem space is simplified to one dimension, and pbest i = p i , gbest = g i , then the speed and position update formula is

[0153]

[0154] Through the iterative relationship, combined with formula (15), we get

[0155]

[0156] where: is the speed and position of the i-th particle in the d-dimension in the t-th iteration process; is the speed and position of the i-th particle in the d-dimension in the t+1-th iteration process; is the speed and position of the i-th particle in the d-dimension in the t+2-th iteration process;

[0157] Combined with formula (15) and formula (16), the second-order difference equation about the particle speed can be obtained:

[0158]

[0159] Stability analysis of the equation can obtain the following adaptive factor:

[0160]

[0161] where: is the value of learning factor c1, c2 in the t-th iteration process;

[0162] Step 24, after the improvement of the weight factor in step 22 and the learning factor in step 23, the optimal solution B is obtained.

[0163] In one embodiment, the step of obtaining the log-likelihood probability by the model is:

[0164] Step 41, local state initialization:

[0165] δ t (i) = π i b i (O1), i = 1, 2, 3,... N (20)

[0166] In the formula, δ t (i) is the maximum probability in all transition paths of hidden state i at time t.

[0167]

[0168] In the formula, δ is the previous node hidden state of the maximum probability in all transition paths of hidden state i at time t;

[0169] Step 42, performing dynamic programming recursion of local state at time t (t = 1, 2, 3,... T);

[0170]

[0171] Step 43, solving the maximum probability of δ T (i) at time T, and solving the maximum δ T (i) at time T.

[0172] P * = max 1≤j≤N δ T (i)

[0173]

[0174] Step 44, starting backtracking by , that is,

[0175]

[0176] obtaining the hidden state sequence with the maximum probability

[0177]

[0178] The recursive expression of δ is

[0179]

[0180] The recursive expression of is

[0181]

[0182] The maximum probability value P* calculated at the last moment T is the log-likelihood probability.

[0183] The pipeline leakage aperture identification method based on APSO-HMM involved in the above embodiment is described in detail below with reference to the accompanying drawings of laboratory simulation.

[0184] Figure 1 The experimental site map is as follows: The experimental pipeline system is a U-shaped PE pipe, the pipe is Φ110×10mm, and the total length is 3599cm, where the air pressure is 0.1MPa and the flow rate is 16m / s; in the straight section on the left side of the pipeline, an inlet ball valve, an upstream infrasonic sensor, a pressure sensor, and a turbine flowmeter are arranged in sequence from the pipeline inlet end; in the straight section on the right side of the pipeline, an outlet ball valve, a downstream infrasonic sensor, a pressure sensor, a turbine flowmeter and a leakage port B are arranged in sequence from the pipeline outlet end; in the curved section of the pipeline (the middle horizontal section), a pressure sensor and a leakage port A are arranged in sequence along the flow direction; adjustable leakage ports A and B with apertures of 1mm, 2mm and 3mm respectively are used to simulate small leaks in urban pipelines, and an air compressor is used to simulate the operation process of urban gas pipelines. Upstream and downstream infrasonic sensors are installed at both ends of the pipeline, and connected to the digital network transmission instrument through a cable, and then transmitted to the analysis software on the PC through the network cable, see Figure 2 .

[0185] A pipeline leakage aperture identification method based on APSO-HMM includes the following steps:

[0186] S1: Under the condition of 0.3MPa pressure, the leakage pipe is tested. When the leakage pipe is tested, the leakage aperture size is changed, and the pipeline infrasound signal data of different leakage apertures under 0.3MPa pressure are collected. The leakage aperture sizes are 1mm, 2mm, and 3mm respectively. The waveform of the original infrasound signal data collected under the leakage condition is as follows: Figure 5 shown.

[0187] S2: performing Gaussian low-pass filtering and denoising on the original signal of the pipeline leakage infrasound wave to obtain a denoised signal; Figure 5 The waveform in the denoising is as follows Figure 6 shown.

[0188] S3: The EEMD algorithm is used to completely decompose the denoised signal in step S2 into independent IMF signals of different scales. The size of the energy entropy can reflect the change of the signal energy with the frequency distribution. Therefore, the IMF energy entropy of the infrasound signal is used as the characteristic value to reflect the change of infrasound caused by pipeline leakage.

[0189] In adopting Figure 1The experimental device collects the leakage signal process, the leakage signal of different aperture uses EEMD analysis method to obtain IMF1 to IMF9 nine IMF components and residual IMF10, the correlation coefficient of each IMF component is calculated according to formula (6), and the results are shown in Figure 7 Table 1.

[0190] Table 1 is the data of correlation coefficient of each IMF component

[0191]

[0192] According to the correlation coefficient relationship between each IMF component and the original signal, the calculation results are analyzed as follows:

[0193] (1) The correlation coefficient of each component of IMF1 to IMF10 and the original signal shows a decreasing rule;

[0194] (2) Through statistics, the energy value of each component of IMF1 to IMF6 accounts for more than 90% of the total energy value of the original signal.

[0195] Based on the above two points, referring to formula (7), the root mean square value vector [E1, E2, E3, E4, E5, E6] of IMF1 to IMF6 is used as the input key feature value of APSO-HMM pipeline leakage model training and model judgment.

[0196] S4: Use adaptive particle swarm optimization (APSO) algorithm to optimize the initial value of each state observation value, that is, the initial value of B. The comparison chart before and after optimization is shown in Figure 8 .

[0197] S5: Fitness value refers to the value of objective function. The following is the difference comparison of model fitness convergence value and convergence speed between APSO-HMM pipeline leakage model and conventional HMM pipeline leakage model. From the results of Figures 9-11 , three different aperture leakage states can significantly improve the convergence speed of HMM by increasing the parameter optimization of adaptive particle swarm global optimization, and avoid falling into local minimum.

[0198] Table 2 is the output log likelihood probability and identification result of aperture of 1mm, 2mm and 3mm

[0199]

[0200] In Table 2, taking the running number 1 (the actual leakage aperture is 1 mm) as an example, the feature vector is input into each trained APSO-HMM model, and the maximum likelihood probability in the model of 1 mm leakage aperture is -19.1923 > -31.0048 > -65.2358, the identification effect is more obvious, which proves the feasibility of the application of the APSO-HMM to the pipeline leakage identification.

[0201] In order to compare with the conventional HMM pipeline leakage model, the same training samples as the APSO-HMM model are selected for leakage identification by using the HMM. The labels of the leakage aperture of 1, 2 and 3 mm under the pressure of 0.1 MPa are set as

[123] , wherein "1" represents the leakage aperture of 1 mm, "2" represents the leakage aperture of 2 mm, and "3" represents the leakage aperture of 3 mm. Correspondingly, the data of the pressure of 0.1 MPa and the leakage aperture of 1, 2 and 3 mm are tested in turn, and the result label is

[112] , that is, only the leakage aperture of 1 mm is identified, and the accuracy is only 33.33%, which is less than the identification result of the APSO-HMM.

[0202] Through multiple sample tests, the comparison of the identification results of the APSO-HMM model and the HMM is shown in Table 3. As shown in the table, the APSO-HMM model has good identification effect.

[0203] Table 3 Comparison of identification results of the APSO-HMM model and the HMM model

[0204]

[0205] Based on the above ideal embodiments according to the present application, through the above description, relevant personnel can make various changes and modifications without deviating from the technical idea of the present application. The technical scope of the present application is not limited to the contents in the specification, and must be determined according to the scope of the claims.

Claims

1. A method for identifying a pipeline leakage aperture based on APSO-HMM, comprising: constructing a plurality of training sets, each corresponding to a pipeline leakage aperture; wherein the construction of each training set comprises: obtaining a plurality of infrasonic original leakage signals x(t) corresponding to the pipeline leakage aperture; extracting a feature vector T of each infrasonic original leakage signal x(t), and taking the set of the feature vectors T of the plurality of infrasonic original leakage signals x(t) as the training set of the HMM model for the pipeline leakage aperture; inputting each training set into an HMM model for training iteration, and using an APSO algorithm to optimize the initial value of the observation value under each state, i.e. the initial value of B, during the training process; bringing the optimal initial value obtained by the APSO algorithm into the HMM model, and using a Baum-Welch algorithm to iteratively calculate the HMM model to obtain a trained APSO-HMM model; extracting a feature vector T of a to-be-identified infrasonic original leakage signal x(t), and inputting the feature vector T into each trained APSO-HMM model to output the pipeline leakage aperture corresponding to the APSO-HMM model with the maximum log-likelihood probability as the identification result; wherein the step of extracting the feature vector T of the infrasonic original leakage signal x(t) comprises: performing Gaussian low-pass filtering denoising processing on the infrasonic original leakage signal x(t) to obtain a denoised signal x'(t) ; decomposing the denoised signal x'(t) using an EEMD algorithm, and selecting effective IMF components from the IMF components obtained by the decomposition based on a correlation coefficient; and taking the energy entropy of the effective IMF components as elements to form the feature vector T of the infrasonic original leakage signal x(t). 2.The method of claim 1, wherein the EEMD algorithm is used to decompose the denoised signal x'(t), and the decomposition comprises: step 11, introducing Gaussian white noise into the denoised signal x'(t), and the calculation formula is: step 12, performing EMD decomposition on the signal after the white noise is added according to the calculation method in EMD to obtain a group of IMF components: step 13, returning to step 11 and repeating N times; and step 14, calculating the mean of each IMF component to obtain the final decomposition result, and the calculation formula is: 3.The method of claim 1, wherein the effective IMF components are selected from the IMF components obtained by the decomposition based on a correlation coefficient, and the selection comprises: calculating the correlation coefficient between each IMF component and the denoised signal x'(t) using a correlation coefficient formula; sorting each IMF component in descending order of the correlation coefficient, and taking the first J IMF components as the effective IMF components, wherein the sum of the energy values of the first J-1 IMF components is less than a preset percentage of the total energy value of the denoised signal x'(t), and the energy values of the first J IMF components are not less than a preset percentage of the total energy value of the denoised signal x'(t). ​ ​ ​ ​ ​ ​ e i (t) = x'(t) + n i (t) where e i (t) is the signal after adding the ith noise, n i (t) is the ith white noise signal; ​ where m ij (t) is the signal e i (t) is the jth IMF component, r ij (t) is the jth residue component resulting from the decomposition of the signal e ​ ​ In the formula, m j (t) is the jth IMF component obtained by EEMD decomposition, and r(t) is the residual component obtained by EEMD decomposition. ​ ​ ​ where cov(x, y) is the covariance of x, y; D(x) and D(y) are the variances of x, y, respectively; p xy The value range is [0, 1], and x and y respectively represent the IMF component to be calculated and the denoising signal x'(t); ​ 4. The APSO-HMM-based pipeline leak aperture identification method according to claim 1, characterized in that, The calculation formula of energy entropy is: where m j (t) is the jth effective IMF component, E j is the energy entropy of the jth effective IMF component.

5. The APSO-HMM-based pipeline leak aperture identification method according to claim 1, characterized in that, Further comprising, normalizing the feature vector T; wherein, In the training process, after normalizing each feature vector T in the training set, input the HMM model training; In the identification process, after normalizing the feature vector T of the infrasound wave original signal leak signal x(t) to be identified, input each trained APSO-HMM model; T = [E1, E2,..., E J ]; T' = [E1 / E, E2 / E,..., En / E] (1) J / E] where E1, E2,..., E J denotes the energy entropy of the 1st to Jth effective IMF components, and T' denotes the normalized vector.

6. The APSO-HMM-based pipeline leak aperture identification method according to claim 1, characterized in that, The initial value of the observation value under each state is optimized by using the APSO algorithm, that is, the initial value of B; specifically, Step 21, first define the population number as M, each particle searches and optimizes in D-dimensional space with a certain speed, the optimal position of the nth particle is defined as pbest n ; the optimal position of the population is defined as gbest, in each iteration process, each particle and the population find their optimal positions, and the speed and position are updated according to the formula, thereby completing an iteration process; The velocity update formula is: The position update formula is: wherein: is the velocity in d-dimension and position of the nth particle in the tth iteration process; is the velocity in d-dimension and position of the nth particle in the t-1th iteration process; is the optimal position of the nth particle in the t-1th iteration process; is the optimal position of the swarm in the t-1th iteration process; w is the inertia weight, c1, c2 are learning factors, r1, r2 ∈ [0, 1] are random numbers; the following variables are defined: In the formula: is the current global optimal value, is the last global optimal value, α t is the population evolution speed: the ratio of the current global optimal value to the last global optimal value, α t The smaller the value of w, the faster the population evolution speed, and w should be increased to improve the global search performance of the population, and vice versa. is the current total average fitness, β t is the population convergence speed factor at the tth iteration, representing the population convergence degree: the ratio of the global optimal value to the current population average fitness. The larger the value, the more likely the population is to fall into a local optimum, and w should be increased to make the particles diverge, so that the population is easy to jump out of the local area, and vice versa; γ t is the population premature convergence factor at the tth iteration, representing γ t The smaller the value, the more likely the population is to tend to premature convergence; Step 22, when w shows linear or nonlinear decline, the algorithm has higher execution efficiency, and the mapping relationship between w and γ is designed as where: w t is the inertia weight for the tth iteration; When the iteration number is less than the preset value T rec , that is, t≤T rec , w always keeps a larger value to make the particle search in a larger range, and the expansion method is as follows: When the preset iteration value is reached, the inertia weight factor should be reduced adaptively, that is, t > T rec when, where w min and w max are the upper and lower bounds of the inertia weight factor, q1, q2 are empirical constants, T max is the maximum number of iterations, and τ is a threshold value. Step 23, to facilitate analysis and presentation, the problem space is reduced to one dimension, defining pbest i = p i , gbest = g i , then the velocity and position update equations are Through the iterative relationship, we get wherein: is the velocity in d-dimension and position of the i-th particle during the t-th iteration; is the velocity in d-dimension and position of the i-th particle during the t+1-th iteration; is the velocity in d-dimension and position of the i-th particle during the t+2-th iteration; In combination, the second-order difference equation for particle velocity is obtained: By analyzing the stability of the equation, we can get the following adaptive factor: In the formula, are the values of the learning factors c1, c2 in the tth iteration. Step 24, after the improvement of the weight factor in step 22 and the learning factor in step 23, the optimal solution B is obtained.

7. The APSO-HMM-based pipeline leak aperture identification method according to claim 1, characterized in that, The step of obtaining the log-likelihood probability by the model is: Step 41, local state initialization: δ t (i) = π i b i (O1),i = 1,2,3,...N where δ t (i) is the maximum probability of all transition paths of hidden state i at time t; wherein is the previous node hidden state of the maximum probability among all transition paths of hidden state i at time t. Step 42, local state of dynamic programming recursion time t (t = 1, 2, 3,..., T); Step 43, solve for T time delta T (i) probability of maximum occurrence, and solve for T time maximum delta T (i) P * = max 1≤j≤N δ T (i) Step 44, by starting backtracking, i.e. The hidden state sequence with the maximum probability is obtained The recursive expression of δ is obtained Estimate The recursive expression for Finally, the maximum probability value P* of the last time T is calculated, which is the log-likelihood probability.

Images