GPR image clutter removal method based on robust non-negative matrix factorization

Abstract

The invention relates to a GPR clutter removal method based on robust non-negative matrix factorization (RNMF), which comprises the following steps of: representing a GPR image by a rectangular matrix X with the size of M * N, and decomposing a data matrix into a sparse error matrix S and the sum of products of two non-negative matrixes W and H by using the RNMF; iteratively updating W, H and S according to the iteration updating equation; and the RNMF is compared with the most advanced clutter removal methods such as MCA and RPCA in addition to the traditional PCA, and an analog data set is compared with a real data set. According to the method, the visualization effect can be improved to a great extent, a clearer background can be displayed, no prior information is needed, and the method has practical application significance.

Description

Technical Field

[0001] The invention relates to the technical field of GPR image processing, and in particular to a GPR image clutter removal method based on robust non-negative matrix decomposition. Background Art

[0002] Clutter caused by ground reflection, coupling between transmitting and receiving antennas, and reflection from underground inhomogeneities seriously deteriorates the performance of underground target detection algorithms. Therefore, clutter removal is an active research direction in the field of Ground Penetrating Radar (GPR). In addition to traditional methods based on subspace technology, such as Singular Value Decomposition (SVD), Principal Component Analysis (PCA) or Independent Component Analysis (ICA), many methods based on multi-resolution analysis or matrix decomposition have been proposed in recent years. Subspace technology decomposes the GPR image into several components, the most important of which corresponds to clutter, followed by the target, and the remaining part corresponds to noise. However, target information can be contained in several components; similarly, multiple targets may not be represented by a single component. Another difficulty of ICA is to determine the components that represent target or clutter information, because ICA does not provide components arranged in order.

[0003] Bao et al. proposed to decompose the original GPR image into sub-images containing both multi-scale information and multi-directional information. The sub-image containing clutter information was removed, and then the remaining sub-image was inversely transformed to obtain a clutter-free image. In the edge-preserving filter for the multi-scale part, better performance was shown in target detection. However, these methods provide diffuse target components in the generated images, which reduces the visual quality.

[0004] Another approach is based on Morphological Component Analysis (MCA), which decomposes GPR images into target and clutter components, each of which is sparsely modeled using an appropriate dictionary learning method. Although this method has significant performance, it suffers from high complexity and inappropriate dictionary selection.

[0005] In recent years, low-rank sparse representation, which represents a given matrix as the sum of a low-rank matrix and a sparse matrix, has found many applications in image processing, such as denoising, restoration, classification, change detection, and background subtraction.

[0006] Low Rank and Sparse Matrix Decomposition is a matrix decomposition technique that combines the ideas of low rank approximation and sparse representation. It decomposes a matrix into the product of a low rank matrix and a sparse matrix to better capture the potential structure and outliers in the data. This technique is widely used in data dimensionality reduction, signal processing, image processing and other fields. Non-negative Matrix Factorization (NMF) is a matrix decomposition technique that decomposes a non-negative matrix into the product of two or more non-negative matrices. This technique has been widely used in data mining, feature extraction, image processing and other fields.

[0007] In the context of ground penetrating radar (GPR), the target component can be considered sparse and the clutter part can be represented by a low-rank matrix. Nowadays, robust principal component analysis (RPCA) is a widely used tool that can perform low-rank and sparse decomposition for clutter removal in underground and through-wall imaging applications. Summary of the invention

[0008] The present invention provides a GPR (Ground Penetrating Radar, GPR) image clutter removal method based on Robust Nonnegative Matrix Factorization (RNMF). Different from the method of obtaining low rank using only Nonnegative Matrix Factorization (NMF), the method proposed in the present invention decomposes the GPR image into low rank and sparse components at the same time, representing clutter and target respectively. The method can greatly improve the visualization effect, show a clearer background, and does not require any prior information, which is more practical.

[0009] The technical solution adopted by the present invention is:

[0010] A GPR image clutter removal method based on robust non-negative matrix factorization is proposed. The GPR image is represented as a data matrix of size M×N. The data matrix is ​​decomposed into a sparse error matrix S and the sum of the products of two non-negative matrices W and H through robust non-negative matrix factorization (RNMF). W, H, and S are iteratively updated according to the iterative update equation to remove GPR image clutter.

[0011] The GPR image is represented by a data matrix X of size M×N, where M and N are the depth index and the number of antenna positions, respectively;

[0012] X=X 目标 +X 杂波 (1);

[0013] In formula (1), X represents the matrix representation of the input GPR image; X 目标 represents the matrix representation of the target image; X 杂波 Represents the matrix expression of the clutter component.

[0014] Robust Non-negative Matrix Factorization (RNMF) decomposes the data matrix X into a sparse error matrix S and the sum of the products of two non-negative matrices W and H:

[0015] X≈WH+S (6);

[0016] Where X∈R M×N , W∈R M×k , H∈R k×N , S∈R M×N , k is the rank of NMF decomposition; R represents the correlation matrix.

[0017] For robust non-negative matrix factorization (RNMF), the optimization problem can be expressed as:

[0018]

[0019] In the formula: min means to find the minimum value; ||·|| F is the F-norm of the matrix; λ>0, λ is a regularization parameter, which is used to control the complexity of the model, prevent overfitting, and possibly improve the generalization ability of the model, and plays a vital role in GPR image decomposition. m represents the total number of rows in the matrix S, and n represents the total number of columns in the matrix S.

[0020] In order to find the appropriate W and H, the iterative update equation is as follows:

[0021]

[0022] In formula (8a) and formula (8b), W and H are called basis matrix and coefficient matrix respectively, both of which are non-negative matrices; S represents the sparse error matrix; |·| represents the determinant; W T , H T They are represented as the transposed matrices of the basis matrix W and the coefficient matrix H respectively; ij represents the basis matrix after W iterations; H ij represents the coefficient matrix after H iteration; ← represents the W ij and H ij Assignment.

[0023] W and H can be normalized as follows:

[0024]

[0025] In formula (9a) and formula (9b): W ki represents the kth row and the ith column; W kj represents the kth row and jth column; n represents the matrix order.

[0026] In order to optimize the sparsity parameter, S, W and H need to be fixed. The update rule for S is as follows:

[0027]

[0028] In formula (10): If Then all elements in S are equal to zero.

[0029] Iterate W, H, and S according to equations (8a) to (10).

[0030] For a given predetermined threshold or maximum number of iterations, the iterative process can be terminated; the resulting S matrix corresponds to the target portion of the GPR image.

[0031] The regularization parameter λ controls the trade-off between the low-rank component L, corresponding to the RNMF, and the sparse component S; the λ value is chosen to be the Peak Signal to Noise Ratio (PSNR) value, and then it is fixed for all experiments. The optimal λ value is given by the peak of the PSNR graph. When the PSNR value is lower than this value, the PSNR value drops sharply because more false information is passed to the S image, and the PSNR value drops gently.

[0032] The present invention provides a GPR image clutter removal method based on robust non-negative matrix decomposition, which has the advantages of:

[0033] 1) Decompose the GPR image into low-rank and sparse components simultaneously, representing the clutter component and the target component respectively. This method is faster and has better results.

[0034] 2) It can greatly improve the visualization effect, show a clearer background, and does not require any prior information, which has great application value. BRIEF DESCRIPTION OF THE DRAWINGS

[0035] The present invention will be further described below in conjunction with the accompanying drawings and examples;

[0036] Figure 1 Schematic diagram of the change of λ value with PSNR.

[0037] Figure 2 Schematic diagram of the experimental setup for the simulated dataset.

[0038] Figure 3 The results of RNMF clutter removal at different levels.

[0039] Figure 4 The simulation data results.

[0040] Figure 5 Experimental settings for real datasets.

[0041] Figure 6 The results are from the real dataset. DETAILED DESCRIPTION

[0042] The GPR clutter removal method based on robust non-negative matrix decomposition includes the following steps:

[0043] Step 1: The GPR image is represented by a rectangular matrix X of size M×N. Robust Non-negative Matrix Factorization (RNMF) decomposes the data matrix into a sparse error matrix S and the sum of the products of two non-negative matrices W and H.

[0044] Step 2: Iteratively update W, H, and S according to the iterative update equation;

[0045] Step 3: In addition to RNMF and traditional PCA, the most advanced clutter removal methods such as MCA and RPCA are compared on simulated and real data sets.

[0046] The raw GPR data is represented as the sum of low-rank and sparse matrices, corresponding to the clutter and target components, respectively. The low-rank and sparse decomposition uses a robust non-negative matrix factorization RNMF. Although similar to the robust principle component analysis (RPCA) widely used in image processing applications and GPR recently, the proposed method is faster and more effective.

[0047] Non-negative Matrix Factorization (NMF) usually uses an iterative optimization algorithm to solve W and H. Its specific steps include:

[0048] (1) Initialization: Randomly initialize the matrices W and H, or use some heuristic method to initialize them;

[0049] (2) Iterative update;

[0050] (3) Termination condition: Stop when the maximum number of iterations is reached, the error is lower than a certain threshold or other stopping criteria.

[0051] Non-negative matrix factorization provides an effective method for data dimensionality reduction and feature extraction, which is applicable to various types of data analysis tasks, such as image processing, text mining, bioinformatics, etc. The advantage of non-negative matrix factorization NMF is that it forces the matrix elements to be non-negative, making the decomposition physically or biologically interpretable, and can also achieve data compression and dimensionality reduction through low-rank approximation.

[0052] The GPR image is represented by a matrix X of size M×N, where M and N are the depth index and number of antenna positions (downward), respectively. It is given as the sum of the target and clutter components:

[0053] X=X 目标 +X 杂波 (1);

[0054] In formula (1), X represents the matrix representation of the input GPR image, X 目标 Represents the matrix representation of the target image, X 杂波 Represents the matrix expression of the clutter component.

[0055] Using non-negative matrix factorization (NMF), the matrix X of the GPR image is decomposed into the product of two matrices W and H, where both matrices are constrained to be non-negative and are expressed as:

[0056] X≈WH stW,H≥0 (2);

[0057] In formula (2), W and H are the basis matrix and coefficient matrix respectively, X∈R M×N , R M×N represents the M×N matrix R, k is the rank of the matrix. Generally, the value of k is k <min(M,N)。

[0058] represents an M×k positive definite matrix, represents a k×N positive definite matrix.

[0059] Formula (2) can be transformed into an optimization problem:

[0060]

[0061] In formula (3), D represents the divergence function, which is used to find a convenient W and H so that the divergence between the data matrix X and the non-negative matrices W and H is minimized; D(X; W, H) represents the divergence between the target matrix X and W and H.

[0062] Using the Kulbeck-Leibler divergence, the iterative update rule can be derived as:

[0063]

[0064] In Equation (4): 1 represents the identity matrix; W and H are respectively called the basis matrix and the coefficient matrix; W T , H T represent the transpose matrices of the basis matrix W and the coefficient matrix H respectively; ← represents further assigning a value to W.

[0065] When k < N, a low-rank approximation of the data matrix X is obtained by choosing k:

[0066] Since clutter dominates the target signal or any other discontinuous part in the soil, the most significant decomposition component provides the clutter when k = 1;

[0067] X 杂波 = W M×1 H 1×N (5);

[0068] In Equation (5), X 杂波 represents the matrix expression form of the clutter component; W M×1 represents an M×1 basis matrix; H 1×N represents a 1×N coefficient matrix. X 目标 is directly obtained from Equation (1).

[0069] Robust Non-negative Matrix Factorization (RNMF) is a variant of NMF, aiming to improve the robustness to noise and outliers in data. In practical applications, data is often disturbed by various factors, such as measurement errors, outliers, or data missing, etc. These factors may affect the performance of traditional NMF algorithms. RNMF can effectively cope with these interferences and improve the robustness of the model by introducing robustness constraints or improving the loss function, etc. RNMF has a wide range of applications in various fields, especially useful when dealing with data with noise and outliers in the real world.

[0070] In Robust Non-negative Matrix Factorization (RNMF), assuming that the data matrix X is sparsely corrupted, RNMF decomposes the data matrix X into the sum of a sparse error matrix S and the product of two non-negative matrices W and H:

[0071] X ≈ WH + S (6);

[0072] where, X ∈ R M×N , W ∈ R M×k , H ∈ R k×N , S ∈ R M×N, k is the rank of NMF decomposition; R represents the matrix.

[0073] For robust non-negative matrix factorization (RNMF), the optimization problem can be expressed as:

[0074]

[0075] In the formula: min means to find the minimum value; ||·|| F is the F-norm of the matrix; λ>0, λ is a regularization parameter, which is used to control the complexity of the model, prevent overfitting, and possibly improve the generalization ability of the model, and plays a vital role in GPR image decomposition. m represents the total number of rows in the matrix S, and n represents the total number of columns in the matrix S.

[0076] In order to find the appropriate W and H, the iterative update equation is as follows:

[0077]

[0078] In formula (8a) and formula (8b), W and H are called basis matrix and coefficient matrix respectively, both of which are non-negative matrices; S represents the sparse error matrix; |·| represents the determinant; W T , H T They are represented as the transposed matrices of the basis matrix W and the coefficient matrix H respectively; ij represents the basis matrix after W iterations; H ij represents the coefficient matrix after H iteration; ← represents the W ij and H ij Assignment.

[0079] W and H can be normalized as follows:

[0080]

[0081] In formula (9a) and formula (9b): W ki represents the kth row and the ith column; W kj represents the kth row and jth column; n represents the matrix order.

[0082] In order to optimize the sparsity parameter, S, W and H need to be fixed. The update rule for S is as follows:

[0083]

[0084] In formula (10): If Then all elements in S are equal to zero. For large regularization parameters, RNMF is the same as the traditional NMF method.

[0085] The GPR image X can be represented as the sum of a low-rank matrix (clutter component), while the target component can be considered sparse relative to the entire GPR image. Therefore, similar to RPCA, RNMF can be used to decompose the GPR image into its low-rank: k < < min(M, N) and sparse components, as described in Equation (1), where X 目标 =WH,X 杂波 =S.

[0086] W, H, S are iterated according to equations (8a) to (10). The iteration process can be terminated for a given predetermined threshold or maximum number of iterations. The resulting S matrix corresponds to the target portion of the GPR image.

[0087] The choice of regularization parameter λ is crucial for RNMF. As known from RPCA applications, λ controls the trade-off between the low-rank component L (corresponding to WH of RNMF) and the sparse component S. As the value of λ decreases, more sparse component information (targets) are transferred from the original GPR image to S. On the other hand, if the value of λ is below a certain value, clutter and noise information will also be transferred along with the target information, making the target indistinguishable. Although the regularization parameter is suggested to be λ in many RPCA applications, But it can vary depending on the application. In the method of the present invention, for the sample GPR images selected from the dataset of the present invention, the λ value is selected as the Peak Signal to Noise Ratio (PSNR) value, and then it is fixed for all experiments. The variation of λ value with PSNR can be found in Figure 1 The optimal λ value is given by the peak of the PSNR graph. When the PSNR value is lower than this value, the PSNR value drops sharply because more false information is transferred to the S image, and the PSNR value drops gently. Table 1 shows the clutter removal algorithm based on RNM in GPR.

[0088] Table 1 Clutter removal algorithm based on RNM in GPR

[0089]

[0090] The proposed RNMF method is compared with the recently introduced MCA, RPCA and NMF as well as the traditional PCA method on both simulated and real datasets. For the convenience of comparison, both visualization and quantitative results are given. The quantitative results are based on PSNR.

[0091] (1) Simulation data set results:

[0092] The simulation was built using gprMax. The antenna type was selected from its built-in library as a geophysical survey system, operating at a frequency of 1.5 GHz. Six different soil types were used in the simulation: dry sand, wet sand, wet sand, dry clay, wet clay, and dry loam; and two different target materials: plastic and aluminum. The experimental setup of the scenario is shown in Figure 2 shown.

[0093] PSNR is given as follows:

[0094]

[0095] In formula (11), PSNR is the peak signal-to-noise ratio, the unit is decibel (dB), log() represents the logarithmic function, and MSE is the mean square error.

[0096] The mean square error (MSE) can be expressed as:

[0097]

[0098] In formula (12), M is the number of rows of the input matrix X, and N is the number of columns of the input matrix X.

[0099] Here, the λ value of RNMF is 3×10 -4 , the λ value of RPCA is 1.8×10 -2 Rank k is another important parameter of RNMF. The effects of different k values ​​are as follows: Figure 3 As shown in Figure 1, the original data and the results obtained when k = 1, k = 2, and k = 5 are given. Since the clutter is much stronger than the reflected signal from the target, k = 1 is selected to directly give the target, while the full rank gives the original GPR image. It can be observed that as the k value increases, the target image deteriorates due to the part-based nature of NMF. The number of iterations is selected as 1×10 4 .

[0100] Figure 4 Sub-images (a) and (b) show the original GPR image and reference image of a simulated scenario where the target is aluminum, the burial environment is dry sand, and the burial depth is 1 cm. Figure 4 Sub-images (c) to (g) show the PCA, NMF, MCA, RPCA, RNMF clutter and RNMF target images respectively. All methods, except PCA, which presents a vertical line at the target position, can obtain the target image. NMF models the clutter as a low-rank representation, and the target is obtained as the difference between the original GPR image and the low-rank representation.

[0101] Since the classic NMF decomposition does not provide any sparse components like RNMF, the resulting target image contains some clutter residues, which are observed as horizontal lines at the depth of the target. Looking closely at the results of MCA, some distortions appear above the target due to the inappropriate dictionary. RPCA and RNMF present similar visualization results.

[0102] Table 2 Peak signal-to-noise ratio results of aluminum and plastic at different burial depths

[0103]

[0104]

[0105] Table 3 Peak signal-to-noise ratio (PSNR) results of aluminum and plastic targets in different soils

[0106]

[0107] Table 4 Running time of clutter removal methods in GPR

[0108]

[0109] Tables 2 and 3 give the quantitative results for different target materials, burial depths and soil types. Table 3 shows the PSNR values ​​of aluminum targets buried in dry sand and plastic targets buried in dry clay at different burial depths using PCA, NMF, MCA, RPCA and the RNMF method proposed in the present invention. For aluminum targets at all burial depths, RNMF achieved the best PSNR value, followed by RPCA and MCA. For plastic targets, the performance of MCA, RPCA and the RNMF proposed in the present invention varied with the burial depth.

[0110] Although the proposed method still lags behind MCA at burial depths of 0 and 1 cm, MCA has higher complexity and requires a suitable dictionary. Only at a burial depth of 4 cm does the proposed method outperform RPCA, followed by MCA. For burial depths of 2 and 10 cm, the proposed method has the best score. The PSNR value of NMF, which only contains the low-rank part, is still much lower than the low-rank sparse decomposition methods RPCA and RNMF. Among all the methods, PCA shows the worst results.

[0111] Table 3 gives the PSNR results of aluminum and plastic targets at 2 cm depth in different soil types. Again, RNMF showed the best results, while MCA had similar results to RPCA, followed by NMF and PCA, except for the plastic target buried in moist sand, where MCA showed approximately the same performance as the proposed method.

[0112] (2) Results of real datasets

[0113] In order to verify the effectiveness of the algorithm of the present invention, real GPR data is used for verification. The experimental design is as follows: Figure 5 As shown, the original data is Figure 6 The upper and lower rows correspond to the dry and wet conditions, respectively. In these scenarios, PMA-3, PMA-1, rock and copper targets are buried at a depth of 5 cm, where PMA-3 and PMA-1 refer to plastic-cased anti-personnel mines. The soil type is moist clay with small rocks. The surface is uneven with a maximum difference of 10 cm between the highest and lowest points. The antenna is placed 5 cm above the highest point on the surface. The GPR settings for this scan are: averaging 10 times per A-Scan, 25 ps sampling interval, 512 samples, 1 GHz antenna, and pulse width 1 ns.

[0114] The scanning range was dx=50 cm×dy=196 cm, and the scanning step length in each direction was 1 cm. After leaving the untouched object for 23 days, the same scan was repeated to obtain the raw data of the drying condition.

[0115] Figure 6 Sub-figures (b) to (f) of the figure show the target images obtained by PCA, NMF, MCA, RPCA, and RNMF in the dry condition. MCA, RPCA, and RNMF outperform PCA and NMF, where RNMF gives the clearest background. In PCA and NMF, the background contains a lot of distortion. The results in the real dataset also match those in the simulated dataset. Figure 6 In the wet cases of sub-images (b) to (f), the performance of the compared methods has deteriorated to the expected extent, especially the residual clutter between targets is more obvious. However, the results of RNMF seem to be more robust to the moisture effect.

[0116] The running time is shown in Table 4. The scenario is as follows Figure 6 As shown in the sub-figure (a) of . PCA and NMF have the least computation time, while MCA has the highest computation time. However, the quantitative and visual results provided by PCA and NMF are not satisfactory. Compared with RPCA, the proposed method RNMF not only outperforms RPCA in both visual and quantitative results, but also has a speed 25 times that of RPCA under the selected λ value. Therefore, in addition to having high visual performance, it is also suitable for working in real time like PCA and NMF.

[0117] This paper proposes a new clutter removal method based on low-rank sparse decomposition of GPR images, called RNMF. It not only provides a low-rank decomposition of the input data like the traditional one, but also a sparse component, which is used to decompose the original GPR image. In the obtained decomposition, the low-rank part corresponds to the clutter, while the sparse part directly obtains the clutter-free target image. Visualization and quantitative results on simulated and real datasets verify the superiority of the proposed RNMF method. RNMF provides a clearer background than MCA, it does not require any prior information such as a proper dictionary, and is as fast as traditional clutter removal methods.

Claims

1. A GPR image clutter removal method based on robust non-negative matrix decomposition, characterized by: The GPR image is represented as a data matrix of size M×N. The data matrix is ​​decomposed into a sparse error matrix S and the sum of the products of two non-negative matrices W and H through robust non-negative matrix factorization (RNMF). W, H, and S are iteratively updated according to the iterative update equation to remove GPR image clutter.

2. The GPR image clutter removal method based on robust non-negative matrix decomposition according to claim 1, characterized in that: The GPR image is represented by a data matrix X of size M×N, where M and N are the depth index and the number of antenna positions, respectively; X=X 目标 +X 杂波 (1); In formula (1), X represents the matrix representation of the input GPR image; X 目标 represents the matrix representation of the target image; X 杂波 Matrix expression representing the clutter component; Robust non-negative matrix factorization (RNMF) decomposes the data matrix X into a sparse error matrix S and the sum of the products of two non-negative matrices W and H: X≈WH+S (6); Where X∈R M×N , W∈R M×k , H∈R k×N , S∈R M×N , k is the rank of NMF decomposition; R represents the correlation matrix; For robust non-negative matrix factorization (RNMF), the optimization problem can be expressed as: In the formula: min means to find the minimum value; ||i|| F is the F-norm of the matrix; λ>0, λ is a regularization parameter used to control the complexity of the model; m represents the total number of rows in the matrix S, and n represents the total number of columns in the matrix S; In order to find the appropriate W and H, the iterative update equation is as follows: In formula (8a) and formula (8b), W and H are called basis matrix and coefficient matrix respectively, both of which are non-negative matrices; S represents the sparse error matrix; |i| represents the determinant; W T , H T They are represented as the transposed matrices of the basis matrix W and the coefficient matrix H respectively; ij represents the basis matrix after W iterations; H ij represents the coefficient matrix after H iteration; ← represents the W ij and H ij Assignment; W and H can be normalized as follows: In formula (9a) and formula (9b): W ki represents the kth row and the ith column; W kj represents the kth row and jth column; n represents the matrix order; In order to optimize the sparsity parameter, S, W and H need to be fixed; the update rule for S is as follows: In formula (10): If Then all elements in S are equal to zero; Iterate W, H, and S according to equations (8a) to (10).

3. The GPR image clutter removal method based on robust non-negative matrix decomposition according to claim 2, characterized in that: For a given predetermined threshold or maximum number of iterations, the iterative process can be terminated; the resulting S matrix corresponds to the target portion of the GPR image.

4. The GPR image clutter removal method based on robust non-negative matrix decomposition according to claim 2, characterized in that: The regularization parameter λ controls the trade-off between the low-rank component L, corresponding to the RNMF, and the sparse component S; the λ value is chosen to be the peak signal-to-noise ratio (PSNR) value, and then it is fixed for all experiments; the optimal λ value is given by the peak value of the PSNR; when the PSNR value is lower than this value, the PSNR value drops sharply because more false information is passed into the S image, and the PSNR value drops gently.

Images