Nothing Special   »   [go: up one dir, main page]

CN109785242B - Hyperspectral image unmixing method based on band-by-band generalized bilinear model - Google Patents

Hyperspectral image unmixing method based on band-by-band generalized bilinear model Download PDF

Info

Publication number
CN109785242B
CN109785242B CN201811097454.XA CN201811097454A CN109785242B CN 109785242 B CN109785242 B CN 109785242B CN 201811097454 A CN201811097454 A CN 201811097454A CN 109785242 B CN109785242 B CN 109785242B
Authority
CN
China
Prior art keywords
matrix
representing
band
hyperspectral image
formula
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201811097454.XA
Other languages
Chinese (zh)
Other versions
CN109785242A (en
Inventor
李畅
刘羽
成娟
宋仁成
陈强
彭虎
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Hefei University of Technology
Original Assignee
Hefei University of Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Hefei University of Technology filed Critical Hefei University of Technology
Priority to CN201811097454.XA priority Critical patent/CN109785242B/en
Publication of CN109785242A publication Critical patent/CN109785242A/en
Application granted granted Critical
Publication of CN109785242B publication Critical patent/CN109785242B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Image Processing (AREA)

Abstract

The invention discloses a hyperspectral image unmixing method based on a band-by-band generalized bilinear model, which comprises the following steps of: 1. establishing a band-by-band generalized bilinear model; 2. establishing an optimization model corresponding to the unmixing method based on a Bayes maximum posterior criterion and a regularization theory; 3. and solving the optimization model by adopting an alternative method multiplier method. The method not only can take the Gaussian noises of different Gaussian levels of different wave bands of the hyperspectral image into account, but also can take mixed noises widely existing in the real hyperspectral image into account, namely, the hyperspectral image not only contains the Gaussian noises, but also contains impulse noises, stripes, dead pixels, dead lines and the like, so the method has the advantage of robustness to the mixed noises in the real hyperspectral image and the noises of different Gaussian levels of different wave bands.

Description

Hyperspectral image unmixing method based on band-by-band generalized bilinear model
Technical Field
The invention relates to the field of hyperspectral image unmixing, in particular to a band-by-band generalized bilinear hyperspectral image unmixing model and method.
Background
The hyperspectral image is acquired by an imager in hundreds of narrow continuous spectral bands. Due to the high spectral resolution, this inevitably leads to the problem of mixing pixels, so that different substances occupy the same pixel. The presence of mixed pixels has a large impact on many applications such as object recognition, sub-pixel mapping and classification, etc. Hyperspectral image unmixing decomposes mixed pixels into a series of pure substances (i.e., end members) and corresponding proportions (i.e., abundances).
Linear hybrid model (LMM) is a simple and widely used model that assumes that each incident ray interacts with only one substance, and thus each pixel is a linear combination of end-members. However, LMMs fail when there is intimate mixing, topographical factors, or multiple scattering effects. Nonlinear hybrid models (NLMMs) provide some methods to overcome the above problems, and can be broadly classified into two categories. The first category includes some flexible models based on signal processing, including post-nonlinear models, neural network models, and kernel models, among others. The second category includes some physics-based models including a tight-mix model, a bilinear-mix model (BMM), and a multilinear-mix model, among others. Wherein the BMM only takes into account second order mixing effects and does not take into account higher order mixing effects. This is because the higher order mixing effect not only contributes little to improving the unmixing accuracy, but also greatly increases the computational complexity. Several representative BMM models have been proposed. The Nascimento model is an extended LMM model with virtual end members, the Fan Model (FM) is a truncated Taylor expansion of a nonlinear mixing function, and the Generalized Bilinear Model (GBM) can be regarded as popularization of the LMM and the FM. Different algorithms have been proposed for GBM unmixing, and Halimi et al propose Bayesian algorithms for estimating the abundance and noise variance of the GBM model, and they also propose a Gradient Descent (GDA) -based pixel-by-pixel unmixing method. Furthermore, Yokoya et al propose a semi-non-negative matrix factorization (semi-NMF) based optimization method for GBM unmixing, and Li et al propose a boundary projection optimal gradient method for GBM unmixing. Most GBM-based unmixing methods imply the assumption that hyperspectral images are only contaminated by Additive White Gaussian Noise (AWGN) and that the noise levels of different bands are assumed to be the same.
However, there are still two challenges to the actual GBM-based hyperspectral image unmixing. On one hand, each wave band of the hyperspectral image contains AWGN with different intensities, and on the other hand, the mixed noise problem widely existing in the real hyperspectral image is solved, namely, the hyperspectral image not only contains Gaussian noise, but also contains impulse noise, stripes, dead pixels, dead lines and the like.
Disclosure of Invention
In order to overcome the corresponding technical defects, the invention provides a method for unmixing a hyperspectral image based on a band-by-band generalized bilinear model, which can solve the influence of Gaussian noises with different intensities of different bands and widely-existing mixed noises on the unmixing of the hyperspectral image, so that the unmixing can adapt to the problem of complex noises in the actual hyperspectral image, and the accuracy and the robustness of the unmixing are improved.
The invention adopts the following technical scheme for solving the technical problems:
the invention relates to a hyperspectral image unmixing method based on a band-by-band generalized bilinear model, which is characterized by comprising the following steps of:
step 1, establishing a band-by-band generalized bilinear model by using a formula (1):
Yi=(EA)i+(FB)i+Si+Ni(1)
in the formula (1), YiA band corresponding to the ith row in the pixel matrix Y representing the hyperspectral image, i ═ 1,2, …, D, and Y ∈ RD×PD and P respectively represent the total number of bands in the spectral dimension and the total number of pixels in the spatial dimension of the hyperspectral image, E ═ E1,e1,…,ej,…,eM]∈RD×MAn end member matrix representing the hyperspectral image, wherein ejRepresents the jth end member in the end member matrix E, j is 1,2, …, M represents the total number of end members in the hyperspectral image, A is [ a ]1,a2,…,ak,…,aP]∈RM×PRepresenting an abundance matrix of the hyperspectral image, wherein akAn abundance vector representing the kth pixel in the abundance matrix a, k ═ 1,2, …, P, F ═ e1⊙e2,...,e1⊙eM,e2⊙e3,...,e2⊙eM,...,eM-1⊙eM]∈RD ×M(M-1)/2Represents a quadratic interaction end-member matrix, wherein ⊙ represents the Hadamard product, B ∈ RM(M-1)/2×PRepresenting a quadratic interaction abundance matrix, SiA band corresponding to the ith row of a sparse noise matrix S representing the hyperspectral image, and S ∈ RD×P,NiThe wave band corresponding to the ith row of the dense noise matrix N of the hyperspectral image is represented and meets the requirement
Figure GDA0002589923900000021
The gaussian noise representing the band corresponding to the ith row follows a gaussian distribution of different strengths with zero mean,
Figure GDA0002589923900000022
representing the variance of the Gaussian noise of the band corresponding to the ith row, IpRepresenting an identity matrix having p elements in the diagonal N ∈ RD×P
Step 2, based on Bayes maximum posterior criterion and regularization theory, establishing an optimization model corresponding to the unmixing method by using formula (2):
Figure GDA0002589923900000023
in the formula (2), W is a diagonal matrix, and the diagonal elements
Figure GDA0002589923900000024
min is a minimization operator, | · | countFRepresenting the Frobenius norm of the matrix, | S | | luminance1=∑i.j|Si,jI represents the sum of absolute values of the ith row and the jth column in the sparse noise matrix S, lambda represents a regularization parameter, s.t. represents a constraint condition, and C represents an upper bound matrix of the secondary interaction abundance matrix B;
step 3, solving the optimization model by adopting an alternative method multiplier method to obtain the abundance matrix A, the secondary interaction abundance matrix B and the sparse noise matrix S:
step 3.1, three auxiliary variables V are introduced1、V2And V3Rewriting the optimization model to obtain a rewritten optimization model shown as a formula (3):
Figure GDA0002589923900000031
in formula (3), X is equal to V1,V2Or V3
Figure GDA0002589923900000032
Representing a non-negative limit R+Illustrative function of, Xi,jThe ith row and jth column element of variable X when X isi,jBelonging to a non-negative range R+At the time of the above-mentioned operation,
Figure GDA0002589923900000033
is 0, otherwise is + ∞;
Figure GDA0002589923900000034
is the interval [0, C]When X is an indicative function ofi,jIs of [0, Ci,j]At the time of the above-mentioned operation,
Figure GDA0002589923900000035
is 0, otherwise is + ∞; ci,jThe ith row and the jth column of elements of the upper bound matrix of the secondary interaction abundance matrix B;
obtaining an augmented Lagrangian function corresponding to the rewritten optimization model by using a formula (4):
Figure GDA0002589923900000036
in the formula (4), the reaction mixture is,
Figure GDA0002589923900000037
mu represents a penalty factor which is given by,
Figure GDA0002589923900000038
the lagrange multiplier is represented by a number of lagrange multipliers,
Figure GDA0002589923900000039
scaling matrix representing stacked Lagrange multipliers, Λ1、Λ2And Λ3Respectively representing scaling matrixes of Lagrange multipliers corresponding to the sparse noise matrix S, the abundance matrix A and the secondary interaction abundance matrix B;
Figure GDA00025899239000000310
i denotes a unit matrix of the cell,
Figure GDA00025899239000000311
g (V, Q) represents the objective function of the rewritten optimization model and has:
Figure GDA00025899239000000312
step 3.2, defining the current iteration number as k, and initializing k to be 0; initialization Ak、Bk、Sk、V1 k、V2 k、V3 k
Figure GDA00025899239000000313
Figure GDA00025899239000000314
Step 3.3: updating abundance matrix A of k +1 th iteration by using formula (6)k+1
Figure GDA0002589923900000041
Step 3.4: updating the quadratic interaction abundance matrix B of the (k + 1) th iteration by using the formula (7)k+1
Figure GDA0002589923900000042
Step 3.5: updating sparse noise matrix S for the (k + 1) th iteration using equation (8)k+1
Figure GDA0002589923900000043
In the formula (8)
Figure GDA0002589923900000044
Figure GDA0002589923900000045
Is a soft contraction operator, and
Figure GDA0002589923900000046
representing a threshold, sgn (x) representing a sign function of x, and max (·) representing a function of taking a larger value;
step 3.6: updating the first auxiliary variable V of the (k + 1) th iteration using equation (9)1 k+1
Figure GDA0002589923900000047
Step 3.7: updating the second auxiliary variable of the (k + 1) th iteration using equation (10)
Figure GDA0002589923900000048
Figure GDA0002589923900000049
Step 3.8: updating the third auxiliary variable for the (k + 1) th iteration using equation (11)
Figure GDA00025899239000000410
Figure GDA00025899239000000411
In formula (11), min (-) represents a smaller value function;
step 3.9: updating sparse noise matrix S, abundance matrix A and quadratic cross of the (k + 1) th iteration by using formula (12)Scaling matrix of Lagrange multiplier corresponding to mutual abundance matrix B
Figure GDA00025899239000000412
And
Figure GDA00025899239000000413
Figure GDA00025899239000000414
step 3.10: updating the original error r of the (k + 1) th iteration by using equation (13)k+1And the dual error d of the (k + 1) th iterationk+1
Figure GDA00025899239000000415
Step 3.11: judging a convergence condition:
if it is
Figure GDA00025899239000000416
And is
Figure GDA00025899239000000417
And if not, assigning k +1 to k, and executing the step 3.3 in a rotating manner.
Compared with the prior art, the invention has the beneficial effects that:
1. the invention establishes a band-by-band generalized bilinear model which can take Gaussian noise with different intensities of different bands and various types of noise into consideration. In addition, the optimization model corresponding to the unmixing method is established, and the optimization model is solved by adopting an alternative method multiplier method. The model and the unmixing method can be suitable for the problem of complex noise in practice, so that the unmixing robustness and accuracy are higher.
2. Firstly, the method not only can take the Gaussian noises of different Gaussian levels of different wave bands of the hyperspectral image into consideration, but also can take the mixed noise widely existing in the real hyperspectral image into consideration, namely, the hyperspectral image not only contains the Gaussian noise, but also contains impulse noise, stripes, dead pixels, dead lines and the like. The noise is divided into dense noise and sparse noise according to the characteristics of the noise, and a band-by-band generalized bilinear model is established on the basis of the assumption that the noise levels of different bands are different. Secondly, according to the Bayes maximum posterior criterion, on the basis of the assumption that each wave band contains independent Gaussian noise with different intensities, a unified frame based on band-by-band generalized bilinear de-mixing is established, and then the sparse characteristic of the sparse noise is taken into consideration according to the regularization theory, so that the de-mixing method established by the invention has the advantage of robustness to the mixed noise in the real hyperspectral image and the different Gaussian level noise of different wave bands. Finally, since the alternative method multiplier method has been widely applied to solving the optimization problem under the constraint condition and obtains satisfactory effect, the alternative method multiplier method is also adopted to solve the optimization model.
Drawings
FIG. 1 is a flow chart of an embodiment of the present invention.
Detailed Description
In this embodiment, as shown in fig. 1, a method for unmixing a hyperspectral image based on a band-by-band generalized bilinear model mainly includes 3 steps: 1. establishing a band-by-band generalized bilinear model; 2. establishing an optimization model corresponding to the unmixing method based on a Bayes maximum posterior criterion and a regularization theory; 3. solving the optimization model by adopting an alternative method multiplier method; specifically, the method comprises the following steps:
step 1, in order to overcome the problem that the traditional GBM model is sensitive to different levels of Gaussian noise of different wave bands and various types of noise, the noise in a real hyperspectral image is divided into two categories, namely dense noise and sparse noise. Dense noise means that most of the hyperspectral image is contaminated by noise, which mainly comprises gaussian noise. Sparse noise means that a small part of a hyperspectral image is polluted by sparse noise, and mainly comprises impulse noise, dead pixels, dead lines and stripes. And assuming that the Gaussian noise levels of different bands are different, establishing a band-by-band generalized bilinear model by using the formula (1):
Yi=(EA)i+(FB)i+Si+Ni(1)
in the formula (1), YiA band corresponding to the ith row in the pixel matrix Y representing the hyperspectral image, i ═ 1,2, …, D, and Y ∈ RD×PD and P denote the total number of bands in the spectral dimension and the total number of pixels in the spatial dimension of the hyperspectral image, respectively, E ═ E1,e1,…,ej,…,eM]∈RD×MAn end-member matrix representing a hyperspectral image, wherein ejRepresents the jth end member in the end member matrix E, j is 1,2, …, M represents the total number of end members in the hyperspectral image, A is [ a ]1,a2,…,ak,…,aP]∈RM×PAn abundance matrix representing a hyperspectral image, wherein akAn abundance vector representing the kth pixel in the abundance matrix a, k ═ 1,2, …, P, F ═ e1⊙e2,...,e1⊙eM,e2⊙e3,...,e2⊙eM,...,eM-1⊙eM]∈RD×M(M-1)/2Represents a quadratic interaction end-member matrix, wherein ⊙ represents the Hadamard product, B ∈ RM(M-1)/2×PRepresenting a quadratic interaction abundance matrix, SiRepresents the wave band corresponding to the ith row of the sparse noise matrix S of the hyperspectral image, and S ∈ RD×P,NiThe wave band corresponding to the ith row of the dense noise matrix N of the hyperspectral image is represented and meets the requirement
Figure GDA0002589923900000061
The gaussian noise representing the band corresponding to the ith row follows a gaussian distribution of different strengths with zero mean,
Figure GDA0002589923900000062
representing the variance of the Gaussian noise of the band corresponding to the ith row, IpRepresenting an identity matrix having p elements in the diagonal N ∈ RD×P
Step 2, can get:
Figure GDA0002589923900000063
assuming that the white gaussian noise of each band is independent of each other, there are:
Figure GDA0002589923900000064
in the formula (2), m represents a constant, W is a diagonal matrix, and diagonal elements
Figure GDA0002589923900000065
As can be seen from equation (2), the greater the variance of gaussian noise of each band, the smaller the weight of the band. Based on the Bayes maximum posterior criterion, the abundance matrix A, the secondary interaction abundance matrix B and the sparse noise matrix S can be converted into the following optimization problems:
Figure GDA0002589923900000066
in equation (3), p (EA + FB + S) represents the prior distribution of EA + FB + S. Since sparse noise has systematic characteristics, the noise can be represented by0Norm characterization, but based on l0The norm optimization problem is usually NP hard, and l is used for this purpose1Norm to replace l0Norm, and the traditional GBM model satisfies the constraint condition in the formula (2), so the optimization model corresponding to the unmixing method is established by using the formula (4):
Figure GDA0002589923900000071
in formula (4), min is the minimization operator, | | · the luminanceFRepresenting the Frobenius norm of the matrix, | S | | luminance1=∑i.j|Si,jI represents the sum of absolute values of the ith row and the jth column in the sparse noise matrix S, lambda represents a regularization parameter, s.t. represents a constraint condition, and C represents an upper bound matrix of the secondary interaction abundance matrix B;
and 3, because the alternative method multiplier method is widely applied to solving the optimization problem under the constraint condition, a satisfactory effect is obtained. Therefore, the invention adopts an alternative method multiplier method to solve an optimization model, and obtains an abundance matrix A, a secondary interaction abundance matrix B and a sparse noise matrix S:
step 3.1, three auxiliary variables V are introduced1、V2And V3And rewriting the optimization model to obtain a rewritten optimization model shown as a formula (5):
Figure GDA0002589923900000072
in formula (5), X is equal to V1,V2Or V3
Figure GDA0002589923900000073
Representing a non-negative limit R+Illustrative function of, Xi,jThe ith row and jth column element of variable X when X isi,jBelonging to a non-negative range R+At the time of the above-mentioned operation,
Figure GDA0002589923900000074
is 0, otherwise is + ∞;
Figure GDA0002589923900000075
is the interval [0, C]When X is an indicative function ofi,jIs of [0, Ci,j]At the time of the above-mentioned operation,
Figure GDA0002589923900000076
is 0, otherwise is + ∞; ci,jThe ith row and the jth column of elements of the upper bound matrix of the secondary interaction abundance matrix B;
and (3) obtaining an augmented Lagrange function corresponding to the rewritten optimization model by using the formula (6):
Figure GDA0002589923900000077
in the formula (4), the reaction mixture is,
Figure GDA0002589923900000078
mu represents a penalty factor which is given by,
Figure GDA0002589923900000079
the lagrange multiplier is represented by a number of lagrange multipliers,
Figure GDA00025899239000000710
scaling matrix representing stacked Lagrange multipliers, Λ1、Λ2And Λ3Respectively representing scaling matrixes of Lagrange multipliers corresponding to the sparse noise matrix S, the abundance matrix A and the quadratic interaction abundance matrix B,
Figure GDA0002589923900000081
i denotes a unit matrix of the cell,
Figure GDA0002589923900000082
g (V, Q) represents the objective function of the rewritten optimization model and has:
Figure GDA0002589923900000083
step 3.2, defining the current iteration number as k, and initializing k to be 0; initialization Ak、Bk、Sk、V1 k、V2 k、V3 k
Figure GDA0002589923900000084
Figure GDA0002589923900000085
Step 3.3: updating abundance matrix A of k +1 th iteration by using formula (8)k+1
Figure GDA0002589923900000086
Step 3.4: updating the quadratic interaction abundance matrix B of the (k + 1) th iteration by using the formula (9)k+1
Figure GDA0002589923900000087
Step 3.5: updating sparse noise matrix S for the (k + 1) th iteration using equation (10)k+1
Figure GDA0002589923900000088
In the formula (10)
Figure GDA0002589923900000089
Figure GDA00025899239000000810
Is a soft contraction operator, and
Figure GDA00025899239000000811
representing a threshold, sgn (x) representing a sign function of x, and max (·) representing a function of taking a larger value;
step 3.6: updating the first auxiliary variable V of the (k + 1) th iteration using equation (11)1 k+1
Figure GDA00025899239000000812
Step 3.7: updating the second auxiliary variable of the (k + 1) th iteration using equation (12)
Figure GDA00025899239000000813
Figure GDA0002589923900000091
Step 3.8: updating the third auxiliary variable V of the (k + 1) th iteration using equation (13)3 k+1
Figure GDA0002589923900000092
In formula (11), min (-) represents a smaller value function;
step 3.9: updating sparse noise matrix S, abundance matrix A and quadratic cross of the (k + 1) th iteration by using formula (14)Scaling matrix of Lagrange multiplier corresponding to mutual abundance matrix B
Figure GDA0002589923900000093
And
Figure GDA0002589923900000094
Figure GDA0002589923900000095
step 3.10: updating the original error r of the (k + 1) th iteration by using equation (15)k+1And the dual error d of the (k + 1) th iterationk+1
Figure GDA0002589923900000096
Step 3.11: judging a convergence condition:
if it is
Figure GDA0002589923900000097
And is
Figure GDA0002589923900000098
And if not, assigning k +1 to k, and executing the step 3.3 in a rotating manner. In addition, the choice of μ has a large influence on the convergence speed, and μ is updated so that the ratio of the norm of the original error and the dual error is kept within a certain range and eventually converges to 0. In the iterative process, the diagonal matrix W also needs to be known in advance, and for this reason, the HySime algorithm is adopted by the invention to estimate the noise level of each band. The core idea of the method is to assume that the bands of the hyperspectral image have strong correlation and can be obtained by solving through a multi-regression theory.
In specific implementation, simulation experiments are carried out to verify the effectiveness of the algorithm (NU-BGBM) provided by the invention, the comparison algorithm comprises a fully constrained least square method (FCLS), a gradient descent method (GDA), a semi-nonnegative matrix factorization method (semi-NMF) and a Boundary Projection Optimal Gradient Method (BPOGM), and the accuracy of unmixing is measured by adopting a Root Mean Square Error (RMSE). In general, the smaller the RMSE, the higher the accuracy of the unmixing.
A spectral database from the United States Geological Survey (USGS) will be used, which has 224 bands ranging from 0.38 μm to 2.5 μm. The method comprises the steps of randomly selecting 6 spectrums as end members, enabling 64 x 64 pixels of a synthesized hyperspectral image to be free of pure pixels, dividing the hyperspectral image into 8 x 8 small blocks, filling the pixels in each small block with one of the 6 end members randomly, then adopting a 9 x 9 space low-pass filter, and enabling the pixels with abundance of more than 80% to be replaced by the average of all the end members, so that the synthesized hyperspectral image does not contain the pure pixels. And generating a bilinear abundance matrix according to the band-by-band generalized bilinear model. Then different types of noise (1) gaussian noise are added: adding zero-mean Gaussian white noise to HSI of all wave bands, wherein the signal-to-noise ratio of each wave band is a random value from 10dB to 50 dB; (2) impulse noise. 30% of impulse noise is added into the wave band of 60-70; (3) and (3) dead line: band 120-130 joins the dead line.
In order to simulate various noises actually existing in the hyperspectral image as much as possible, experiments are carried out on one kind of noise, two kinds of noise and three kinds of noise. Table 1 gives the RMSE for different unmixing methods under different noise situations. It can be seen from table 1 that when the hyperspectral image is only polluted by gaussian noise, the unmixing effect of the algorithm provided by the invention is superior to that of other comparison algorithms, because only the method provided by the invention can take the gaussian noise of different wave bands and different levels of the hyperspectral image into account. When the hyperspectral image is only polluted by impulse noise or dead lines, the unmixing effect obtained by the method provided by the invention is obviously better than that obtained by other comparison methods, because the traditional unmixing algorithm is sensitive to sparse noise, only the method provided by the invention can take the sparse characteristic of the sparse noise into consideration. When the hyperspectral image is polluted by two or three kinds of noise, the effect of the unmixing method provided by the invention is better than that of other comparison methods, because only the method provided by the invention can take the Gaussian noise with different wave bands and different intensities and various mixed noises into consideration. Therefore, compared with a comparison method, the method provided by the invention is more robust to unmixing under various complex noise conditions, and can improve the accuracy of unmixing.
Table 1 compares the RMSE (10) for different unmixing methods under different noise conditions-2)
Type of noise FCLS GDA semi-NMF BPOGM NU-BGBM
Gaussian noise 7.103 6.053 5.520 5.157 0.990
Impulse noise 7.123 6.395 6.161 5.403 0.167
Dead line 6.812 5.773 5.796 5.436 0.171
Gaussian noise and impulse noise 8.411 7.781 7.651 7.065 1.004
Gaussian noise, dead line 8.084 7.185 7.197 6.996 1.003
Impulse noise, dead line 7.941 7.254 7.469 6.962 0.296
Gaussian noise, impulse noise, dead line 9.010 8.395 8.609 8.216 1.021

Claims (1)

1. A hyperspectral image unmixing method based on a band-by-band generalized bilinear model is characterized by comprising the following steps:
step 1, establishing a band-by-band generalized bilinear model by using a formula (1):
Yi=(EA)i+(FB)i+Si+Ni(1)
in the formula (1), YiA band corresponding to the ith row in the pixel matrix Y representing the hyperspectral image, i ═ 1,2, …, D, and Y ∈ RD ×PD and P respectively represent the total number of bands in the spectral dimension and the total number of pixels in the spatial dimension of the hyperspectral image, and E ═ E1,e1,…,ej,…,eM]∈RD×MAn end member matrix representing the hyperspectral image, wherein ejRepresents the jth end member in the end member matrix E, j is 1,2, …, M represents the total number of end members in the hyperspectral image, A is [ a ]1,a2,…,ak,…,aP]∈RM×PRepresenting an abundance matrix of the hyperspectral image, wherein akAn abundance vector representing the kth pixel in the abundance matrix a, k ═ 1,2, …, P, F ═ e1⊙e2,...,e1⊙eM,e2⊙e3,...,e2⊙eM,...,eM-1⊙eM]∈RD×M(M-1)/2Represents a quadratic interaction end-member matrix, wherein ⊙ represents the Hadamard product, B ∈ RM(M-1)/2×PRepresenting a quadratic interaction abundance matrix, SiA band corresponding to the ith row of a sparse noise matrix S representing the hyperspectral image, and S ∈ RD×P,NiThe wave band corresponding to the ith row of the dense noise matrix N of the hyperspectral image is represented and meets the requirement
Figure FDA0002589923890000013
The gaussian noise representing the band corresponding to the ith row follows a gaussian distribution of different strengths with zero mean,
Figure FDA0002589923890000014
represents the ith row pairVariance of Gaussian noise of the corresponding band, IpRepresenting an identity matrix having p elements in the diagonal N ∈ RD×P
Step 2, based on Bayes maximum posterior criterion and regularization theory, establishing an optimization model corresponding to the unmixing method by using formula (2):
Figure FDA0002589923890000011
in the formula (2), W is a diagonal matrix, and the diagonal elements
Figure FDA0002589923890000012
min is a minimization operator, | · | countFRepresenting the Frobenius norm of the matrix, | S | | luminance1=∑i.j|Si,jI represents the sum of absolute values of the ith row and the jth column in the sparse noise matrix S, lambda represents a regularization parameter, s.t. represents a constraint condition, and C represents an upper bound matrix of the secondary interaction abundance matrix B;
step 3, solving the optimization model by adopting an alternative method multiplier method to obtain the abundance matrix A, the secondary interaction abundance matrix B and the sparse noise matrix S:
step 3.1, three auxiliary variables V are introduced1、V2And V3Rewriting the optimization model to obtain a rewritten optimization model shown as a formula (3):
Figure FDA0002589923890000021
in formula (3), X is equal to V1,V2Or V3
Figure FDA0002589923890000022
Representing a non-negative limit R+Illustrative function of, Xi,jThe ith row and jth column element of variable X when X isi,jBelonging to a non-negative range R+At the time of the above-mentioned operation,
Figure FDA0002589923890000023
is 0, otherwise is + ∞;
Figure FDA0002589923890000024
is the interval [0, C]When X is an indicative function ofi,jIs of [0, Ci,j]At the time of the above-mentioned operation,
Figure FDA0002589923890000025
is 0, otherwise is + ∞; ci,jThe ith row and the jth column of elements of the upper bound matrix of the secondary interaction abundance matrix B;
obtaining an augmented Lagrangian function corresponding to the rewritten optimization model by using a formula (4):
Figure FDA0002589923890000026
in the formula (4), the reaction mixture is,
Figure FDA0002589923890000027
mu represents a penalty factor which is given by,
Figure FDA0002589923890000028
the lagrange multiplier is represented by a number of lagrange multipliers,
Figure FDA0002589923890000029
scaling matrix representing stacked Lagrange multipliers, Λ1、Λ2And Λ3Respectively representing scaling matrixes of Lagrange multipliers corresponding to the sparse noise matrix S, the abundance matrix A and the secondary interaction abundance matrix B;
Figure FDA00025899238900000210
i denotes a unit matrix of the cell,
Figure FDA00025899238900000211
g (V, Q) represents the objective function of the rewritten optimization model and has:
Figure FDA00025899238900000212
step 3.2, defining the current iteration number as k, and initializing k to be 0; initialization Ak、Bk、Sk、V1 k
Figure FDA00025899238900000213
Figure FDA00025899238900000214
Step 3.3: updating abundance matrix A of k +1 th iteration by using formula (6)k+1
Figure FDA00025899238900000215
Step 3.4: updating the quadratic interaction abundance matrix B of the (k + 1) th iteration by using the formula (7)k+1
Figure FDA00025899238900000216
Step 3.5: updating sparse noise matrix S for the (k + 1) th iteration using equation (8)k+1
Figure FDA0002589923890000031
In the formula (8)
Figure FDA0002589923890000032
Sτ[x]Sgn (x) max (| x | - τ,0) is the soft-shrink operator, and
Figure FDA0002589923890000033
representing a threshold, sgn (x) representing a sign function of x, and max (·) representing a function of taking a larger value;
step 3.6: update using equation (9)First auxiliary variable V of k +1 iterations1 k+1
Figure FDA0002589923890000034
Step 3.7: updating the second auxiliary variable of the (k + 1) th iteration using equation (10)
Figure FDA0002589923890000035
Figure FDA0002589923890000036
Step 3.8: updating the third auxiliary variable V of the (k + 1) th iteration using equation (11)3 k+1
Figure FDA0002589923890000037
In formula (11), min (-) represents a smaller value function;
step 3.9: updating scaling matrixes of Lagrange multipliers corresponding to the sparse noise matrix S, the abundance matrix A and the quadratic interaction abundance matrix B of the (k + 1) th iteration by using the formula (12)
Figure FDA0002589923890000038
And
Figure FDA0002589923890000039
Figure FDA00025899238900000310
step 3.10: updating the original error r of the (k + 1) th iteration by using equation (13)k+1And the dual error d of the (k + 1) th iterationk +1
Figure FDA00025899238900000311
Step 3.11: judging a convergence condition:
if it is
Figure FDA00025899238900000313
And is
Figure FDA00025899238900000312
And if not, assigning k +1 to k, and executing the step 3.3 in a rotating manner.
CN201811097454.XA 2018-09-19 2018-09-19 Hyperspectral image unmixing method based on band-by-band generalized bilinear model Active CN109785242B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201811097454.XA CN109785242B (en) 2018-09-19 2018-09-19 Hyperspectral image unmixing method based on band-by-band generalized bilinear model

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201811097454.XA CN109785242B (en) 2018-09-19 2018-09-19 Hyperspectral image unmixing method based on band-by-band generalized bilinear model

Publications (2)

Publication Number Publication Date
CN109785242A CN109785242A (en) 2019-05-21
CN109785242B true CN109785242B (en) 2020-08-28

Family

ID=66496301

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201811097454.XA Active CN109785242B (en) 2018-09-19 2018-09-19 Hyperspectral image unmixing method based on band-by-band generalized bilinear model

Country Status (1)

Country Link
CN (1) CN109785242B (en)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110378356B (en) * 2019-07-16 2021-07-02 北京中科研究院 Fine-grained image identification method based on multi-target Lagrangian regularization
CN110490386B (en) * 2019-08-26 2022-12-06 苏州树森信息科技有限公司 Comprehensive energy scheduling method and comprehensive energy scheduling system
CN113591581B (en) * 2021-06-24 2024-02-13 南京航空航天大学 Spatio-temporal sub-pixel mapping based on consideration of variation differences between images
CN114596482B (en) * 2022-02-10 2023-04-07 复旦大学 Hyperspectral image nonlinear unmixing method based on extended multi-linear mixed model

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104899850A (en) * 2014-03-03 2015-09-09 五邑大学 High-spectrum image unmixing method based on weighted joint sparse regression
CN104952050A (en) * 2015-07-07 2015-09-30 西安电子科技大学 Self-adaptive hyperspectral image unmixing method based on region segmentation
CN105138860A (en) * 2015-10-23 2015-12-09 武汉大学 Hyperspectral nonlinear demixing method based on boundary projection optimal gradient
CN106056524A (en) * 2016-05-25 2016-10-26 天津商业大学 Hyper-spectral image nonlinear de-mixing method based on differential search
CN106778530A (en) * 2016-11-28 2017-05-31 复旦大学 A kind of hyperspectral image nonlinear solution mixing method based on bilinearity mixed model

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104899850A (en) * 2014-03-03 2015-09-09 五邑大学 High-spectrum image unmixing method based on weighted joint sparse regression
CN104952050A (en) * 2015-07-07 2015-09-30 西安电子科技大学 Self-adaptive hyperspectral image unmixing method based on region segmentation
CN105138860A (en) * 2015-10-23 2015-12-09 武汉大学 Hyperspectral nonlinear demixing method based on boundary projection optimal gradient
CN106056524A (en) * 2016-05-25 2016-10-26 天津商业大学 Hyper-spectral image nonlinear de-mixing method based on differential search
CN106778530A (en) * 2016-11-28 2017-05-31 复旦大学 A kind of hyperspectral image nonlinear solution mixing method based on bilinearity mixed model

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
Hyperspectral image denoising with segmentation-based low rank representation;Jiayi Ma 等;《2016 Visual Communications and Image Processing》;20170105;1-4 *
Hyperspectral Unmixing with Bandwise Generalized Bilinear Model;Chang Li 等;《remote sensing》;20181009;1-19 *
Robust Sparse Hyperspectral Unmixing With ℓ2,1 Norm;Yong Ma 等;《IEEE Transactions on Geoscience and Remote Sensing》;20161215;1227-1239 *
Sparse unmixing of hyperspectral data with noise level estimation;Chang Li 等;《Remote Sensing》;20171113;1-28 *
高光谱遥感影像处理中的若干关键技术研究;李畅;《中国博士学位论文全文数据库 工程科技Ⅱ辑》;20190115;C028-25 *

Also Published As

Publication number Publication date
CN109785242A (en) 2019-05-21

Similar Documents

Publication Publication Date Title
CN109785242B (en) Hyperspectral image unmixing method based on band-by-band generalized bilinear model
Yao et al. Nonconvex-sparsity and nonlocal-smoothness-based blind hyperspectral unmixing
Rabin et al. Adaptive color transfer with relaxed optimal transport
CN110458777B (en) Hyperspectral image denoising method, system and medium based on adaptive rank correction
US9317929B2 (en) Decomposition apparatus and method for refining composition of mixed pixels in remote sensing images
CN109241843B (en) Space-spectrum combined multi-constraint optimization non-negative matrix unmixing method
US20180324465A1 (en) Edge-aware spatio-temporal filtering and optical flow estimation in real time
US8908989B2 (en) Recursive conditional means image denoising
CN110428454B (en) Hyperspectral unmixing method and device, electronic equipment and storage medium
Jiang et al. Fog density estimation and image defogging based on surrogate modeling for optical depth
CN113763271B (en) High-quality spectrum denoising method based on physical noise model
Singh et al. Illumination estimation for nature preserving low-light image enhancement
CN103778598B (en) Disparity map ameliorative way and device
CN104484884A (en) Intrinsic image decomposition method based on multi-scale L0 sparse constraint
CN102646267B (en) Degraded image restoration method and system
CN114998750B (en) Method, device, equipment and medium for removing random length strips of remote sensing image
CN110070506B (en) Video rain removing method based on multi-scale mixed index model
CN106778530B (en) Hyperspectral image nonlinear unmixing method based on bilinear hybrid model
Jha et al. l2‐norm‐based prior for haze‐removal from single image
CN111340697A (en) Clustering regression-based image super-resolution method
JP7463186B2 (en) Information processing device, information processing method, and program
CN113781375B (en) Vehicle-mounted vision enhancement method based on multi-exposure fusion
CN106296749B (en) RGB-D image eigen decomposition method based on L1 norm constraint
CN111402173A (en) Hybrid noise removing method and device, electronic equipment and storage medium
CN112784747B (en) Multi-scale eigen decomposition method for hyperspectral remote sensing image

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant