PaperlayerView exact revision on arXiv

Source-linked AI summary

Monotone Piecewise-Linear Harmonization of Multitemporal Multispectral Surface Reflectance With Unsupervised Deep Multiview Pseudoinvariant Features Learning

Lino Garda Denaro, Chia-Hsiang Lin

arXiv:2610.21183v1eess.IVcs.CV

TL;DR

Residual interdate reflectance biases can distort multitemporal analysis, motivating a method that identifies deep pseudoinvariant features from satellite observations. DeepMGCCA combines multiview learning with robust reweighting and reconstruction, and a quantile-derived monotone piecewise-linear transform achieved the strongest reported harmonization results.

  • Problem

    Residual biases between acquisition dates can make reflectance differences less reliable for multitemporal analysis, while linear methods do not capture nonlinear spatial–spectral relationships.

    Source: §I · pp. 1–2

  • Method

    DeepMGCCA learns deep pseudoinvariant features using a shared convolutional encoder, robust reweighting, and reconstruction regularization, then applies a monotone piecewise-linear transform fitted on DIF quantiles.

    Source: §I · p. 2, 2) Deep learning approach · p. 5, §IV · p. 12

  • Results

    The proposed transform achieved the best RMSE, PSNR, and SSIM in every ROI and the best between-date spectral angle and Pearson’s ρ in three of four, lowering ROI-mean RMSE by 32%.

    Source: §IV · p. 13

  • Takeaways & Limitations

    The findings support deep invariant-feature learning for multitemporal surface-reflectance bias removal while preserving temporal and spectral information.

    Source: §I · p. 3, §IV · p. 13

Abstract

from arXiv · show

Surface reflectance products derived from optical satellite may exhibit systematic inconsistencies across acquisition dates, because atmospheric correction is generally performed independently for each image. Such residual biases distort the temporal behavior of invariant ground targets and reduce the reliability of multitemporal surface reflectance analysis. This study proposes Deep Multitemporal Generalized Canonical Correlation Analysis (DeepMGCCA), an unsupervised deep multiview pseudoinvariant feature learning framework that integrates Deep Generalized Canonical Correlation Analysis (GCCA) with Multivariate Alteration Detection (MAD) to identify deep-pseudoinvariant features (DIFs) for the removal of surface reflectance biases across multitemporal satellite observations. Unlike conventional MAD-based approaches that rely on linear GCCA, the proposed framework replaces date-specific linear projections with a shared convolutional encoder. The framework learns deep invariant features (DIFs) directly from pixel-weighted canonical variates via iterative reweighting strategy, while a reconstruction regularizer preserves surface reflectance information and prevents degenerate collapse from suppressing genuine temporal variability. Experiments were conducted on four annual Landsat-8/9 and Sentinel-2 acquisitions across four regions of interest (ROIs) to evaluate four DIF-based harmonization transforms. A monotone piecewise-linear transform derived from DIFs quantiles consistently achieved the largest bias reduction, lowering the between-date spectral angle by up to 2.54 and RMSE by 0.015, and raising PSNR by up to 4.04~dB. These findings demonstrate the potential of deep invariant feature learning for multitemporal surface reflectance bias removal while highlighting the importance of preserving physically meaningful temporal and spectral information during nonlinear representation learning.

I. INTRODUCTION

The paper addresses residual interdate reflectance bias and the need for statistically principled pseudoinvariant-feature selection. It introduces DeepMGCCA and a globally fitted monotone piecewise-linear transform for learning invariant features and harmonizing multitemporal imagery.

I. INTRODUCTION

  • Residual interdate bias remains a multitemporal harmonization problem even after cross-sensor harmonization, while normalization reliability depends on selecting reliable PIFs.PIF extraction is therefore an independently important problem for multitemporal reflectance analysis.Source: §I · p. 1
  • Linear GCCA and MAD provide ordered invariance-related components, but fixed pixelwise projections cannot exploit spatial context; kernel methods face scene-scale memory limits.Convolutional encoders offer a way to learn nonlinear spatial–spectral relationships without explicit pairwise kernel storage.Source: §I · p. 2
  • DeepMGCCA uses a shared nonlinear encoder with date-specific heads, while reconstruction regularization and robust weighting address collapse that could erase change-sensitive information.The resulting features support pseudoinvariant-feature extraction for change detection and harmonization.Source: §I · pp. 2–3

II. METHODOLOGY

  • The paper frames nonlinear DIF-based transformation as a way to harmonize while preserving spectral signatures needed for downstream identification and classification.The introduction also describes linear GCCA–MAD as the methodological backbone and DeepMGCCA as its deep-learning instantiation.Source: §I · p. 3, §II · p. 3

I. INTRODUCTION

  • The work contributes a globally fitted monotone piecewise-linear transform estimated from invariant-pixel quantiles across date pairs, alongside DeepMGCCA for unsupervised deep invariant-feature learning.The transform is designed to preserve temporal and spectral consistency, with curvature regularization and no privileged reference date.Source: §I · p. 3

A. Linear GCCA approach

The linear GCCA–MAD formulation identifies invariant features by jointly relating multiple image dates under linear dependence. Its difference variates yield a normalized change statistic that can be thresholded to select PIFs.

A. Linear GCCA approach

  • Linear GCCA characterizes multiple observations under linear dependence and identifies invariant features through a joint constrained optimization.For connected dates, maximizing covariance is equivalent to minimizing the variance of their difference variate, which is the MAD variate.Source: A. Linear GCCA approach · p. 3
  • The NMAD statistic sums squared standardized MAD components, with the number of spectral bands defining the chi-square degrees of freedom.The statistic is thresholded to separate PIFs from non-PIFs.Source: A. Linear GCCA approach · pp. 3–4
  • Pixel weights replace unweighted population moments in the GCCA criterion to down-weight observations likely to have changed.The iterative procedure updates weights and repeats until convergence.Source: A. Linear GCCA approach · p. 4

B. Nonlinear GCCA approach

Nonlinear GCCA extends the linear formulation to higher-order dependencies using kernel-based or deep-learning strategies. Kernel GCCA uses implicit feature maps but incurs scaling costs that motivate the deep alternative.

B. Nonlinear GCCA approach

  • Nonlinear GCCA aims to capture higher-order dependencies among dates, and the paper investigates kernel-based and deep-learning-based strategies.The kernel approach replaces linear projections with implicit nonlinear feature-space projections represented through Gram matrices.Source: B. Nonlinear GCCA approach · p. 4, 1) Kernel approach · p. 4

1) Kernel approach:

  • Kernel GCCA models nonlinear correlations while retaining the generalized eigenvalue structure without explicitly evaluating the feature map.A ridge term is added to regularize the kernel-matrix inversion.Source: 1) Kernel approach · p. 4
  • Kernel matrices scale as n × n, making kernel GCCA impractical for full-scene imagery with millions of pixels.The paper presents kernel GCCA primarily as a comparison baseline rather than a scalable solution.Source: 1) Kernel approach · p. 4

B. Nonlinear GCCA approach

  • The deep-learning alternative replaces fixed linear projections with learned nonlinear representations from a shared encoder and date-specific heads.The overview connects those representations to consensus-aligned variates and the downstream harmonization pipeline.Source: 1) Kernel approach · p. 4, Figure 1 · p. 4

2) Deep learning approach:

DeepMGCCA learns date-specific feature outputs from a shared convolutional encoder, aligns them through consensus, and uses correlation and reconstruction objectives. Its weighted iterative workflow produces thresholded deep-pseudoinvariant features.

2) Deep learning approach:

  • A shared convolutional encoder learns spatial–spectral structure across dates, while date-specific heads model residual date behavior in r-dimensional representations.The representations are flattened, weighted, and mean-centered during training.Source: 2) Deep learning approach · p. 5
  • Weighted center-scaling enforces the unit-variance constraint for each date’s head outputs during training.The resulting representations are then consensus-aligned to form canonical variates.Source: 2) Deep learning approach · pp. 5–6
  • The correlation objective is paired with a reconstruction regularizer so learned outputs remain decodable to their inputs rather than collapsing to trivial agreement.The reconstruction error sums squared differences between each input and its decoder reconstruction, weighted by λrecon.Source: 2) Deep learning approach · p. 5
  • Consensus-aligned variates are differenced across dates, aggregated and standardized into NMAD scores, then thresholded to produce the final DIF image.The stated training procedure repeats pixel reweighting until convergence and outputs the variates and DIFs.Source: 2) Deep learning approach · p. 6
  • Initial pixel weights combine a cloud-robust temporal reference with latent-space spectral-angle factors, and later updates use the NMAD rule on trained canonical variates.The initial temporal reference uses the lower of the mean and median across dates; the per-date latent spectral-angle score aggregates connected partners.Source: 2) Deep learning approach · p. 6

C. Image harmonization

Each date and band receives a monotone piecewise-linear transform fitted to pooled-DIF quantiles, with agreement, identity anchoring, spectral consistency, and curvature terms shaping the fit. The resulting convex, bound-constrained least-squares problem applies the fitted transforms to every pixel.

C. Image harmonization

  • Each date-band transform uses pooled-DIF quantile knots and a hinge basis, replacing a two-parameter affine fit with a (K+1)-parameter piecewise-linear transform.The implementation uses K = 3 knots; the first and last ramps extrapolate linearly beyond the DIF range.Source: C. Image harmonization · pp. 6–7
  • Monotonicity is enforced through coefficient bounds that impose a minimum segment slope, preventing flat segments from collapsing bright pixels to one value.Pixels outside the DIF range are linearly extrapolated using the end slopes.Source: C. Image harmonization · p. 7
  • The transforms are fitted to DIF-population quantiles, and their agreement term aligns transformed quantile functions across dates without privileging a reference date.The identity anchor prevents the dates from shifting jointly to an arbitrary common level.Source: C. Image harmonization · p. 7
  • The spectral term couples hinge coefficients across adjacent bands, while the curvature term penalizes coefficient second differences across knots to preserve spectral signatures and smooth transforms.The curvature penalty has no counterpart in the affine formulation.Source: C. Image harmonization · p. 7
  • The four terms form a convex bound-constrained weighted least-squares problem because each transform is linear in its coefficients, despite being nonlinear in reflectance.The fitted transform is applied to every pixel in each date-band image, not only to DIF pixels.Source: C. Image harmonization · p. 7

A. Datasets

The study evaluates four ROIs spanning urban, forest, desert, and coastal settings. It uses six spectral bands common to Landsat 8/9 OLI and Sentinel-2 MSI.

A. Datasets

  • The four study ROIs are Taipei, Java, Turkey, and Saudi Arabia, selected to represent varied environmental and landscape conditions.Taipei is predominantly urban, Java forest-dominated, Turkey combines urban development and a strait, and Saudi Arabia is arid with urban characteristics.Source: A. Datasets · pp. 7–8
  • The analysis is restricted to blue, green, red, NIR, SWIR-1, and SWIR-2, the six bands with near-equivalent counterparts on both sensors.The common-band restriction supports cross-sensor analysis and comparability with the Harmonized Landsat–Sentinel-2 product and prior studies.Source: A. Datasets · p. 8

B. Experimental result and comparison

Across four ROIs, the DeepMGCCA-DIF transform removes visible interdate seams and reduces profile shifts while retaining spatial detail and genuine variation. The evaluation uses a shared residual encoder and fits the harmonization transform on DIF quantiles.

B. Experimental result and comparison

  • The three-block residual shared encoder has nonlinear per-date heads, and training runs for ten epochs under iterative reweighting.The harmonization transform is fitted on the global DIF quantile population and applied to every scene pixel.Source: B. Experimental result and comparison · p. 8
  • The harmonized mosaics remove visible interdate brightness and color seams while preserving fine spatial texture such as roads, buildings, and field boundaries.The same ground features retain their color and brightness across date boundaries.Source: B. Experimental result and comparison · p. 8, Figure 4 · p. 9
  • Before harmonization, per-band reflectance profiles shift between date segments; afterward, the shifts are strongly reduced and the four dates overlap along the row axis.The profiles indicate reduced interdate drift on invariant terrain without suppressing genuine spatial or temporal variation.Source: B. Experimental result and comparison · p. 9
  • In the Turkey ROI, the lower-row profile decline reflects the strait’s low reflectance, and the sea-dominated segment was least consistent before harmonization.The reported qualitative results indicate that the DIF-guided transform removes interdate bias without flattening spatial texture or altering genuinely changed surfaces.Source: B. Experimental result and comparison · p. 9

C. Performance evaluation

Across four ROIs, the proposed transform improves several harmonization metrics, while trading slightly higher self-distortion for lower between-date disagreement. Profile comparisons and a Taipei exception qualify the overall pattern.

C. Performance evaluation

  • Across four ROI means, the proposed transform lowers RMSE from 0.0295 to 0.0200 (−32%) and SAMh from 5.30° to 3.46° (−35%), while raising PSNR from 32.60 to 35.64 dB and SSIM from 0.8310 to 0.8693.Table III compares variants fitted on the same DIF population and scored on identical PIF masks.Source: C. Performance evaluation · p. 9, Table III · p. 11
  • The proposed transform achieved the best RMSE, PSNR, and SSIM in all four ROIs, while GSC had lower self-distortion and the proposed transform lower SAMh.The comparison is a systematic trade-off: the proposed transform’s averaged SAMh is 3.46° versus 3.96° for GSC, while GSC’s averaged self-distortion is 2.06°.Source: C. Performance evaluation · pp. 9–10
  • The proposed transform’s higher SAMsc does not indicate band misalignment: per-band profiles remain aligned, and its strictly positive minimum slope preserves spectral ordering.The authors characterize GSC and the proposed transform as defensible points on a Pareto front between spectral self-preservation and greater interdate bias removal.Source: C. Performance evaluation · p. 10, Figure 5 · p. 10
  • In Taipei, reference-based fits minimized SAMh, and the proposed transform’s Pearson’s ρ fell from 0.7416 to 0.7106, unlike the other methods.The passage attributes this pattern to high-frequency variation in Taipei’s 10 m pansharpened data and notes that Pearson’s ρ captures only linear association.Source: C. Performance evaluation · p. 10
  • Reference-affine and sequential stepwise methods align dates to a reference, so fit quality decays with radiometric distance; sequential errors can accumulate until correction is worse than none.The reported examples are 5.24° self-distortion for RA in Turkey and SSA PSNR of 33.72 dB versus 33.74 dB unprocessed in Java.Source: C. Performance evaluation · pp. 10–11

1) DeepMGCCA architecture and training:

DeepMGCCA uses shared nonlinear features and date-specific heads to learn cross-date canonical representations without supervision. Reconstruction regularization and robust weighting are used to keep the features informative and the invariant-pixel selection reliable.

1) DeepMGCCA architecture and training:

  • A shared encoder learns features common to all dates, while date-specific 1×1 heads model residual behavior without spatial reach, preserving the GCCA structure and requiring agreement to arise from shared features.The correlation itself supplies the learning signal, so the pipeline is fully unsupervised.Source: 1) DeepMGCCA architecture and training · p. 11
  • A per-date decoder reconstructs each input from its features to counter collapse; λrecon = 1 balances reconstruction and correlation terms, with DIF quality reported as insensitive over [0.5, 2].At zero reconstruction weight collapse is unchecked, while very large values make the model an autoencoder.Source: 1) DeepMGCCA architecture and training · p. 11

2) Objective, weighting, and reliability:

  • The maximum-variance GCCA criterion is solved through a closed-form eigenproblem, while the agreement rows, identity anchors, band smoothness, curvature penalties, and positive-slope bounds form the transform fit.The bound-constrained least-squares problem is linear in its coefficients and therefore convex; fitted transforms are applied to every pixel.Source: 2) Objective, weighting, and reliability · pp. 11–12
  • IR-MAD weights define the training covariances, constraint, and detection statistic; objective and reflectance-RMSE ratios between weight deciles serve as reliability diagnostics.The weights are seeded with temporal and correlation-scaled spectral-angle factors, then updated with the NMAD rule.Source: 2) Objective, weighting, and reliability · p. 12
  • The eigenproblem is independent of pixel count, with only per-date covariances growing at O(Nr^2), unlike kernel GCCA’s O(n^2) memory and O(n^3) time.The authors state this permits million-pixel scenes to be processed at modest cost.Source: 2) Objective, weighting, and reliability · p. 12

3) Thresholding and harmonization:

The study combines DeepMGCCA with a constrained monotone piecewise-linear transform fitted on DIF quantiles. Across four environmental ROIs, the transform improves key harmonization metrics, while the authors identify jointly optimizing the deep learner and transform as future work.

3) Thresholding and harmonization:

  • The DIF population is selected by thresholding standardized canonical differences against a chi-square null; the fixed tail cutoff adapts to each scene’s score distribution.The cutoff admits a wider no-change region in noisy scenes and a narrower one in stable scenes; ties are broken by a seeded random draw.Source: 3) Thresholding and harmonization · p. 12
  • The transform fits a globally aligned common level from DIF quantiles using a convex, bound-constrained least-squares objective, with monotonicity enforced by a positive minimum slope.Its agreement, identity-anchor, band-smoothness, and knot-curvature terms constrain the fit without privileging a reference date.Source: 3) Thresholding and harmonization · p. 12

IV. CONCLUSIONS AND FUTURE WORKS

  • DeepMGCCA learns DIFs from pixel-weighted canonical variates using a shared encoder and date-specific heads, while reconstruction regularization guards against collapse.The conclusion proposes future work on jointly optimizing the convex transform and deep invariant-feature learner.Source: §IV · pp. 12–13
  • Across four ROIs, the proposed transform achieved the best RMSE, PSNR, and SSIM in every ROI, and the best between-date spectral angle and Pearson’s ρ in three of four.ROI-mean RMSE fell by 32% and between-date spectral angle by 35%; the transform accepts slightly higher self-distortion than GSC.Source: §IV · p. 13

Related papers

Loading 2610.21183v1…