Hyperspectral anomaly detection method based on low-rank and sparse prior constrained autoencoder

Abstract

The present invention discloses a hyperspectral anomaly detection method based on a low-rank and sparse prior constrained autoencoder, which mainly solves the problem that the existing detection scheme has poor detection effect in complex application scenarios. It includes: 1) constructing an autoencoder with two shared encoders; 2) establishing an objective function to decompose the hyperspectral image into low-rank components, sparse components and residual components; 3) using low-rank components and sparse components as labels to calculate the reconstruction losses of the two autoencoders respectively; 4) integrating hyperspectral image decomposition and autoencoder training into a unified framework, and jointly optimizing in an end-to-end manner; 5) using the RX detector to process and fuse the reconstruction results to obtain the final detection results. The present invention combines the advantages of a linear low-rank and sparse model with a nonlinear autoencoder, and performs joint optimization in an end-to-end manner to avoid the model from falling into a local optimum, thereby effectively improving the performance of hyperspectral anomaly detection.

Description

Technical Field

[0001] The present invention belongs to the technical field of hyperspectral remote sensing image processing, and further relates to hyperspectral anomaly detection, specifically a hyperspectral anomaly detection method based on low-rank and sparse prior constrained autoencoder, which can be used for mineral exploration, fire detection, marine vessel detection, etc. Background Art

[0002] Hyperspectral anomaly detection aims to distinguish background from anomalies without prior knowledge. Since it does not require any prior knowledge, hyperspectral anomaly detection technology is widely used in mineral exploration, fire detection, marine ship detection and other fields.

[0003] In recent years, in the field of hyperspectral anomaly detection, linear low-rank and sparse models and nonlinear autoencoders have developed rapidly and have attracted widespread attention from researchers. For example, Zhang YuXiang et al. proposed a hyperspectral anomaly detection method based on low-rank and sparse matrix decomposition in the paper "A Low-Rank and Sparse Matrix Decomposition-Based Mahalanobis Distance Method for Hyperspectral Anomaly Detection". Xiang Pei et al. proposed a hyperspectral anomaly detection method based on guided autoencoders in the paper "Hyperspectral Anomaly Detection with Guided Autoencoder", which aims to reduce the feature representation of abnormal targets. In the above methods, the linear low-rank and sparse model is interpretable, but it cannot characterize complex scenes; on the contrary, the nonlinear autoencoder can extract the distinguishing features of background and anomalies in complex scenes, but it is not interpretable. Therefore, how to combine the advantages of the two to improve the accuracy of hyperspectral anomaly detection is a problem worth exploring. Summary of the invention

[0004] The purpose of the present invention is to address the deficiencies of the above-mentioned prior art and propose a hyperspectral anomaly detection method based on a low-rank and sparse prior constrained autoencoder, which solves the problem that the existing detection scheme has poor detection effect in complex application scenarios. The present invention combines the advantages of a nonlinear autoencoder with a linear low-rank and sparse model, and optimizes the model in an end-to-end manner, thereby avoiding the model from falling into a local optimum, thereby effectively improving the accuracy of hyperspectral anomaly detection.

[0005] The idea of ​​implementing the present invention is: first, construct two autoencoders, where the two autoencoders share an encoder, and establish an objective function to decompose the hyperspectral image into low-rank components, sparse components and residual components, and use the obtained low-rank components and sparse components as labels to calculate the reconstruction losses of the two autoencoders respectively; secondly, integrate the hyperspectral image decomposition and the training of the autoencoder into a unified framework to achieve end-to-end joint optimization; finally, send the hyperspectral image to the trained autoencoder to obtain the reconstructed result, and use the RX detector to process the reconstructed result and perform nonlinear fusion to obtain the final detection result.

[0006] The specific steps of the present invention to achieve the above-mentioned purpose include the following:

[0007] (1) Constructing an autoencoder with low rank and sparse prior constraints:

[0008] The first autoencoder LR' is composed of the encoder_E and the low-rank decoder_LR, and the second autoencoder S' is composed of the encoder_E and the sparse decoder_S, and the first autoencoder LR' and the second autoencoder S' share the encoder_E, so as to obtain the first autoencoder LR' and the second autoencoder S' with low-rank and sparse prior constraints;

[0009] (2) End-to-end joint optimization:

[0010] (2.1) Decompose the hyperspectral image X into a low-rank component L, a sparse component S and a residual component E;

[0011] (2.2) Using the low-rank component L and the sparse component S as labels, we establish the reconstruction loss l of the first autoencoder LR' and the second autoencoder S' respectively. LR and l S , and train them;

[0012] (2.3) The decomposition of the hyperspectral image X and the training of the autoencoder are integrated into a unified framework to achieve joint optimization in an end-to-end manner, and the following joint optimization objective function is obtained:

[0013]

[0014] stX=L+S+E

[0015] Among them, λ, β, γ, α and ν represent the coefficients of weighing each regularization term, which are the first weighing coefficient, the second weighing coefficient, the third weighing coefficient, the fourth weighing coefficient and the fifth weighing coefficient respectively; ||·||*, ||·|| 1 and ||·|| 2,1 They represent the nuclear norm, L 1 Norm and L 2,1 norm; st represents the constraint condition;

[0016] (2.4) The joint optimization objective function is solved using the augmented Lagrangian function, and back propagation and Adam optimizer are used for optimization to obtain the optimal parameters Θ of encoder_E, low-rank decoder_LR and sparse decoder_S en '、 and Get the trained first autoencoder LR' and second autoencoder S';

[0017] (3) Hyperspectral anomaly detection:

[0018] (3.1) The hyperspectral image X is fed into the trained first autoencoder LR' and the second autoencoder S' to obtain the reconstructed low-rank part L re and the reconstructed sparse part S re ;

[0019] (3.2) Calculate residual data using RX detector The test results of LR And reconstruct the sparse part S re The test results of S , and the test result d is calculated according to the following formula LR and d S Perform nonlinear fusion to obtain the final detection result d f :

[0020]

[0021] Among them, d f represents the final detection result, and k represents the fusion coefficient.

[0022] Compared with the prior art, the present invention has the following advantages:

[0023] First, the present invention effectively combines the advantages of linear-based low-rank and sparse models with those of nonlinear-based autoencoders, which can not only effectively characterize complex scenes, but also extract discriminative features of background and anomalies in complex scenes, and has strong interpretability.

[0024] Second, since the present invention integrates the decomposition of hyperspectral images and the training of autoencoders into a unified framework and performs optimization in an end-to-end manner, it avoids the model from falling into the local optimum, thereby effectively improving the hyperspectral anomaly detection performance. BRIEF DESCRIPTION OF THE DRAWINGS

[0025] Figure 1 Flow chart for realizing the method of the present invention;

[0026] Figure 2 It is a schematic diagram of the application process of the method of the present invention in a specific example;

[0027] Figure 3 Schematic diagram of the data set used in the present invention; wherein (I) represents the pseudo-color image and its label of the Texas Coast data set; (II) represents the pseudo-color image and its label of the Los Angeles data set; (III) represents the pseudo-color image and its label of the Salinas data set;

[0028] Figure 4 : It is a comparison chart of the detection results of the method of the present invention and the mainstream method on three data sets; wherein (Ia) represents the label of the Texas Coast data set; (Ib) to (Ii) respectively represent the detection images of RX, CRD, LSMAD, LSDM-MoG, RGAE, GAED, NJCR and the method of the present invention on the Texas Coast data set; (II-a) represents the label of the Los Angeles data set; (II-b) to (II-i) respectively represent the detection images of RX, CRD, LSMAD, LSDM-MoG, RGAE, GAED, NJCR and the method of the present invention on the Los Angeles data set; (III-a) represents the label of the Salinas data set; (III-b) to (III-i) respectively represent the detection images of RX, CRD, LSMAD, LSDM-MoG, RGAE, GAED, NJCR and the method of the present invention on the Salinas data set;

[0029] Figure 5 The results of the proposed method and the mainstream method on three data sets are shown in Figure 2. Curve graph; (I), (II), (III) respectively represent the results of the method of the present invention and the mainstream method on the Texas Coast, Los Angeles and Salinas data sets Curve, where P D represents the true positive rate, P F represents the false positive rate;

[0030] Figure 6 The results of the proposed method and the mainstream method on three data sets are shown in Figure 2. Curve graph; (I), (II), (III) respectively represent the results of the method of the present invention and the mainstream method on the Texas Coast, Los Angeles and Salinas data sets Curve, where P F represents the false positive rate, τ represents the threshold; DETAILED DESCRIPTION

[0031] The present invention will be further described in detail below in conjunction with the accompanying drawings and specific embodiments.

[0032] Example 1: Reference Figure 1 and Figure 2 The present invention proposes a hyperspectral anomaly detection method based on a low-rank and sparse prior constrained autoencoder, which is implemented in the following steps:

[0033] Step 1. Build an autoencoder with low-rank and sparse prior constraints:

[0034] The first autoencoder LR' is composed of encoder_E and low-rank decoder_LR, the second autoencoder S' is composed of encoder_E and sparse decoder_S, and the first autoencoder LR' and the second autoencoder S' share encoder_E, thereby obtaining the first autoencoder LR' and the second autoencoder S' with low-rank and sparse prior constraints. The number of nodes of encoder_E, low-rank decoder_LR and sparse decoder_S in this embodiment are {d, round(d / 2), round(d / 4), n}, {n, round(d / 4), round(d / 2), d} and {n, round(d / 4), round(d / 2), d}, respectively, where round represents a rounding operation, d represents the number of hyperspectral image bands, and n represents the number of hidden layer nodes.

[0035] Step 2. End-to-end joint optimization:

[0036] (2.1) The hyperspectral image X is decomposed into a low-rank component L, a sparse component S, and a residual component E. Specifically, the decomposition is performed according to the following objective function:

[0037]

[0038] (2.2) Using the low-rank component L and the sparse component S as labels, we establish the reconstruction loss l of the first autoencoder LR' and the second autoencoder S' respectively. LR and l S , and train it; the specific implementation is as follows:

[0039]

[0040]

[0041] Where b represents the batch size, and x i represent the low-rank component L, the sparse component S and the i-th pixel in the hyperspectral image X, respectively, and f en (·;·), and denote the encoding operation of encoder_E, the decoding operation of low-rank decoder_LR and the decoding operation of sparse decoder_S, Θ en , and denote the parameters of encoder_E, low-rank decoder_LR and sparse decoder_S respectively, Indicates L 2 The square of the norm;

[0042] During the training process, the batch size b is set to the total number of pixels N of the hyperspectral image X, and the loss functions of the first autoencoder LR' and the second autoencoder S' are transformed into:

[0043]

[0044]

[0045] in, Represents the square of the Frobenius norm.

[0046] (2.3) The decomposition of the hyperspectral image X and the training of the autoencoder are integrated into a unified framework to achieve joint optimization in an end-to-end manner, and the following joint optimization objective function is obtained:

[0047]

[0048] stX=L+S+E

[0049] Among them, λ, β, γ, α and ν represent the coefficients of weighing each regularization term, which are the first weighing coefficient, the second weighing coefficient, the third weighing coefficient, the fourth weighing coefficient and the fifth weighing coefficient respectively; ‖·‖ * 、||·|| 1 and ||·|| 2,1 They represent the nuclear norm, L 1 Norm and L 2,1 norm; st represents the constraint condition;

[0050] The loss function l in the first autoencoder LR' and the second autoencoder S' LR and l S Ignore coefficient The joint optimization objective function is converted into the following form:

[0051]

[0052] stX=L+S+E.

[0053] (2.4) The joint optimization objective function is solved using the augmented Lagrangian function, and back propagation and Adam optimizer are used for optimization to obtain the optimal parameters Θ of encoder_E, low-rank decoder_LR and sparse decoder_S en '、 and The trained first autoencoder LR' and second autoencoder S' are obtained.

[0054] Step 3. Hyperspectral anomaly detection:

[0055] (3.1) The hyperspectral image X is fed into the trained first autoencoder LR' and the second autoencoder S' to obtain the reconstructed low-rank part L re and the reconstructed sparse part S re ;

[0056] (3.2) Calculate residual data using RX detector The test results of LR And reconstruct the sparse part S re The test results of S , and the test result d is calculated according to the following formula LR and d S Perform nonlinear fusion to obtain the final detection result d f :

[0057]

[0058] Among them, d f represents the final detection result, and k represents the fusion coefficient.

[0059] Embodiment 2: The overall implementation steps of the hyperspectral anomaly detection method proposed in this embodiment are the same as those in Embodiment 1. Now, the joint optimization objective function is solved by using the augmented Lagrangian function, which is further described as follows:

[0060] First, construct the augmented Lagrangian function l:

[0061]

[0062] Where tr[] represents the trace of the matrix, Y represents the Lagrange multiplier, and μ represents the penalty coefficient;

[0063] Then, the alternating direction multiplier method is used to solve the optimal variable, that is, one variable is updated each time and the remaining variables are fixed:

[0064] i) Update L and fix the rest of the variables:

[0065]

[0066]

[0067] in, Ψ represents the singular value thresholding method.

[0068] ii) Update S and fix the remaining variables:

[0069]

[0070] in, Φ represents the soft thresholding method.

[0071] iii) Update E and fix the rest of the variables:

[0072]

[0073] Where Π represents L 2,1 norm minimization operator;

[0074] iv) Update Y:

[0075] Y=Y+μ(XLSE)

[0076] v) Update μ:

[0077] μ=min(ρμ,μ max )

[0078] vi) Update Θ en , and Fix the remaining variables:

[0079]

[0080] Finally, back propagation and Adam optimizer are used to optimize the parameter Θ en , and Get the optimal parameters Θ of encoder_E, low-rank decoder_LR and sparse decoder_S en '、 and

[0081] In this embodiment, it is preferred to set the initial values ​​of the low-rank component L, the sparse component S, and the residual component E to 0, and the initial value of μ to 10 -4 , μ max For 10 10 , ρ is 1.1, and the maximum number of iterations is 200; of course, the initial value of μ can also be set to 10 -3 , 10 -5 Or otherwise, the maximum number of iterations may also be preset to be no less than 100 times, so as to achieve the optimization purpose.

[0082] The effect of the present invention will be further described below in conjunction with experiments.

[0083] 1. Experimental conditions:

[0084] The experiments of the present invention were carried out in a hardware environment with a CPU main frequency of 2.00 GHz, a memory of 16 GB, and a software environment of Windows 10 operating system and Python 3.6.13.

[0085] 2. Experimental content:

[0086] This experiment evaluates the performance of the method of the present invention on the Texas Coast, Los Angeles and Salinas datasets. On each dataset, the detection effect of the method of the present invention is compared with that of seven mainstream methods from two perspectives: qualitative (including detection graphs and two ROC curves) and quantitative (including two AUC values).

[0087] The mainstream methods mainly include RX [paper: Adaptive multiple-band CFAR detection of an optical pattern with unknown spectral distribution], CRD [paper: Collaborative representation for hyperspectral anomaly detection], LSMAD [paper: A low-rank and sparse matrix decomposition-based mahalanobis distance method for hyperspectral anomaly detection], LSDM-MoG [paper: Low-rank and sparsedecomposition with mixture of Gaussian for hyperspectral anomaly detection], RGAE [paper: Hyperspectral anomaly detection with robust graph autoencoders], GAED [paper: Hyperspectral anomaly detection with guided autoencoder] and NJCR [paper: Nonnegative-constrained joint collaborative representation with union dictionary for hyperspectral anomaly detection].

[0088] 3. Simulation results and analysis:

[0089] Figure 4 The figure shows the comparison of the detection results of the method of the present invention and the mainstream methods on three data sets; (I)-(III) correspond to the Texas Coast, Los Angeles and Salinas data sets respectively, (a) represents the label, (b)-(h) correspond to RX, CRD, LSMAD, LSDM-MoG, RGAE, GAED, NJCR, a total of 7 existing methods, and (i) corresponds to the method of the present invention; specifically: Figure 4 In the figure, (Ia) represents the label of the Texas Coast dataset, (Ib) to (Ii) respectively represent the detection results of RX, CRD, LSMAD, LSDM-MoG, RGAE, GAED, NJCR and the method of the present invention on the Texas Coast dataset; (II-a) represents the label of the Los Angeles dataset, (II-b) to (II-i) respectively represent the detection results of RX, CRD, LSMAD, LSDM-MoG, RGAE, GAED, NJCR and the method of the present invention on the Los Angeles dataset; (III-a) represents the label of the Salinas dataset, (III-b) to (III-i) respectively represent the detection results of RX, CRD, LSMAD, LSDM-MoG, RGAE, GAED, NJCR and the method of the present invention on the Salinas dataset.

[0090] By comparing the existing methods, the method of the present invention and the label, it can be seen that the method of the present invention has the best effect and is significantly better than the mainstream methods in the existing 7 prior arts.

[0091] Figure 5 and Figure 6 The results of the proposed method and the mainstream method on three data sets are Curve and Curve; where (I) represents the Texas Coast dataset, (II) represents the Los Angeles dataset, and (III) represents the Salinas dataset.

[0092] The larger the area under the curve, the better the detection performance. Figure 5 It can be seen that the detection effect of the method of the present invention is the best; The smaller the area under the curve, the lower the false alarm rate. Figure 6 It can be seen that the false alarm rate of the method of the present invention is the lowest and the performance is the best.

[0093] The following table (Table 1) lists the AUC values ​​of the method of the present invention and seven mainstream methods on three data sets (i.e. and ). The AUC value is between 0 and 1. The closer its value is to 1, the better the detection effect is. The closer the value is to 0, the lower the false alarm rate is, that is, the better the performance is.

[0094] Table 1 Comparison of AUC values ​​of the proposed method and seven mainstream methods on three data sets

[0095]

[0096] As can be seen from the table above, on three different data sets (Texas Coast, Los Angeles, and Salinas), the method of the present invention has the largest and the smallest It is further shown that the method of the present invention has the highest detection accuracy and the lowest false alarm rate compared with the existing methods, that is, the performance is the best.

[0097] The above experimental results prove the correctness and effectiveness of the method proposed in the present invention.

[0098] Parts of the present invention that are not described in detail belong to common knowledge among those skilled in the art.

[0099] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Obviously, for professionals in this field, after understanding the content and principles of the present invention, they may make various modifications and changes in form and details without departing from the principles and structures of the present invention. However, these modifications and changes based on the ideas of the present invention are still within the scope of protection of the claims of the present invention.

Claims

1. A hyperspectral anomaly detection method based on low-rank and sparse prior constrained autoencoder, characterized in that: The steps are as follows: (1) Constructing an autoencoder with low rank and sparse prior constraints: The first autoencoder LR' is composed of the encoder_E and the low-rank decoder_LR, and the second autoencoder S' is composed of the encoder_E and the sparse decoder_S, and the first autoencoder LR' and the second autoencoder S' share the encoder_E, so as to obtain the first autoencoder LR' and the second autoencoder S' with low-rank and sparse prior constraints; (2) End-to-end joint optimization: (2.1) Decompose the hyperspectral image X into a low-rank component L, a sparse component S and a residual component E; (2.2) Using the low-rank component L and the sparse component S as labels, we establish the reconstruction loss l of the first autoencoder LR' and the second autoencoder S' respectively. LR and l S , and train them; (2.3) The decomposition of the hyperspectral image X and the training of the autoencoder are integrated into a unified framework to achieve joint optimization in an end-to-end manner, and the following joint optimization objective function is obtained: stX=L+S+E Among them, λ, β, γ, α and ν represent the coefficients of weighing each regularization term, which are the first weighing coefficient, the second weighing coefficient, the third weighing coefficient, the fourth weighing coefficient and the fifth weighing coefficient respectively; ||·|| * ,||·||1and||·|| 2,1 Represent the nuclear norm, L1 norm and L 2,1 norm; st represents the constraint condition; (2.4) The joint optimization objective function is solved using the augmented Lagrangian function, and back propagation and Adam optimizer are used for optimization to obtain the optimal parameters Θ of encoder_E, low-rank decoder_LR and sparse decoder_S. en '、 and Get the trained first autoencoder LR' and second autoencoder S'; (3) Hyperspectral anomaly detection: (3.1) The hyperspectral image X is fed into the trained first autoencoder LR' and the second autoencoder S' to obtain the reconstructed low-rank part L re and the reconstructed sparse part S re ; (3.2) Calculate residual data using RX detector The test results of LR And reconstruct the sparse part S re The test results of S , and the test result d is calculated according to the following formula LR and d S Perform nonlinear fusion to obtain the final detection result d f : Among them, d f represents the final detection result, and k represents the fusion coefficient.

2. The method according to claim 1, characterized in that: The number of nodes of the encoder_E, low-rank decoder_LR and sparse decoder_S in step (1) are {d, round(d / 2), round(d / 4), n}, {n, round(d / 4), round(d / 2), d} and {n, round(d / 4), round(d / 2), d} respectively, where round represents rounding operation, d represents the number of hyperspectral image bands, and n represents the number of hidden layer nodes.

3. The method according to claim 1, characterized in that: In step (2.1), the hyperspectral image X is decomposed into a low-rank component L, a sparse component S, and a residual component E according to the following objective function:

4. The method according to claim 1, characterized in that: In step (2.2), the reconstruction loss l of the first autoencoder LR' and the second autoencoder S' is established LR and l S , and train it as follows: Where b represents the batch size, and x i represent the low-rank component L, the sparse component S and the i-th pixel in the hyperspectral image X, respectively, and f en (·;·), and denote the encoding operation of encoder_E, the decoding operation of low-rank decoder_LR and the decoding operation of sparse decoder_S, Θ respectively. en , and denote the parameters of encoder_E, low-rank decoder_LR and sparse decoder_S respectively, represents the square of L2 norm; During the training process, the batch size b is set to the total number of pixels N of the hyperspectral image X, and the loss functions of the first autoencoder LR' and the second autoencoder S' are transformed into: in, Represents the square of the Frobenius norm.

5. The method according to claim 4, characterized in that: The loss function l in the first autoencoder LR' and the second autoencoder S' LR and l S Ignore coefficient The joint optimization objective function described in step (2.3) is converted into the following form: stX=L+S+E.

6. The method according to claim 1, characterized in that: The joint optimization objective function is solved by using the augmented Lagrangian function as described in step (2.4), which is implemented as follows: (2.4.1) Construct the augmented Lagrangian function l: Among them, tr[] represents the trace of the matrix, Y represents the Lagrange multiplier, and μ represents the penalty coefficient; (2.4.2) The alternating direction multiplier method is used to solve the optimal variable, that is, one variable is updated each time and the other variables are fixed: i) Update L and fix the rest of the variables: in, Ψ represents the singular value thresholding method, ii) Update S and fix the remaining variables: in, Φ represents the soft thresholding method, iii) Update E and fix the rest of the variables: Where Π represents L 2,1 norm minimization operator; iv) Update Θ en , and Fix the remaining variables: (2.4.3) Optimize the parameter Θ using back propagation and Adam optimizer en , and Get the optimal parameters Θ of encoder_E, low-rank decoder_LR and sparse decoder_S en '、 and

Images