Hyperspectral image sparse unmixing method based on sigma norm

Abstract

The invention provides a hyperspectral image sparse unmixing method based on a sigma norm, and the method comprises the steps: obtaining a hyperspectral image, and carrying out the wave band selection of the hyperspectral image through a neighborhood grouping method; s2, obtaining a USGS spectrum library, and reserving spectral characteristics in the USGS spectrum library, the spectral characteristics having the same spectral band range as the spectral characteristics obtained in the step S1, so as to obtain a ground object spectrum library; the method comprises the following steps: segmenting a hyperspectral image into a plurality of superpixels to obtain a hyperspectral image matrix, and calculating a homogeneity coefficient of each superpixel; based on the super-pixel homogeneous coefficient, establishing a target function of a coarse mixing solution stage based on a sigma norm; learning an abundance matrix of the rough mixing and solving stage based on the target function of the rough mixing and solving stage; obtaining a pruning spectrum library based on the abundance matrix and the ground object spectrum library in the rough unmixing stage; and performing fine unmixing based on the abundance matrix in the coarse unmixing stage and the pruning spectrum library to output a final abundance matrix. The spectral image abundance estimation method has the characteristics of noise robustness and easy optimization, and can effectively improve the accuracy of spectral image abundance estimation.

Description

Technical Field

[0001] The present invention relates to the technical field of hyperspectral image unmixing, and in particular to a hyperspectral image sparse unmixing method based on σ norm. Background Art

[0002] Hyperspectral images are usually taken by hyperspectral imagers carried on satellite platforms or high-altitude airborne, and they have the characteristics of integrating spatial and spectral information. On the one hand, at any spatial position point, the spectral characteristics of hyperspectral images can intuitively represent the radiation values ​​of different objects within the band range. On the other hand, any band of hyperspectral images can also present the two-dimensional spatial information of the distribution of objects. At present, hyperspectral imagers have high spectral resolution and can collect radiation characteristics on hundreds of spectral bands within the surface area. However, due to the limitations of sensor technology, its spatial resolution is often low, resulting in a large number of mixed pixels in the image. These mixed pixels are composed of the spectral characteristics of multiple substances superimposed in a certain pattern, which makes the analysis and interpretation of hyperspectral images complicated.

[0003] The purpose of hyperspectral image unmixing is to decompose the mixed pixels into pure spectral features of ground objects and their abundance, so as to more accurately reflect the real situation of the ground objects, which is of great significance for ground object classification and detection, environmental monitoring and agricultural production. Most of the existing unmixing methods learn the abundance matrix by minimizing the reconstruction error based on the square of the F norm or the L2,1 norm based on the preset spectral mixing model, and combining other reasonable abundance matrix constraint regularization terms. However, the reconstruction error objective functions based on the square of the F norm and the L2,1 norm have different defects. For example, when there is high-intensity noise in the matrix, the square of the F norm will seriously deviate from the true value, that is, it is extremely sensitive to noise, which can easily lead to serious limitations on the accuracy of abundance estimation. The defect of the L2,1 norm is that it is not differentiable at zero, which greatly increases the difficulty of algorithm optimization.

[0004] The invention patent with publication number CN117636162A discloses a sparse unmixing method for hyperspectral images, including: obtaining a hyperspectral image, dividing the hyperspectral image to obtain multiple homogeneous regions; performing image regularization on each of the homogeneous regions to obtain a first regularization term; obtaining an abundance matrix based on the global spatial information of the hyperspectral image, and performing low-rank regularization on the abundance matrix to obtain a second regularization term; obtaining a basic sparse unmixing model, obtaining a final sparse unmixing model based on the basic sparse unmixing model, the first regularization term and the second regularization term, and solving the sparse unmixing model. This invention improves the accuracy of the sparse unmixing algorithm, thereby improving the accuracy of the mixed solution results. However, the accuracy of abundance estimation of the current unmixing method is severely limited, and there are deficiencies in the difference learning of spectral similarity in local regions of hyperspectral images, and it is usually difficult to capture accurate spatial correlation. Summary of the invention

[0005] In view of the technical problems that the accuracy of abundance estimation of existing methods is severely limited and it is difficult to capture accurate spatial correlation, the present invention proposes a sparse unmixing method for hyperspectral images based on the σ norm, which adopts image segmentation and homogeneity learning to more accurately obtain the local similarity of hyperspectral images, and adopts the reconstruction error constrained unmixing target based on the σ norm, which can effectively avoid the defects of the F norm and L2,1 mixing norm, and has the characteristics of robustness to noise and easy optimization.

[0006] In order to achieve the above object, the technical solution of the present invention is achieved as follows:

[0007] A sparse unmixing method for hyperspectral images based on σ norm comprises the following steps:

[0008] S1: Obtain a hyperspectral image, select the band of the hyperspectral image using the neighborhood grouping method, and obtain the spectral band range;

[0009] S2: Obtain the USGS spectral library, retain the spectral features in the USGS spectral library that have the same spectral band range as that obtained in step S1, and obtain the ground feature spectral library A;

[0010] S3: Perform local area spectral feature homogeneity learning: divide the hyperspectral image into multiple superpixels, obtain the hyperspectral image matrix Y, and calculate the homogeneity coefficient of each superpixel;

[0011] S4: Based on the superpixel homogeneity coefficient, establish the objective function of the rough mixed solution stage based on the σ norm;

[0012] S5: learning the abundance matrix of the coarse unmixing stage based on the objective function of the coarse unmixing stage;

[0013] S6: Based on the abundance matrix of the rough unmixing stage and the ground object spectral library A, a pruned spectral library B is obtained;

[0014] S7: Perform fine unmixing based on the abundance matrix of the rough unmixing stage and the pruned spectral library B, and output the final abundance matrix.

[0015] The method for band selection of a hyperspectral image using the neighborhood grouping method described in step S1 is: the hyperspectral image is evenly divided into blocks according to the spectral dimension, each block corresponds to a spectral band range, and then the local density and information entropy of the spectral features of each spectral band range are calculated, and the spectral features of each spectral band range with the largest product of local density and information entropy in the hyperspectral image are retained, and the spectral features of the remaining spectral band ranges are eliminated to complete the band selection.

[0016] The method for segmenting the hyperspectral image into multiple superpixels described in step S3 is as follows: set the number of regions K desired to be segmented, initialize the hyperspectral image segmentation using a regular hexagonal grid, use the pixels in the same hexagonal grid as the initial cluster, randomly select the pixel cluster center in each hexagonal grid, set the search area radius to twice the side length of the hexagonal grid, and then search within a range of twice the side length of the hexagonal grid, and calculate the fusion distance D between the pixel cluster center and the pixels in the search range corresponding to the pixel cluster center. j , according to the fusion distance between the pixels in the search range corresponding to the pixel cluster center and the pixel cluster center, the pixels are redistributed to clusters, and the pixels are divided into the categories to which the pixel cluster center closest to them belongs until all cluster centers converge and remain unchanged, completing the image superpixel segmentation, and obtaining the superpixel matrix Y=Y1,Y2,…,Y after the image superpixel segmentation i ,…Y K , where the i-th superpixel matrix n i is the number of pixels in the i-th superpixel matrix.

[0017] The fusion distance D j The calculation process is:

[0018]

[0019] Among them, arccos() represents the inverse cosine function; Represents the i-th pixel cluster center and pixel cluster centers The jth pixel in the corresponding search range The spectral angular distance, Represents the center of pixel cluster Coordinates and pixels The Euclidean distance between coordinates, The pixel cluster center The coordinates of (m j ,n j ) is the pixel The coordinates of the hexagonal grid are shown in Figure 2. p represents the side length of the hexagonal grid, and v is the spectral angular distance. and Euclidean distance The balance parameters of

[0020] The method for calculating the homogeneity coefficient of each superpixel in step S3 is:

[0021]

[0022] Among them, ε is a small positive number.

[0023] The method for establishing the objective function of the rough mixed solution stage based on the σ norm in step S4 is as follows: Establishing the linear mixed model Y of the rough mixed solution stage i =AX i , define the σ norm, and establish the objective function of the rough mixed solution stage based on the σ norm based on the σ norm and the superpixel homogeneity coefficient; the defined σ norm is:

[0024]

[0025] In which, the definition matrix E i =Y i -AX i , L is the matrix E i The number of row vectors, e j Represents the matrix E i The j-th row vector of , σ is the parameter;

[0026] The objective function of the rough mixed solution stage is:

[0027] Among them, X i is the superpixel matrix Y i The corresponding abundance matrix, ||Y i -AX i || σ represents the reconstruction error based on the σ norm, c i ||X i || 2,0 represents the sparse constraint term of the abundance matrix with superpixel homogeneity coefficient as weight, is the constraint abundance matrix X i is a positive indicator function.

[0028] The method for learning the abundance matrix of the coarse demixing stage based on the objective function of the coarse demixing stage described in step S5 is: the optimization solution of the i-th sub-problem of the objective function of the coarse demixing stage is equivalent to the solution of the Lagrangian function, and the abundance matrix corresponding to all superpixel matrices of the coarse demixing stage is obtained according to the Lagrangian function using an iterative optimization method.

[0029] The i-th sub-problem is:

[0030]

[0031] The Lagrangian function is:

[0032]

[0033] Among them, Λ1, Λ2, Λ3 represent Lagrange multipliers, |||| 2,0 represents the number of non-zero row vectors, μ represents the penalty weight, V i and U i is an auxiliary variable.

[0034] The method based on iterative optimization is:

[0035] a: Initialize the abundance matrix Where I represents the unit matrix, initialization Initialize the Lagrange multiplier Λ1=Λ2=Λ3=0;

[0036] b: According to KKT conditions, through formula X i =(A T A+2I) -1 [A T (Y i -E i -Λ1)+V i +U i +Λ2+Λ3]

[0037] Update the abundance matrix X i ;

[0038] c: Calculate the diagonal matrix D and update the matrix E i =μ(2D+μI) -1 (Y i -AX i -Λ1), the diagonal elements of the diagonal matrix D

[0039]

[0040] d: By minimizing Get the auxiliary variable V i ;

[0041] e: Through the formula Ui =max(X i -Λ3,0) Update the calculation auxiliary variable U i ;

[0042] f: Calculate the Lagrange multipliers Λ1, Λ2, Λ3:

[0043] Λ1←Λ1-(Y i -E i -AX i )

[0044] Λ2←Λ2-(X i -V i )

[0045] Λ3←Λ3-(X i -U i )

[0046] The above process is repeated until the loop stop condition is reached to obtain the abundance matrix X corresponding to all superpixel matrices in the rough unmixing stage i .

[0047] The method for obtaining the pruned spectral library B in step S6 is: according to the abundance matrix X=[X1, X2, ..., X K ], calculate the contribution of each spectral feature in the ground object spectral library A in characterizing the mixed pixels of the hyperspectral image x i represents the i-th row vector of the abundance matrix X in the rough unmixing stage, M represents the number of rows of the abundance matrix X in the rough unmixing stage, retains the spectral features with high contribution, and removes the inactive spectral features in the spectral library to obtain the pruned spectral library B.

[0048] The method described in step S7 for performing fine unmixing based on the abundance matrix of the rough unmixing stage and the pruned spectral library B to output the final abundance matrix is: replacing the original ground feature spectral library A with the pruned spectral library B, using the abundance matrix X of the rough unmixing stage obtained in step S5 as the initialization abundance matrix, and cyclically repeating steps b to f of step S5 until the loop stop condition is reached to obtain the final abundance matrix.

[0049] Beneficial effects of the present invention;

[0050] 1) Using image segmentation and homogeneity learning, we can more accurately obtain the local similarity of hyperspectral images and make reasonable adjustments to the intensity of the sparse characteristics of abundance;

[0051] 2) The reconstruction error constrained unmixing objective based on the σ norm is adopted, which can effectively avoid the defects of the F norm and L2,1 mixing norm, and at the same time has the characteristics of quadratic differentiability and emphasis on band-by-band reconstruction.

[0052] 3) When the hyperspectral image is seriously polluted by noise, the present invention can still estimate a clear abundance close to the true value, that is, the present invention has strong robustness to noise. BRIEF DESCRIPTION OF THE DRAWINGS

[0053] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings required for use in the embodiments or the description of the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

[0054] Figure 1 The figure is a flow chart of the method of the present invention.

[0055] Figure 2 Schematic diagram of superpixel homogeneity differences after image segmentation in the present invention.

[0056] Figure 3 It is a relationship diagram of the σ norm, L2,1 norm and F norm square of the present invention.

[0057] Figure 4 This is a comparison chart of the unmixing results on the first simulated hyperspectral image in the present invention.

[0058] Figure 5 It is a schematic diagram for comparing the results of the second simulation experiment in the present invention. DETAILED DESCRIPTION

[0059] The following will be combined with the drawings in the embodiments of the present invention to clearly and completely describe the technical solutions in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0060] like Figure 1 As shown, a sparse unmixing method for hyperspectral images based on σ norm includes:

[0061] S1: Obtain a hyperspectral image, use the neighborhood grouping method to select the band of the hyperspectral image, and obtain the spectral band range.

[0062] Band selection is performed on the hyperspectral image taken by the spectral imager to select the most representative spectral features of the ground objects. In general, the hyperspectral image is divided into l blocks according to the spectral dimension. The number of blocks l is determined by the number of bands L of the hyperspectral image. Where a∈{2,3,4,5,6,7,8}. Each block corresponds to a spectral band range, and then the local density and information entropy of the spectral features of each spectral band range are calculated, the spectral features of each spectral band range with the largest product of local density and information entropy in the hyperspectral image are retained, and the spectral features of the remaining spectral band ranges are eliminated to complete the band selection.

[0063] S2: Obtain the USGS spectral library, retain the spectral features in the USGS (United States Geological Survey) spectral library with the same spectral band range as that obtained in step S1, and obtain a complete ground feature spectral library A corresponding to the range after the band selection in step S1.

[0064] S3: Perform local area spectral feature homogeneity learning: divide the hyperspectral image into multiple superpixels, obtain the hyperspectral image matrix Y, and calculate the homogeneity coefficient of each superpixel.

[0065] The method for segmenting a hyperspectral image into multiple superpixels is as follows: before performing superpixel segmentation on the hyperspectral image, the number of regions K to be segmented is set. K is determined by the size of the hyperspectral image. If the image size is M*N, then Where b∈{20,25,30,35,40}. Use a regular hexagonal grid to initialize the image segmentation, that is, the pixels in the same grid are used as the initial segmentation result, and each local block in the segmentation result is regarded as a cluster; randomly select the pixel cluster center in each grid, set the search area radius to twice the length of the hexagonal grid side, and then search within the range of twice the length of the hexagonal side, and calculate the fusion distance D between the pixel cluster center and the pixels in the search range j ; Fusion distance D j The calculation process is:

[0066]

[0067] Among them, arccos() represents the inverse cosine function; Represents the i-th pixel cluster center and the jth pixel in the cluster search range The spectral angular distance, represents the Euclidean distance between the two coordinates, The pixel that is the center of the pixel cluster The coordinates of (m j ,n j ) is the pixel The coordinates of the hexagonal grid are shown in Figure 2. p represents the side length of the hexagonal grid, which is determined by the size of the hyperspectral image and the number of regions K that are expected to be segmented. v is the spectral angular distance and the Euclidean distance of spatial coordinates The balance parameter is set to v = 0.7 in this example. In each iteration, the fusion distance D between the pixel and the pixel cluster center is calculated. j The pixels are reassigned to clusters and divided into categories to which the cluster centers closest to them belong until all cluster centers converge and remain unchanged, thus completing the image superpixel segmentation. Figure 2 Shows the pixel similarity differences in four superpixels after segmentation of a hyperspectral image of Cuprite, Nevada, USA. Figure 2 It can be seen that the pixels in the two superpixels corresponding to the first row have a high degree of overlap, indicating that the pixels in the two neighborhoods have a high degree of similarity. However, the pixels in the two superpixels corresponding to the second row are significantly different. This also verifies that the spectral homogeneity of different superpixels is significantly different, so it is not appropriate to constrain different superpixel abundance matrices with the same sparsity.

[0068] After segmentation, the hyperspectral image can be represented as a superpixel matrix Y = Y1, Y2, …, Y i ,…Y K , where K is the number of superpixels, is the matrix composed of all pixels in the i-th superpixel, and the number of pixels in this superpixel is n i , Indicates the nth i Pixel.

[0069] Then the homogeneity coefficient of each superpixel is calculated one by one to characterize the homogeneity of the superpixel spectral characteristics, that is, the similarity of all pixels in the superpixel. The superpixel homogeneity coefficient c i The spectral angle distances of all pixels in the superpixel and their pixel cluster centers are closely related, and the calculation method is as follows:

[0070]

[0071] where n i Represents the number of pixels in the ith superpixel. ε is a small positive number to avoid meaningless fractions. In this example, ε is set to 10. -16 (Usually ε<10 -10 ). From the calculation method of the homogeneity coefficient of the superpixel, it can be seen that the homogeneity coefficient of the superpixel represents the similarity of the pixels in the superpixel. i The larger the value, the higher the homogeneity of the pixels in the superpixel. Therefore, it is reasonable to impose a stronger sparsity constraint on the superpixel.

[0072] S4: Based on the superpixel homogeneity coefficient, establish the objective function of the coarse mixed solution stage based on the σ norm.

[0073] Establish the linear mixed model Y in the rough mixing stage i =AX i , define the σ norm, and establish the objective function of the rough mixed solution stage based on the σ norm based on the σ norm and the superpixel homogeneity coefficient; the defined σ norm is:

[0074]

[0075] In which, the definition matrix E i =Y i -AX i , L is the matrix E i The number of row vectors, e j Represents the matrix E i The j-th row vector of , σ is a parameter, and its value range is 0<σ<+∞; From the definition, we can see that the matrix E i The σ norm of is quadratically differentiable, which can meet the requirement of being easy to optimize as a loss function. Figure 3 It shows that when the parameter σ takes different values, ||E i || σ and ||E i || 2,1 , relationship, through Figure 3 It can be seen that as the parameter σ approaches 0 infinitely, the matrix E i The σ norm of ||E i || σ Infinitely close to the matrix E i L 2,1 norm; when the parameter σ approaches positive infinity, the matrix E i The σ norm of ||E i || σ is infinitely close to the matrix E i The square of the F-norm.

[0076] The objective function of the rough mixed solution stage is:

[0077] Among them, X i is the superpixel matrix Y i The corresponding abundance matrix, ||Y i -AX i || σ represents the reconstruction error based on the σ norm, c i ||X i || 2,0 represents the sparse constraint term of the abundance matrix with superpixel homogeneity coefficient as weight, is the constraint abundance matrix X i is a positive indicator function. It is defined as when the abundance matrix X i When the middle element is positive, otherwise Is infinite.

[0078] According to the objective function of the rough solution mixed model, the minimization problem of the objective function is essentially the minimization of the sum of K sub-problems. Obviously, the optimal solution of the objective function is equivalent to the minimization of each sub-problem, so its minimization problem can be achieved by minimizing the K sub-problems separately. This embodiment will next elaborate on the i-th sub-problem Optimization solution steps.

[0079] S5: Learn the abundance matrix of the coarse unmixing stage based on the objective function of the coarse unmixing stage.

[0080] According to the alternating direction multiplication method, let Y i -AX i =E i , introduce auxiliary variable V i and U i The optimization of the subproblem is equivalent to solving the following Lagrangian function:

[0081]

[0082] Where L() is the Lagrangian function, Λ1, Λ2, Λ3 represent Lagrangian multipliers, |||| 2,0 It represents the number of non-zero row vectors, and μ represents the penalty weight.

[0083] The solution of this function can be obtained through an iterative optimization method. The detailed steps are as follows:

[0084] Step a: Initialize the abundance matrix Where I represents the unit matrix, initialization Initialize the Lagrange multipliers Λ1=Λ2=Λ3=0.

[0085] Step b: Abundance matrix X i The optimization is only related to the fourth, fifth and sixth terms of the Lagrangian function. According to the KKT (Karush-Kuhn-Tucker conditions), the abundance matrix X can be updated by i :

[0086] X i =(A T A+2I) -1 [A T (Y i -E i -Λ1)+V i +U i +Λ2+Λ3].

[0087] Step c: Matrix E iThe optimization of is related to the first and fourth terms of the Lagrangian function. This optimization problem can be solved by iteratively calculating the diagonal matrix D and the matrix difference E i The cycle process obtains the matrix E i The diagonal elements of the diagonal matrix D are Matrix E i The update method is E i =μ(2D+μI) -1 (Y i -AX i -Λ1).

[0088] Step d, auxiliary variable V i The optimization of is related to the second and fifth terms of the Lagrangian function, which can be achieved by minimizing This problem can be solved by setting the auxiliary variable V i =X i -Λ2 and when the auxiliary variable V i The L2 norm of a row vector in is less than When all the elements of the row vector are forced to be 0, the auxiliary variable V of the solution of the subproblem is obtained. i In this step, by applying a weight of c i The sparse constraint term of the abundance matrix can be adjusted. This is because when the superpixel homogeneity coefficient c i When it is larger, the superpixel homogeneity is stronger, which means that the corresponding abundance matrix is ​​sparser and the corresponding abundance matrix has more zero row elements. On the contrary, when the superpixel homogeneity coefficient c is i When it is small, it means that the sparsity of its abundance matrix is ​​weak, so the auxiliary variable V i The number of all-zero rows is small.

[0089] Step e: Auxiliary variable U i The optimization of is only related to the third and sixth terms in the Lagrangian function, which can be expressed by the formula U i =max(X i -Λ3,0) Update calculation U i Auxiliary variable U i For positive time, When the function is 0 and the element is negative, the function is infinite. Since the goal is to minimize the function, it is necessary to ensure that the matrix elements are not negative.

[0090] Step f: Calculate the Lagrange multipliers Λ1, Λ2, and Λ3 by the following formula.

[0091] Λ1←Λ1-(Y i -E i -AX i );

[0092] Λ2←Λ2-(Xi -V i );

[0093] Λ3←Λ3-(X i -U i ).

[0094] The above process is repeated until the loop stop condition is reached (in this embodiment, the number of loops is set to 100 times) to obtain the abundance matrix X corresponding to all superpixel matrices in the rough unmixing stage i (i=1,2,...,K).

[0095] S6: Based on the abundance matrix of the rough unmixing stage and the ground object spectral library A, a pruned spectral library B is obtained.

[0096] According to the abundance matrix X = [X1, X2, ..., X K ], calculate the contribution of each spectral feature in the ground object spectral library A in characterizing the mixed pixels of the hyperspectral image The spectral features with high contribution are retained, and the inactive spectral features in the spectral library are removed to obtain the pruned spectral library B.

[0097] In contribution u i In the definition of i represents the i-th row vector of the abundance matrix X in the rough unmixing stage, and M represents the number of rows of the abundance matrix X in the rough unmixing stage. Obviously, the larger the value of the element in the row vector of the abundance matrix in the rough unmixing stage, the greater the contribution of the corresponding spectral feature in characterizing the mixed pixels of the hyperspectral image. In this embodiment, the inactive spectrum refers to the contribution u in the rough unmixing stage. i Spectral features less than 0.01 (usually set to the range [0.01, 0.02]).

[0098] S7: Based on the abundance matrix of the rough unmixing stage and the pruned spectral library B, fine unmixing is performed to output the final abundance matrix.

[0099] The original ground feature spectral library A is replaced by the pruned spectral library B, the abundance matrix X of the rough unmixing stage obtained in step S5 is used as the initialization abundance matrix, and steps b to f of step S5 are repeated cyclically until the loop stop condition is reached to obtain the final abundance matrix.

[0100] In order to verify the effectiveness of the present invention, this embodiment conducted a comparative experiment on two simulated hyperspectral images. The first simulated hyperspectral image is composed of 5 spectra randomly selected from the ground object spectrum library A, and Gaussian noise with a signal-to-noise ratio in the range of [25, 35] dB is added to each band. Figure 4is the unmixing result on the simulated hyperspectral image. The first row in the figure is the real abundance map of the five spectral features in the simulated hyperspectral image; the second and third rows are the abundance maps estimated by the JSpBLRU and PMoSU methods, respectively. It can be seen that the abundance map is full of random noise, and it is difficult to restore the clarity and smoothness of the real abundance map, which shows that the JSpBLRU and PMoSU methods are greatly affected by Gaussian noise. The fourth row, that is, the abundance map estimated by the present invention, shows a clean and smooth characteristic, and the edges of the blocks are clear, and it is little affected by Gaussian noise. Compared with the abundance maps estimated by the JSpBLRU and PMoSU methods, it is closer to the real abundance map. This experimental result comparison verifies the effectiveness of the present invention, especially its strong robustness to Gaussian noise.

[0101] The second simulation experiment is composed of 9 randomly selected spectra in the spectral library, and 20% of the bands are randomly selected on the simulated hyperspectral image to add 10 strip noises. Figure 5 is the unmixing result of this simulation experiment. Similarly, the first row in the figure shows the real abundance value of the hyperspectral image used in this simulation experiment. It can be seen that the distribution of different substances is irregular and the edges of the substances are irregular, which undoubtedly increases the difficulty of unmixing. The second and third rows in the figure are the abundance maps estimated by the JSpBLRU and PMoSU methods, respectively. It can be seen that the unmixing accuracy of these two unmixing methods is seriously affected by the strip noise, especially the clarity of the material edge and the accuracy of the abundance value are seriously reduced. The fourth row shows the abundance map estimated by the present invention. Compared with the estimated abundance maps in the second and third rows, the present invention is very little affected by the strip noise. In particular, the present invention can still clearly restore the real abundance value corresponding to the spectral feature under the pollution of the strip noise, and the estimated abundance value is closer to the real abundance map.

[0102] In summary, the present invention discloses a sparse unmixing method for hyperspectral images based on the σ norm, and mainly describes the necessity of learning the homogeneity of local regions of hyperspectral images and the optimization of reconstruction errors based on the σ norm and its advantages in the unmixing method. The present invention sets the confidence coefficient (i.e., the superpixel homogeneity coefficient) according to the similarity score of the spectral features of the local region, and uses the confidence coefficient to reasonably relax or strengthen the constraint degree of the sparsity of the abundance matrix, which can effectively utilize the spatial correlation of the spectral features; the problem of abundance matrix estimation is processed by using the two stages of coarse unmixing and fine unmixing, and the inactive spectral features in the spectral library are eliminated according to the contribution of the spectrum, which can effectively improve the efficiency of abundance estimation; a new σ norm loss function is set, which has the characteristics of robustness to noise and easy optimization; the present invention is simple to implement and can effectively improve the accuracy of spectral image abundance estimation.

[0103] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principle of the present invention should be included in the protection scope of the present invention.

Claims

1. A sparse unmixing method for hyperspectral images based on σ norm, characterized in that: Includes steps: S1: Obtain a hyperspectral image, select the band of the hyperspectral image using the neighborhood grouping method, and obtain the spectral band range; S2: Obtain the USGS spectral library, retain the spectral features in the USGS spectral library that have the same spectral band range as that obtained in step S1, and obtain the ground feature spectral library A; S3: Perform local area spectral feature homogeneity learning: divide the hyperspectral image into multiple superpixels, obtain the hyperspectral image matrix Y, and calculate the homogeneity coefficient of each superpixel; S4: Based on the superpixel homogeneity coefficient, establish the objective function of the rough mixed solution stage based on the σ norm; S5: learning the abundance matrix of the coarse unmixing stage based on the objective function of the coarse unmixing stage; S6: Based on the abundance matrix of the rough unmixing stage and the ground object spectral library A, a pruned spectral library B is obtained; S7: Perform fine unmixing based on the abundance matrix of the rough unmixing stage and the pruned spectral library B, and output the final abundance matrix.

2. The sparse unmixing method of hyperspectral images based on σ norm according to claim 1 is characterized in that: The method for band selection of a hyperspectral image using the neighborhood grouping method described in step S1 is: the hyperspectral image is evenly divided into blocks according to the spectral dimension, each block corresponds to a spectral band range, and then the local density and information entropy of the spectral features of each spectral band range are calculated, and the spectral features of each spectral band range with the largest product of local density and information entropy in the hyperspectral image are retained, and the spectral features of the remaining spectral band ranges are eliminated to complete the band selection.

3. The sparse unmixing method for hyperspectral images based on σ norm according to claim 2 is characterized in that: The method for segmenting the hyperspectral image into multiple superpixels described in step S3 is as follows: set the number of regions K desired to be segmented, initialize the hyperspectral image segmentation using a regular hexagonal grid, use the pixels in the same hexagonal grid as the initial cluster, randomly select the pixel cluster center in each hexagonal grid, set the search area radius to twice the side length of the hexagonal grid, and then search within a range of twice the side length of the hexagonal grid, and calculate the fusion distance D between the pixel cluster center and the pixels in the search range corresponding to the pixel cluster center. j , according to the fusion distance between the pixels in the search range corresponding to the pixel cluster center and the pixel cluster center, the pixels are redistributed to clusters, and the pixels are divided into the categories to which the pixel cluster center closest to them belongs until all cluster centers converge and remain unchanged, completing the image superpixel segmentation, and obtaining the superpixel matrix Y=Y1,Y2,…,Y after the image superpixel segmentation i ,…Y K , where the i-th superpixel matrix n i is the number of pixels in the i-th superpixel matrix.

4. The sparse unmixing method for hyperspectral images based on σ norm according to claim 3 is characterized in that: The fusion distance D j The calculation process is: Among them, arccos() represents the inverse cosine function; Represents the i-th pixel cluster center and pixel cluster centers The jth pixel in the corresponding search range The spectral angular distance, Represents the center of pixel cluster Coordinates and pixels The Euclidean distance between coordinates, The pixel cluster center The coordinates of (m j ,n j ) is the pixel The coordinates of the hexagonal grid are shown in Figure 2. p represents the side length of the hexagonal grid, and v is the spectral angular distance. and Euclidean distance The balance parameters of The method for calculating the homogeneity coefficient of each superpixel in step S3 is: Among them, ε is a small positive number.

5. The sparse unmixing method for hyperspectral images based on σ norm according to claim 4 is characterized in that: The method for establishing the objective function of the rough mixed solution stage based on the σ norm in step S4 is as follows: Establishing the linear mixed model Y of the rough mixed solution stage i =AX i , define the σ norm, and establish the objective function of the rough mixed solution stage based on the σ norm based on the σ norm and the superpixel homogeneity coefficient; the defined σ norm is: In which, the definition matrix E i =Y i -AX i , L is the matrix E i The number of row vectors, e j Represents the matrix E i The j-th row vector of σ is a parameter; The objective function of the rough mixed solution stage is: Among them, X i is the superpixel matrix Y i The corresponding abundance matrix, ||Y i -AX i || σ represents the reconstruction error based on the σ norm, c i ||X i || 2,0 represents the sparse constraint term of the abundance matrix with superpixel homogeneity coefficient as weight, is the constraint abundance matrix X i is a positive indicator function.

6. The sparse unmixing method for hyperspectral images based on σ norm according to any one of claims 1 to 5, characterized in that: The method for learning the abundance matrix of the coarse demixing stage based on the objective function of the coarse demixing stage described in step S5 is: the optimization solution of the i-th sub-problem of the objective function of the coarse demixing stage is equivalent to the solution of the Lagrangian function, and the abundance matrix corresponding to all superpixel matrices of the coarse demixing stage is obtained according to the Lagrangian function using an iterative optimization method.

7. The sparse unmixing method for hyperspectral images based on σ norm according to claim 6 is characterized in that: The i-th sub-problem is: The Lagrangian function is: Among them, Λ1, Λ2, Λ3 represent Lagrange multipliers, |||| 2,0 represents the number of non-zero row vectors, μ represents the penalty weight, V i and U i is an auxiliary variable.

8. The sparse unmixing method for hyperspectral images based on σ norm according to claim 7, characterized in that: The method based on iterative optimization is: a: Initialize the abundance matrix Where I represents the unit matrix, initialization Initialize the Lagrange multiplier Λ1=Λ2=Λ3=0; b: According to KKT conditions, through formula X i =(A T A+2I) -1 [A T (Y i -E i -Λ1)+V i +U i +Λ2+Λ3] Update the abundance matrix X i ; c : Calculate the diagonal matrix D and update the matrix E i =μ(2D+μI) -1 (Y i -AX i -Λ1), the diagonal elements of the diagonal matrix D d: By minimizing Get the auxiliary variable V i ; e: Through the formula U i =max(X i -Λ3,0) Update the calculation auxiliary variable U i ; f: Calculate the Lagrange multipliers Λ1, Λ2, Λ3: Λ1←Λ1-(Y i -E i -AX i ) Λ2←Λ2-(X i -V i ) Λ3←Λ3-(X i -U i ) The above process is repeated until the loop stop condition is reached to obtain the abundance matrix X corresponding to all superpixel matrices in the rough unmixing stage i .

9. The sparse unmixing method for hyperspectral images based on σ norm according to claim 8, characterized in that: The method for obtaining the pruned spectral library B in step S6 is: according to the abundance matrix X=[X1, X2, ..., X K ], calculate the contribution of each spectral feature in the ground object spectral library A in characterizing the mixed pixels of the hyperspectral image x i represents the i-th row vector of the abundance matrix X in the rough unmixing stage, M represents the number of rows of the abundance matrix X in the rough unmixing stage, retains the spectral features with high contribution, and removes the inactive spectral features in the spectral library to obtain the pruned spectral library B.

10. The sparse unmixing method for hyperspectral images based on σ norm according to claim 9, characterized in that: The method described in step S7 for performing fine unmixing based on the abundance matrix of the rough unmixing stage and the pruned spectral library B to output the final abundance matrix is: replacing the original ground feature spectral library A with the pruned spectral library B, using the abundance matrix X of the rough unmixing stage obtained in step S5 as the initialization abundance matrix, and cyclically repeating steps b to f of step S5 until the loop stop condition is reached to obtain the final abundance matrix.

Images