Source-linked AI summary
Nonlinear hyperspectral unmixing with robust nonnegative matrix factorization
Cédric Févotte, Nicolas Dobigeon
TL;DR
Hyperspectral unmixing often relies on linear mixing, although some observations contain nonlinear effects that standard models may not capture. The paper introduces robust LMM and its robust NMF formulation, using sparse additive residuals, beta-divergence, and block-coordinate optimization. Synthetic and real-data simulations show that the approach competes with state-of-the-art linear and nonlinear unmixing methods.
Problem
Linear mixing can be inaccurate for hyperspectral applications containing nonlinear effects, while existing nonlinear models require specifying a particular nonlinearity.
Method
The robust LMM adds a nonnegative, group-sparse residual to linear mixing and estimates endmembers, abundances, and residuals with beta-divergence and block-coordinate descent using majorization-minimization updates.
Results
Simulations on synthetic and real data show that the rLMM-based unmixing technique mitigates several nonlinear effects while preserving good estimation performance on linear mixtures.
Takeaways & Limitations
The proposed model captures nonlinear effects as sparse additive outliers without requiring a particular analytical nonlinearity model.
Abstract
from arXiv · showhide
This paper introduces a robust mixing model to describe hyperspectral data resulting from the mixture of several pure spectral signatures. This new model not only generalizes the commonly used linear mixing model, but also allows for possible nonlinear effects to be easily handled, relying on mild assumptions regarding these nonlinearities. The standard nonnegativity and sum-to-one constraints inherent to spectral unmixing are coupled with a group-sparse constraint imposed on the nonlinearity component. This results in a new form of robust nonnegative matrix factorization. The data fidelity term is expressed as a beta-divergence, a continuous family of dissimilarity measures that takes the squared Euclidean distance and the generalized Kullback-Leibler divergence as special cases. The penalized objective is minimized with a block-coordinate descent that involves majorization-minimization updates. Simulation results obtained on synthetic and real data show that the proposed strategy competes with state-of-the-art linear and nonlinear unmixing methods.
I. INTRODUCTION
The paper proposes robust LMM, which augments linear spectral mixing with a nonnegative, group-sparse residual for pixels exhibiting nonlinear effects. This yields a constrained robust NMF formulation that uses beta-divergence and is optimized by block-coordinate updates.
- Motivation: Spectral unmixing decomposes hyperspectral observations into endmember spectra and their pixel-wise abundance proportions.The conventional linear mixing model approximates observations as Y ≈ MA, but may be inaccurate for applications involving nonlinear interactions.
- Proposed model: The robust LMM models observations as Y ≈ MA + R, where nonnegative residual columns represent pixels affected by nonlinearities.Unlike specified nonlinear models, the residual does not require an analytical form for the nonlinearity.
- Proposed model: Group sparsity assumes most residual columns are zero, so the linear mixing model applies to most pixels while nonlinear effects remain localized.The residual energy is penalized at the group level, treating an entire pixel residual as active or inactive.
- Extensions: The paper extends its preliminary conference work by generalizing squared Euclidean fitting to beta-divergence and deriving multiplicative updates through majorization-minimization.It also introduces a rule of thumb for selecting the penalty weight and reports expanded synthetic and real-data experiments.
- Proposed model: The formulation preserves nonnegativity for Y, M, A, and R while constraining abundance coefficients to sum to one.These constraints follow standard assumptions used in hyperspectral data models.
B. Objective function
The objective combines beta-divergence data fidelity with a group-sparse penalty on the nonnegative residual, under nonnegativity and abundance sum-to-one constraints. This defines a robust NMF formulation that generalizes earlier squared-Euclidean robust NMF models.
- Objective function: The optimization minimizes J(M, A, R) = D(Y|MA + R) + λ∥R∥2,1 over the endmembers, abundances, and residual.The penalty weight λ controls the trade-off between reconstruction fidelity and residual sparsity.
- Constraints: The robust NMF problem imposes nonnegativity on the abundance coefficients and retains the abundance sum-to-one constraint.The residual is also nonnegative, consistent with the paper’s robust mixing model.
- Novelty: The formulation is novel in combining a nonnegative residual, an ℓ2,1 penalty, abundance sum-to-one constraints, and beta-divergence fitting.Earlier robust NMF formulations considered squared Euclidean distance rather than the more general beta-divergence.
- Data fidelity: Beta-divergence provides a continuous fitting family whose parameter β governs the assumed observation-noise statistics.Its special cases include Itakura-Saito divergence at β = 0, Kullback-Leibler divergence at β = 1, and squared Euclidean distance at β = 2.
III. BLOCK-COORDINATE DESCENT ALGORITHM
The rNMF objective is minimized by cyclically updating M, A, and R, with MM-derived multiplicative updates for M and R. The algorithm preserves nonnegativity and has linear per-iteration complexity.
- Overall algorithm: The block-coordinate descent updates M, A, and R in turn while conditionally decreasing the nonconvex objective.Because the objective is nonconvex, initialization can affect the returned local solution.
- Parameter updates: The M and R updates are obtained through majorization-minimization, whereas the A update uses a heuristic scheme for the sum-to-one constrained abundances.The abundance constraint is handled through a change of variables because the direct Lagrange-multiplier approach generally fails except at β = 1.
- Implementation: Multiplicative updates preserve nonnegativity automatically and require linear complexity O(LKP) per iteration.The implementation can use γ(β) = 1 in practice while still reducing the objective and accelerating convergence.
- Update of M: Majorization-minimization constructs a tight upper bound and minimizes it to produce a valid descent update for M.The bound satisfies equality at the current iterate, so each MM step does not increase the objective.
- Update of M: The β-divergence is decomposed into convex and concave parts, enabling Jensen and tangent inequalities to form separable auxiliary functions.The auxiliary function breaks the sum over endmembers, making optimization separable with respect to individual m_lk entries.
B. Update of the outlier term R
The outlier term R captures nonlinear effects as nonnegative, group-sparse deviations from the low-rank component MA. Its update majorizes both the β-divergence data fit and the group-sparsity penalty.
- Outlier formulation: R is estimated by minimizing D(Y|MA + R) + λ∥R∥2,1 subject to R ≥ 0.The ℓ2,1 penalty promotes group sparsity in the outlier component.
- Majorization: The data-fitting term is majorized using Jensen and tangent inequalities applied to the convex-concave decomposition of the β-divergence.The low-rank component is [MA], while the current approximation is formed by adding the current R iterate.
- Majorization: The group-sparsity penalty is itself majorized with a quadratic tight upper bound, decoupling spectral bands within the square root.Additional monomial majorization matches the leading degree of the data-fit term, with separate constructions for β > 2 and β ≤ 2.
- Update properties: The resulting multiplicative update preserves nonnegativity when initialized nonnegative and was observed to decrease the objective for the tested β values.The update direction follows the sign of the gradient, moving the variable toward lower objective values.
- Algorithm: The overall procedure is implemented in matrix form as Algorithm 1, which updates the outlier term within group robust NMF.The algorithm initializes M, A, and R, updates parameters iteratively, and checks relative objective decrease for convergence.
D. Setting the value of λ
The penalty weight λ is selected by a method-of-moments rule that matches empirical and prior expectations. The resulting estimate is a practical heuristic without statistical guarantees.
- Role of λ: λ controls the trade-off between the β-divergence data fit D(Y|MA + R) and the group-sparsity penalty ∥R∥2,1.Choosing its value is difficult because it governs the balance between reconstruction fidelity and sparsity.
- Selection rule: The proposed rule of thumb interprets the objective as a joint likelihood and matches the empirical data mean with its prior expectation.The prior model assumes E[[MA]lp] = ρ, while the empirical expectation is ˆµ = (LP)^−1 P ylp.
- Caveat: The estimate is a gross practical approximation with no statistical guarantee and strong dependence on the prior expectation ρ.Because ρ is lower bounded by zero, the estimate has lower bound λ0 = C/ˆµ.
- Evaluation choice: The evaluations use λ = λ0, which the authors report produced satisfactory results.This choice corresponds to a plausible minimum degree of sparsity under the stated construction.
IV. EXPERIMENTS WITH SYNTHETIC DATA
The paper evaluates the relevance of the proposed rLMM and the accuracy of its rNMF algorithm using synthetic data.
- Experimental scope: Synthetic-data experiments evaluate both the proposed rLMM and the accuracy of its corresponding rNMF algorithm.
A. Data generation
Synthetic hyperspectral images were generated from linear and nonlinear mixing scenarios, with controlled abundances, nonlinear interactions, and Gaussian noise. A second experiment tested robustness when pure pixels were absent.
- Data generation: Four 64 × 64 images with K = 3 or 6 pure spectral components represented one linear and three nonlinear mixing scenarios.The nonlinear images contained LMM pixels alongside pixels generated with nonlinear component interactions.
- Data generation: Three fourths of pixels followed the LMM, while the remaining fourth followed a nonlinear interaction model.This construction made nonlinear effects sparse at the pixel level.
- Data generation: The nonlinear scenarios used the Nascimento, Fan bilinear, and generalized bilinear models.The generalized bilinear model included interaction coefficients γijp controlling pairwise endmember interactions.
- Data generation: Abundances were sampled uniformly from admissible abundance sets, with a separate experiment truncating those sets to remove pure pixels.The second experiment evaluated robustness to the absence of pure pixels.
- Data generation: Generalized bilinear interaction coefficients were uniformly sampled from (0, 1), and white Gaussian noise was added at SNR = 30dB.The selected noise level was described as admissible for most real imaging spectrometers.
B. Compared methods
The comparison combined rNMF with standard and nonlinear model-specific baselines. Unlike the two-step baselines, rNMF jointly estimated endmembers, abundances, and the nonlinear component.
- Compared methods: rNMF was compared with state-of-the-art algorithms designed for the corresponding linear and nonlinear mixing models.The baselines first estimated endmembers and then estimated abundances using a separate inversion step.
- Compared methods: VCA followed by FCLS provided the standard baseline for linear mixing models.Both methods are standard endmember extraction and inversion techniques for the LMM.
- Compared methods: Heylen’s endmember extraction method was paired with NM, FM, GBM, and polynomial post-nonlinear inversion methods.These inversion methods targeted different nonlinear mixing models.
- Compared methods: rNMF used β = 2 and λ = λ0, with initialization from either VCA or Heylen’s algorithm.Optimization stopped when the relative objective change fell below 10−5.
C. Performance measures
Unmixing performance was assessed separately for endmember estimation and abundance estimation. The evaluation used average spectral angle mapper and global mean square error, respectively.
- Performance measures: Endmember estimation accuracy was measured using the average spectral angle mapper (aSAM).The metric evaluates spectral similarity between estimated and reference endmembers.
- Performance measures: Abundance estimation accuracy was measured using the global mean square error (GMSE).The metric evaluates errors between estimated and reference abundance coefficients.
D. Results and discussion
Experiments on synthetic and real hyperspectral data evaluated rNMF against competing methods and examined β-divergence choices for missing-data reconstruction. Results showed competitive unmixing and best reconstruction at β = 1 or 1.5, depending on the image.
- D. Results and discussion: Synthetic results showed that rNMF generally improved endmember accuracy from VCA or Heylen initialization, with or without pure pixels.The same experiments found flexible abundance estimation across multiple mixing scenarios.
- D. Results and discussion: The rLMM-based method mitigated several nonlinear effects while preserving good estimation performance on linear mixtures.This conclusion was based on the synthetic evaluation across linear and nonlinear scenarios.
- D. Results and discussion: The real-data study used 50×50-pixel Moffett and Madonna sub-images with partial ground truth.Moffett contains lake and coast regions, while Madonna is mainly composed of oak and chestnut trees plus an additional unknown endmember.
- D. Results and discussion: Table II reports endmember aSAM and abundance GMSE2, with best and second-best scores identified for rNMF initialized by VCA or Heylen’s method.The table compares estimation performance across the evaluated methods.
- D. Results and discussion: Missing-data reconstruction removed 25, 50, or 75% of pixels and fitted MA to observed entries while omitting R because missing-entry outliers were not identifiable.The reconstructed missing pixels were evaluated using aSAM over 10 runs for β values from −1 to 3.
- D. Results and discussion: β mattered for reconstruction: β = 1 was best for Moffett and β = 1.5 for Madonna, while β in [0, 2] gave sensibly similar performance.The figure reports average SAM values with standard deviations over 10 runs and observation rates from 25% to 75%.
C. Robust unmixing results
The proposed rLMM extends linear unmixing with a sparse additive residual for nonlinear effects and is evaluated on two real hyperspectral images. The results align with prior work while revealing pixels not explained by the standard LMM.
- On the Moffett and Madonna images, β = 1 and β = 2 produced visually similar abundance maps and endmember estimates.
- Residual-energy maps identified mostly LMM-consistent pixels in Moffett, with outliers near the lake shore, and residual patterns near oak trees, an unknown endmember, and sensor artifacts in Madonna.
- The reported real-data unmixing results agree with previous results while additionally identifying pixels unexplained by the standard LMM.
- The proposed rLMM extends the standard LMM with a residual term that captures nonlinear effects as sparsely active additive outliers.
- Unlike approaches requiring a specified nonlinearity, rLMM treats nonlinear effects as outliers without imposing a particular analytical model.
- The resulting robust NMF problem is minimized using a block-coordinate descent algorithm with multiplicative updates.