Upgrade to Pro — share decks privately, control downloads, hide ads and more …

Frédéric Schmidt

S³ Seminar
September 25, 2020

Frédéric Schmidt

(Géosciences Paris-Saclay, Université Paris-Saclay)

Title — Démélange spectral sur Mars

Abstract — Le démélange linéaire spectral supervisé et non-supervisé a fait l’objet de nombreuses études en mathématique appliquée. Parmi les récentes avancées, la prise en compte de la positivité et de la parcimonie a suscité un engouement important dans ce domaine. La présentation détaillera plusieurs types d’applications sur des problématiques de Planétologie pour l’étude de la surface et de l’atmosphère de Mars, utilisant des données spatiales (spectroscopie, spectro-imagerie). La dernière partie insistera sur les problèmes ouverts du domaine.

Biography — Frédéric Schmidt est actuellement Professeur au Département des Sciences de la Terre d’Orsay (Université Paris-Saclay) et responsable du Master inter-établissement « Planétologie et Exploration Spatiale ». Il est membre du laboratoire Géosciences Paris-Saclay (GEOPS) depuis 2009. Ses recherches portent sur la thématique des glaces et des interactions surface-atmosphère des corps planétaires. Il est scientifique associé à plusieurs instruments embarqués sur les missions spatiales Mars Express et ExoMars Trace Gas Orbiter de l’Agence Spatiale Européenne.

S³ Seminar

September 25, 2020
Tweet

More Decks by S³ Seminar

Other Decks in Research

Transcript

  1. Spectral unmixing on Mars
    GEOPS
    Frédéric Schmidt
    Université Paris-Saclay, CNRS, GEOPS, 91405, Orsay, France
    S3 seminar - Central SUPELEC 25th September 2020
    samedi 26 septembre 20

    View Slide

  2. Plan
    • Different use of linear unmixing for detection
    • blind linear unmixing (what and where?)
    • supervised linear unmixing (where?)
    F. Schmidt Plan
    samedi 26 septembre 20

    View Slide

  3. Plan
    • Different use of linear unmixing for detection
    •blind linear unmixing (what and where?)
    • supervised linear unmixing (where?)
    F. Schmidt Plan
    samedi 26 septembre 20

    View Slide

  4. • Unsupervised automatic detection
    • Summary (quick look) of the dataset
    • Based on advance signal treatment and
    machine learning
    • General tool but application to NOMAD/
    Exomars TGO
    • NO PRECISE QUANTIFICATION
    Goal
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  5. Theory
    • Estimation of endmember spectra (source) and
    abundances under constraints:
    observation
    mittance space:
    T⇤
    (

    ) =
    T
    (

    )
    C
    (

    )
    the spectra into absorbance:
    X
    (

    ) = 1
    T⇤
    (

    )
    is the linear mixture :
    X
    (

    )

    NS
    X
    i=1
    Si(

    )
    .Ai
    e source spectra and Ai
    the spectral abunda
    subject to
    positivity
    Number of
    endmember/source
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  6. Non Negative Matrix
    Factorization (NNMF)
    • fast
    • not robust
    Theory
    • Estimation of endmember spectra (source) and
    abundances under constraints:
    observation
    X
    (

    ) = 1
    T⇤
    (

    )
    p is the linear mixture :
    X
    (

    )

    NS
    X
    i=1
    Si(

    )
    .Ai
    he source spectra and Ai
    the spectral abund
    hysical meaning of Si(

    )
    and Ai
    is lost but th
    4
    out of the collection, with the highest abundance of a s
    In the following, we will use the continuum estimation
    c least square (Eilers and Boelens, 2005), with parameter
    = 1 10
    2, 10 number of iterations.
    2. Non negative matrix factorization
    For a collection of spectra, eq. 5 can be written in matri
    th i the source index (from 1 to NS
    ) , j the observation
    d k the wavenumber index (from 1 to N⌫
    ). Thus, one ha
    , by minimizing the objective function:
    F
    =
    kX S.Ak2
    with k.k, the Frobenius norm (usual L2
    norm).
    Several algorithms have been proposed to solve this
    MU Gillis, N. & Glineur, F. 2012, http://dx.doi.org/10.1162/neco_a_00256
    subject to
    positivity
    Method
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  7. Sparse Non Negative Matrix Factorization (sNNMF)
    • fast
    • sparsity
    • how to determine hyperparameter ?
    Theory
    • Estimation of endmember spectra (source) and
    abundances under constraints:
    observation
    X
    (

    ) = 1
    T⇤
    (

    )
    p is the linear mixture :
    X
    (

    )

    NS
    X
    i=1
    Si(

    )
    .Ai
    he source spectra and Ai
    the spectral abund
    hysical meaning of Si(

    )
    and Ai
    is lost but th
    4
    subject to
    positivity
    Method
    proach to reduce the computation time is to select only
    the dataset (Moussaoui et al., 2008), but then the stat
    chmidt et al., 2010). Thanks to the advances of comp
    opose to treat the full dataset. This kind of algorithm is
    e formulation is Bayesian, it converge toward an unique
    NMF. In order to regularize the problem of eq. 6, on
    nalization term to enforce sparsity on A (only few non
    im and Park, 2007) :
    F
    =
    kX S.Ak2
    +
    kAk
    1
    With k.k
    1
    , the L1
    norm. The first term is called data a
    ual squared difference). The second is called regularizat
    Li, Y., Ngom, A.: The non-negative matrix factorization toolbox for biological data mining.
    Source Code Biol. Med. 8(1), 10 (2013)
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  8. • each component as gamma
    • MCMC approach
    • no tuning parameter
    • slow computation
    Method
    Bayesian Prior Source Separation
    BPSS Moussaoui, S.; et al. 2008, http://dx.doi.org/10.1016/j.neucom.2007.07.034
    Schmidt, F.; et al., 2010,http://dx.doi.org/10.1109/TGRS.2010.2062190
    N OF NON-NEGATIVE MIXTURE OF NON-NEGATIVE SOURCES
    By assuming the mutual independence of the source s
    and the mixing coefficients, their associated prior densit
    then expressed by
    where the vectors and
    contain the parameters of the Gamma distributions.
    C. Posterior Density and Estimation Issues
    T⇤
    (

    ) =
    T
    (

    )
    C
    (

    )
    (3)
    ert the spectra into absorbance:
    X
    (

    ) = 1
    T⇤
    (

    )
    (4)
    tep is the linear mixture :
    X
    (

    )

    NS
    X
    i=1
    Si(

    )
    .Ai
    (5)
    the source spectra and Ai
    the spectral abundance. In this de-
    physical meaning of Si(

    )
    and Ai
    is lost but the apparent SNR
    4
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  9. Probabilistic Sparse Non-negative Matrix Factorization
    psNMF
    • each component as truncated normal/exponential
    • no explicit sparsity
    • variational Bayesian approach
    • no tuning parameter (soft prior with low )
    Method
    Hinrich, J. L. & Mørup, M., 2018, http://dx.doi.org/10.1007/978-3-319-93764-9_45
    490 J. L. Hinrich and M. Mørup
    λd
    ∼ Gamma (αλ, βλ
    ) , (1)
    wid
    | λd
    ∼ T N 0, λ−1
    d
    , 0, ∞ or wid
    |λd
    ∼ Exponential λ−1
    d
    , (2)
    hdj
    | λd
    ∼ T N 0, λ−1
    d
    , 0, ∞ or hdj
    |λd
    ∼ Exponential λ−1
    d
    , (3)
    τ ∼ Gamma (ατ , βτ
    ) , x
    j
    | W, h
    j, τ ∼ N x
    j
    | Wh
    j, τ−1II . (4)
    Where the reconstruction error x
    j
    − Wh
    j
    follows a normal distribution with
    Bayesian modeling ·
    1 Introduction
    Non-negative matrix fact
    ization [1] has become a
    pretable part based repr
    XI×J ≈ W I×DHD×J in w
    i.e. wid
    ≥ 0 ∀i, d and hdj
    in NMF have been propo
    this includes multiplicativ
    component-wise updating
    Unfortunately, NMF
    does not sufficiently span
    tions may equally well re
    λd
    ∼ Gamma (αλ, βλ
    ) , (1)
    wid
    | λd
    ∼ T N 0, λ−1
    d
    , 0, ∞ or wid
    |λd
    ∼ Exponential λ−1
    d
    , (2)
    hdj
    | λd
    ∼ T N 0, λ−1
    d
    , 0, ∞ or hdj
    |λd
    ∼ Exponential λ−1
    d
    , (3)
    τ ∼ Gamma (ατ , βτ
    ) , x
    j
    | W, h
    j, τ ∼ N x
    j
    | Wh
    j, τ−1II . (4)
    Where the reconstruction error x
    j
    − Wh
    j
    follows a normal distribution with
    se precision (inverse variance) τ. Each component of W and H (i.e. w
    d, h
    d
    )
    res a common prior λd
    used to infer the scale of a component, see also [14].
    te λd
    represents the precision or rate for the truncated normal and exponential
    tribution, respectively. This prior can be used for model order selection and
    commonly called an automatic relevance determination (ARD) prior. The
    terior distribution of θ = {W, H, τ, {λd
    }d=1,2,...,D
    } is inferred from the data
    while the shape α⋆
    and rate β⋆
    of the Gamma distributions are fixed. In
    ence of a priori information, weak priors are specified (c.f. α⋆
    = β⋆
    = 10−4),
    reby the distribution of τ and λd
    are determined primarily from the data.
    Probabilistic Sparse Non-negative Matrix Factorization
    obabilistic modeling can be used to automatically infer the sparsity pattern
    d level, as supported by the data. The Bayesian framework facilitates sparse
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  10. Results on Mars
    atmospheric data
    • Methane (life) on Mars ?
    • PH3 (life) on Venus ?
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  11. NOMAD/
    ExoMars TGO
    Fig. 8. Evolution of the slit angle during a typical occultation (inclination exag-
    gerated for clarity).
    on in the UV–visible (Left Panel) and the IR (Right Panel) recorded by the UVIS and SO
    1 ppb) for a typical clear Mars atmosphere, considering the Rayleigh scattering and a
    ond. For the IR spectra: The limits of the different diffraction orders covered are shown
    transmittances have been added in light green to clearly indicate where its absorption
    Fig. 6. Examples of simulated transmittances obtained during a typical solar occultation in the UV–visible (Left Panel) and the IR (Right Panel) recorded by the UVIS and SO
    channels respectively. These spectra contain the absorption of CO2
    , H2
    O, O3
    and CH4
    (1 ppb) for a typical clear Mars atmosphere, considering the Rayleigh scattering and a
    dust loading of tau¼0.2. The colour code indicates the altitudes to which they correspond. For the IR spectra: The limits of the different diffraction orders covered are shown
    at the top of the Figure and spectra have been artificially shifted by 0.40 for clarity. CH4
    transmittances have been added in light green to clearly indicate where its absorption
    lines are located.
    A.C. Vandaele et al. / Planetary and Space Science 119 (2015) 233–249
    242
    insulation (MLI), and instrument-to-spacecraft mounting
    hardware (3.7 kg).
    NOMAD has to survive in the environmental conditions
    imposed by the spacecraft. The most severe constraints are the
    Fig. 1. Different observation modes with NOMAD in orbit around
    Mars (1 ˆ nadir, 2 ˆ limb, 3 ˆ SO).
    Table 1. Coalignment Contributors
    Coalignment
    Contributors
    Accuracy
    Limit
    Knowledge
    Accuracy
    Solar LOS to NOMAD
    mechanical axis
    0.15 mrad
    Nadir LOS to NOMAD
    mechanical axis
    10.0 mrad
    NOMAD mechanical
    axis to spacecraft axis
    0.20 mrad 0.05 mrad
    Any NOMAD solar
    LOS to spacecraft axisa
    0.30 mrad
    aIt is assumed that the ACS instrument, which performs SO measurements at
    the same time as NOMAD, has a similar coalignment budget of 0.30 mrad,
    leading to a maximum misalignment between the NOMAD and ACS solar
    lines of sight of 0.60 mrad.
    Table 2. Pointing Budget
    Pointing Error Limit Knowledge
    Accuracy
    Absolute pointing error (APE)
    solar
    ≤1.23 mrad
    Relative pointing error (RPE)
    solar (short + long term)a
    ≤1.23 mrad
    Absolute pointing error (APE)
    nadir
    ≤3.50 mrad
    Relative pointing error (RPE)
    nadir (long term)a
    ≤3.00 mrad
    Relative pointing error (RPE)
    nadir (short term)a
    ≤0.54 mrad
    Overall pointing knowledge solar ≤0.45 mrad
    Overall pointing knowledge nadir ≤0.54 mrad
    aShort term is 1 s, long term is 60 s.
    8496 Vol. 54, No. 28 / October 1 2015 / Applied Optics Research Article
    Vandaele, A.; et al, Planetary and Space Science, Elsevier BV, 2015, 119, 233-249, http://dx.doi.org/10.1016/j.pss.2015.10.003
    Transmittance
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  12. NOMAD/ExoMars TGO
    • Methane (life) on Mars ?
    • detection limits down to 10 ppt = 10 part per
    trillion = 10-11 !
    • how to discover new unexpected species out
    of large dataset ?!
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  13. • Synthetic toy model
    • Simulations
    • Real data
    Results on Mars atmospheric
    data
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  14. •Synthetic toy model
    • Simulations
    • Real data
    Results on Mars atmospheric
    data
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  15. -6 -4 -2 0 2 4
    BD -3
    0
    50
    100
    150
    200
    frequency
    Histogram of BD CH4 @ 3067.2 cm-1
    3060 3065 3070 3075
    -1
    0.002
    0.004
    0.006
    0.008
    0.01
    0.012
    0.014
    0.016
    0.018
    0.02
    Absorbance
    Spectra max BD CH4 @ 3067.2 cm-1
    Toy model
    Figure 1: Synthetic dataset containing 104 spectra with various abundances of H2O and 100
    containing CH4 at 3- level of the noise. In blue the reference spectra SH2O of H2O (coming
    from actual data analysis). In red the reference spectra SCH4
    of CH4 (from theoretical data).
    The final synthetic dataset is represented in Fig. 1.
    221
    In order to check the quality of the estimation, we simply compute the
    222
    correlation coefficient between SCH4
    and the estimated NS
    sources ˙
    S, using:
    223
    Q
    =
    corr nSCH4
    , ˙
    S:i
    o (12)
    The ith source with the maximum correlation is identified to CH4
    contri-
    224
    bution. The value to the maximum correlation is used as metric to assess the
    225
    quality of the retrieval.
    226
    4.1.2. Results
    227
    By plotting the 10000 samples of the dataset, one is able to identify easily
    228
    the H2
    O bands. Nevertheless, we cannot observe the target CH4
    in the average
    229
    spectrum, even at 3- level, because it is lost in the baseline changes.
    230
    • linear mixture
    • 100 spectra of CH4 at 3-sigma
    10 000 spectra
    example, we simulate a linear mixture of NO = 10
    observations
    ⌫ = 320
    spectels (see fig. 1) similar to order 136 of NOMAD-
    rum is a mixture of a spectra of water vapor SH2
    O
    (coming
    source estimated from real data using psNMF) and theoretical
    rom Villanueva et al. (2018), with corresponding abundances
    X
    =
    SH2
    O
    .AH2
    O +
    SCH4
    .ACH4 +
    n (11)
    is assumed to be a Gaussian process with a standard deviation
    no bias: n
    =
    G
    (0
    ,
    )
    . All spectra contain pure water vapor with
    owing AH2
    O = 5
    /
    6
    .
    (1
    ,
    10) + 1
    /
    6
    .U
    (0
    ,
    1)
    , a mixture of beta ( )
    5/6 of the sample and an uniform (U) distribution for 1/6 of
    s process mimics well the water vapor band depth distribution
    on in section 3.3) of the real dataset (see Fig.2). As the baseline
    ero, we also mimic baseline correction errors. In addition 100
    0000 contain methane with ACH4 = 1
    , such that the band depth
    level. Please note that the model to generate the data is not
    m-to-one constraint, but fully fulfilling the positivity constraint.
    d noise and signal level, the noise RMSD is 0.16.
    7
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  16. Toy model
    • linear mixture
    • 100 spectra of CH4
    at 3-sigma
    10 000 spectra
    Figure 1: Synthetic dataset containing 104 spectra with various abundances of H2O and 100
    containing CH4 at 3- level of the noise. In blue the reference spectra SH2O of H2O (coming
    from actual data analysis). In red the reference spectra SCH4
    of CH4 (from theoretical data).
    The final synthetic dataset is represented in Fig. 1.
    221
    In order to check the quality of the estimation, we simply compute the
    222
    correlation coefficient between SCH4
    and the estimated NS
    sources ˙
    S, using:
    223
    Q
    =
    corr nSCH4
    , ˙
    S:i
    o (12)
    The ith source with the maximum correlation is identified to CH4
    contri-
    224
    bution. The value to the maximum correlation is used as metric to assess the
    225
    quality of the retrieval.
    226
    4.1.2. Results
    227
    By plotting the 10000 samples of the dataset, one is able to identify easily
    228
    the H2
    O bands. Nevertheless, we cannot observe the target CH4
    in the average
    229
    spectrum, even at 3- level, because it is lost in the baseline changes.
    230
    is able to detect the hidden compounds.
    example, we simulate a linear mixture of NO = 10
    4 observations
    ⌫ = 320
    spectels (see fig. 1) similar to order 136 of NOMAD-
    trum is a mixture of a spectra of water vapor SH2
    O
    (coming
    source estimated from real data using psNMF) and theoretical
    from Villanueva et al. (2018), with corresponding abundances
    X
    =
    SH2
    O
    .AH2
    O +
    SCH4
    .ACH4 +
    n (11)
    is assumed to be a Gaussian process with a standard deviation
    no bias: n
    =
    G
    (0
    ,
    )
    . All spectra contain pure water vapor with
    owing AH2
    O = 5
    /
    6
    .
    (1
    ,
    10) + 1
    /
    6
    .U
    (0
    ,
    1)
    , a mixture of beta ( )
    5/6 of the sample and an uniform (U) distribution for 1/6 of
    is process mimics well the water vapor band depth distribution
    ion in section 3.3) of the real dataset (see Fig.2). As the baseline
    zero, we also mimic baseline correction errors. In addition 100
    0000 contain methane with ACH4 = 1
    , such that the band depth
    level. Please note that the model to generate the data is not
    m-to-one constraint, but fully fulfilling the positivity constraint.
    ed noise and signal level, the noise RMSD is 0.16.
    -6 -4 -2 0 2 4
    BD 10-3
    0
    50
    100
    150
    200
    frequency
    Histogram of BD CH4 @ 3067.2 cm-1
    3060 3065 3070 3075
    wavenumber [cm-1]
    0.002
    0.004
    0.006
    0.008
    0.01
    0.012
    0.014
    0.016
    0.018
    0.02
    Absorbance
    Spectra max BD CH4 @ 3067.2 cm-1
    Figure 3: (left) Histogram of Band Depth at 3067.2 cm 1 from the dataset containing 100
    CH4 at 3- level out of 104 spectra. (right) 100 spectra with the maximum Band Depth at
    3067.2 cm 1 specific of CH4. Signal is dominated by water and by noise. No specific signature
    of CH4 is visible.
    3060 3065 3070 3075
    wavelength
    0.005
    0.01
    0.015
    0.02
    0.025
    0.03
    0.035
    0.04
    0.045
    Intensity
    psNMF : Endmember spectra
    No1 relevance =36.77%
    No2 relevance =32.6%
    No3 relevance =30.24%
    No4 relevance =0.39%
    Figure 3: (left) Histogram of Band Depth at 3067.2 cm 1 from the dataset containing 100
    CH4 at 3- level out of 104 spectra. (right) 100 spectra with the maximum Band Depth at
    3067.2 cm 1 specific of CH4. Signal is dominated by water and by noise. No specific signature
    of CH4 is visible.
    3060 3065 3070 3075
    wavelength
    0.005
    0.01
    0.015
    0.02
    0.025
    0.03
    0.035
    0.04
    0.045
    Intensity
    psNMF : Endmember spectra
    No1 relevance =36.77%
    No2 relevance =32.6%
    No3 relevance =30.24%
    No4 relevance =0.39%
    omial, that is filtered out.
    this continuum removal rationale is that when the optical
    R is decreased and the noise effect on continuum removal
    at.).
    is rationale, we propose to first correct for the continuum
    ance space:
    T⇤
    (

    ) =
    T
    (

    )
    C
    (

    )
    (3)
    pectra into absorbance:
    X
    (

    ) = 1
    T⇤
    (

    )
    (4)
    e linear mixture :
    X
    (

    )

    NS
    X
    i=1
    Si(

    )
    .Ai
    (5)
    • Results : correlation 0.98
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  17. Toy model
    • linear mixture
    10 000 spectra
    spanned over N⌫ = 320
    spectels (see fig. 1) similar to order 136 of N
    205
    SO. Each spectrum is a mixture of a spectra of water vapor SH2
    O
    206
    from one actual source estimated from real data using psNMF) and th
    207
    methane SCH4
    from Villanueva et al. (2018), with corresponding ab
    208
    AH2
    O
    , ACH4
    :
    209
    X
    =
    SH2
    O
    .AH2
    O +
    SCH4
    .ACH4 +
    n
    The noise n is assumed to be a Gaussian process with a standard d
    210
    of =0.001 and no bias: n
    =
    G
    (0
    ,
    )
    . All spectra contain pure water va
    211
    a coefficient following AH2
    O = 5
    /
    6
    .
    (1
    ,
    10) + 1
    /
    6
    .U
    (0
    ,
    1)
    , a mixture of
    212
    distribution for 5/6 of the sample and an uniform (U) distribution f
    213
    the sample. This process mimics well the water vapor band depth dis
    214
    (BD, see definition in section 3.3) of the real dataset (see Fig.2). As the
    215
    of SH2
    O
    is not zero, we also mimic baseline correction errors. In add
    216
    spectra out of 10000 contain methane with ACH4 = 1
    , such that the ba
    217
    of SCH4
    is at 3- level. Please note that the model to generate the da
    218
    fulfilling the sum-to-one constraint, but fully fulfilling the positivity co
    219
    Given the defined noise and signal level, the noise RMSD is 0.16.
    220
    7
    MU psNMF BPSS2
    Quality Q 0.35±0.12 0.822 ±0.005 0.41±0.06
    RMSD relative error 0.1455±
    2
    .
    10
    6 0.1461±
    5
    .
    10
    6 0.1468±
    3
    .
    10
    4
    Computation time (s) 13±8 46±9 413±21
    Table 1: Results (mean and standard deviation) from 10 realizations of a toy synthetic example
    with NS = 5 (in agreement with next section on synthetic tests), NO = 10000, N⌫ = 320 and
    300 CH4 spectra hidden at a level of 1 std of the noise. Quality is computed as a correlation
    coefficient (see Eq. 12). RMSD is computed from Eq. 8. Computation time is expressed in
    second.
    usual squared difference). The second is called regularization term.
    161
    lem with this approach, is that hyperparameter is not known an
    162
    tuned manually. A recent approach has been proposed to solve this
    163
    the Bayesian framework (Hinrich and Mørup, 2018). The main idea i
    164
    pass all variables and hyperparameters in a unique problem that i
    165
    with variational update principle. We will refer this algorithm to
    166
    sparse NMF (psNMF). This algorithm has the advantage to have
    167
    computation time and no hyperparameter tuning. It also has a reg
    168
    term to avoid strong dependence of the initialization on the final so
    169
    In order to estimate the precision of the reconstruction, we use
    170
    Mean Square Difference RMSD:
    171
    RMSD
    =
    s⌧⇣X ˙
    S. ˙
    A⌘2
    hXi
    With h.i, the mean.
    172
    Once the sources are estimated, we quantify their relevance for
    173
    dataset. From the total reconstruction ˙
    Xkj=
    ˙
    Ski
    . ˙
    Aij
    , for all i, we c
    174
    the contribution of source i0, that is to say: ˙
    Xi
    kj=
    ˙
    Ski0
    . ˙
    Ai0j
    . Thus, th
    175
    of source i is defined as:
    176
    Ri
    =
    D ˙
    Xi ˙
    X E
    D ˙
    XE
    This definition is convenient since the sum of all Ri is one (th
    177
    is only present when sources and abundances are positive) and we
    178
    estimate the % contribution of each source in the final reconstructio
    179
    dataset containing 104 spectra with various abundances of H2O and 100
    level of the noise. In blue the reference spectra SH2O of H2O (coming
    ysis). In red the reference spectra SCH4
    of CH4 (from theoretical data).
    hetic dataset is represented in Fig. 1.
    heck the quality of the estimation, we simply compute the
    ient between SCH4
    and the estimated NS
    sources ˙
    S, using:
    Q
    =
    corr nSCH4
    , ˙
    S:i
    o (12)
    e with the maximum correlation is identified to CH4
    contri-
    e to the maximum correlation is used as metric to assess the
    ieval.
    e 10000 samples of the dataset, one is able to identify easily
    Nevertheless, we cannot observe the target CH4
    in the average
    3- level, because it is lost in the baseline changes.
    • 300 spectra of CH4 at 1-sigma
    • average of 10 realizations
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  18. Toy model
    • linear mixture
    10 000 spectra
    SO. Each spectrum is a mixture of a spectra of water vapor SH2
    O
    206
    from one actual source estimated from real data using psNMF) and th
    207
    methane SCH4
    from Villanueva et al. (2018), with corresponding ab
    208
    AH2
    O
    , ACH4
    :
    209
    X
    =
    SH2
    O
    .AH2
    O +
    SCH4
    .ACH4 +
    n
    The noise n is assumed to be a Gaussian process with a standard d
    210
    of =0.001 and no bias: n
    =
    G
    (0
    ,
    )
    . All spectra contain pure water va
    211
    a coefficient following AH2
    O = 5
    /
    6
    .
    (1
    ,
    10) + 1
    /
    6
    .U
    (0
    ,
    1)
    , a mixture of
    212
    distribution for 5/6 of the sample and an uniform (U) distribution f
    213
    the sample. This process mimics well the water vapor band depth dis
    214
    (BD, see definition in section 3.3) of the real dataset (see Fig.2). As the
    215
    of SH2
    O
    is not zero, we also mimic baseline correction errors. In add
    216
    spectra out of 10000 contain methane with ACH4 = 1
    , such that the ba
    217
    of SCH4
    is at 3- level. Please note that the model to generate the da
    218
    fulfilling the sum-to-one constraint, but fully fulfilling the positivity co
    219
    Given the defined noise and signal level, the noise RMSD is 0.16.
    220
    7
    • plateau at 5-10 endmembers
    figure 5
    3 4 5 6 7 8 9
    Nb endmember
    0
    0.1
    0.2
    0.3
    0.4
    0.5
    0.6
    0.7
    0.8
    0.9
    1
    Value
    Effect of Ns for noise factor 2 with 100 CH
    4
    spectra
    MU mean corr.
    psNMF mean corr.
    MU frac. corr. > 0.5
    psNMF frac. corr. > 0.5
    3 4 5 6 7 8 9
    Nb endmember
    0
    0.1
    0.2
    0.3
    0.4
    0.5
    0.6
    0.7
    0.8
    0.9
    1
    Value
    Effect of Ns for noise factor 3 with 100 CH
    4
    spectra
    MU mean corr.
    psNMF mean corr.
    MU frac. corr. > 0.5
    psNMF frac. corr. > 0.5
    Factor above noise level of 2
    Nb endmember
    Value
    100 hidden CH4 spectra
    Factor above noise level of 3
    Nb endmember
    Value
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  19. Toy model
    • linear mixture
    10 000 spectra
    SO. Each spectrum is a mixture of a spectra of water vapor SH2
    O
    206
    from one actual source estimated from real data using psNMF) and th
    207
    methane SCH4
    from Villanueva et al. (2018), with corresponding ab
    208
    AH2
    O
    , ACH4
    :
    209
    X
    =
    SH2
    O
    .AH2
    O +
    SCH4
    .ACH4 +
    n
    The noise n is assumed to be a Gaussian process with a standard d
    210
    of =0.001 and no bias: n
    =
    G
    (0
    ,
    )
    . All spectra contain pure water va
    211
    a coefficient following AH2
    O = 5
    /
    6
    .
    (1
    ,
    10) + 1
    /
    6
    .U
    (0
    ,
    1)
    , a mixture of
    212
    distribution for 5/6 of the sample and an uniform (U) distribution f
    213
    the sample. This process mimics well the water vapor band depth dis
    214
    (BD, see definition in section 3.3) of the real dataset (see Fig.2). As the
    215
    of SH2
    O
    is not zero, we also mimic baseline correction errors. In add
    216
    spectra out of 10000 contain methane with ACH4 = 1
    , such that the ba
    217
    of SCH4
    is at 3- level. Please note that the model to generate the da
    218
    fulfilling the sum-to-one constraint, but fully fulfilling the positivity co
    219
    Given the defined noise and signal level, the noise RMSD is 0.16.
    220
    7
    1 1.5 2 2.5 3
    Noise factor
    0
    0.1
    0.2
    0.3
    0.4
    0.5
    0.6
    0.7
    0.8
    0.9
    1
    Value
    Effect of noise factor with 100 CH
    4
    spectra
    MU mean corr.
    psNMF mean corr.
    MU frac. corr. > 0.5
    psNMF frac. corr. > 0.5
    1 1.5 2 2.5 3
    Noise factor
    0
    0.1
    0.2
    0.3
    0.4
    0.5
    0.6
    0.7
    0.8
    0.9
    1
    Value
    Effect of noise factor with 50 CH
    4
    spectra
    MU mean corr.
    psNMF mean corr.
    MU frac. corr. > 0.5
    psNMF frac. corr. > 0.5
    figure 6
    100 hidden CH4 spectra
    Factor above noise level
    Value
    Factor above noise level
    Value
    50 hidden CH4 spectra
    • high correlation up to 2-2.5 noise level
    1.5 2 2.5 3
    Noise factor
    of noise factor with 100 CH
    4
    spectra
    5
    1 1.5 2 2.5 3
    Noise factor
    0
    0.1
    0.2
    0.3
    0.4
    0.5
    0.6
    0.7
    0.8
    0.9
    1
    Value
    Effect of noise factor with 50 CH
    4
    spectra
    MU mean corr.
    psNMF mean corr.
    MU frac. corr. > 0.5
    psNMF frac. corr. > 0.5
    of the MU and psNMF algorithm for NS = 5, NO = 10000, N⌫ = 320, as a
    ctor above noise level. The average Q of 10 realizations of the best estimated
    6
    hidden CH4 spectra
    actor above noise level Factor above noise level
    Value
    50 hidden CH4 spectra
    s of the MU and psNMF algorithm for NS = 5, NO = 10000, N⌫ = 320, as a
    ctor above noise level. The average Q of 10 realizations of the best estimated
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  20. • Synthetic toy model
    •Simulations
    • Real data
    Results on Mars atmospheric
    data
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  21. Simulation
    • non-linear : GCM + radiative transfer + instrument
    simulation : 12 486 spectra, 106 occultation, order
    136, May-December 2018
    • parameters : H2O density + CH4 density + random
    addition of CH4 + noise level
    • average of 10 noise realizations
    Daerden, F.,et al., 2019.Icarus 326, 197–224.
    3060 3065 3070 3075 3080
    Wavenumber (cm-1)
    0
    0.01
    0.02
    0.03
    0.04
    0.05
    0.06
    Transmittance
    Simulation spectra
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  22. 10 2 10 3
    0
    0.2
    0.4
    0.6
    0.8
    1
    value
    Fraction of the 4 main peaks detected
    fraction CH
    4
    = 100%
    fraction CH
    4
    = 50%
    fraction CH
    4
    = 10%
    fraction CH
    4
    = 5%
    fraction CH
    4
    = 1%
    10 2 10 3
    0
    0.5
    1
    1.5
    2
    value
    Mean distance to the expected center [in pixel]
    10 2 10 3
    CH
    4
    density [in ppt]
    0
    0.05
    0.1
    0.15
    0.2
    0.25
    value
    Abundance of CH
    4
    in the source
    noise level 0.001
    noise level 0.0001
    10 2 10 3
    0
    0.2
    0.4
    0.6
    0.8
    1
    value
    Fraction of the 4 main peaks detected
    fraction CH
    4
    = 100%
    fraction CH
    4
    = 50%
    fraction CH
    4
    = 10%
    fraction CH
    4
    = 5%
    fraction CH
    4
    = 1%
    10 2 10 3
    0
    0.5
    1
    1.5
    2
    value
    Mean distance to the expected center [in pixel]
    10 2 10 3
    CH
    4
    density [in ppt]
    0
    0.05
    0.1
    0.15
    value
    Abundance of CH
    4
    in the source
    noise level 0.001
    noise level 0.0001
    Figure 7: Results of the psNMF algorithm for NS = 5 on simulation dataset, averaged over 10
    figure 7
    Fraction of the 4 main peaks detected
    Abundance of CH4 in the source
    Mean distance to the expected center [in pixel]
    CH4 density [in ppt] CH4 density [in ppt]
    Value Value Value
    Value
    Value
    Value
    Figure 7: Results of the psNMF algorithm for NS = 5 on simulation dataset, averaged over 10
    10 ppm H2O 100 ppm H2O
    • detection limit of CH4 100-500 ppt
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  23. • Synthetic toy model
    • Simulations
    •Real data
    Results on Mars atmospheric
    data
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  24. Real data pre-treatment
    all SO from order 134
    Pipeline
    • from version 1p0a (up to
    15 January 2020)
    • Resampling to correct for
    detector temperature
    • 0 < Altitude < 50 km
    • SNR > 100
    365 985 spectra
    Neefs, E.; et al, 2015, http://dx.doi.org/10.1364/ao.54.008494
    Vandaele, A. C. et al., 2018 http://dx.doi.org/10.1007/s11214-018-0517-2
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  25. all SO from order 134
    Pipeline
    • Continuum removal
    (ALS)
    • Convert to absorbance
    365 985 spectra
    Real data pre-treatment
    Eilers, P. H. & Boelens, H. F., Leiden University Medical Centre Report, 2005
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  26. all SO from order 134
    Pipeline
    • Shift correction
    (residual of detector
    temperature)
    365 985 spectra
    F. Schmidt Automatic Detection
    Real data pre-treatment
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  27. 3015 3020 3025 3030 3035
    0
    0.01
    0.02
    0.03
    0.04
    0.05
    0.06
    sources (total rms = 0.0012748)
    mean = 16%
    mean = 83%
    • 2 sources
    Results order 134
    ).
    ationale, we propose to first correct for the continuum
    e space:
    T⇤
    (

    ) =
    T
    (

    )
    C
    (

    )
    (3)
    ctra into absorbance:
    X
    (

    ) = 1
    T⇤
    (

    )
    (4)
    inear mixture :
    X
    (

    )

    NS
    X
    i=1
    Si(

    )
    .Ai
    (5)
    spectra and Ai
    the spectral abundance. In this de-
    eaning of Si(

    )
    and Ai
    is lost but the apparent SNR
    sources
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  28. 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
    Abundance
    5
    10
    15
    20
    25
    30
    35
    40
    45
    Altitude
    Abundances (total RMS = 0.0012748)
    mean = 16%
    mean = 83%
    10 000 spectra
    3015 3020 3025 3030 3035
    0
    0.01
    0.02
    0.03
    0.04
    0.05
    0.06
    sources (total rms = 0.0012748)
    mean = 16%
    mean = 83%
    abundances
    cles, scattering by molecules and particles, and also reflection at
    broadband features. Such large features are modeled by a conti
    often taken as a polynomial, that is filtered out.
    The problem with this continuum removal rationale is that w
    depth is large, the SNR is decreased and the noise effect on cont
    amplified (see Sup. Mat.).
    Instead of using this rationale, we propose to first correct for
    C
    (

    )
    in the transmittance space:
    T⇤
    (

    ) =
    T
    (

    )
    C
    (

    )
    Then convert the spectra into absorbance:
    X
    (

    ) = 1
    T⇤
    (

    )
    The final step is the linear mixture :
    X
    (

    )

    NS
    X
    i=1
    Si(

    )
    .Ai
    samedi 26 septembre 20

    View Slide

  29. all SO from order 134
    • 2358 orbits
    • 365 985 spectra
    3015 3020 3025 3030 3035
    Wavenumber (cm-1)
    0.005
    0.01
    0.015
    0.02
    0.025
    Intensity
    psNMF : Endmember spectra
    No1 relevance =31.07%
    No2 relevance =19.58%
    No3 relevance =19.44%
    No4 relevance =15.87%
    No5 relevance =14.04%
    3015 3020 3025 3030 3035
    Wavenumber (cm-1)
    0
    0.02
    0.04
    0.06
    0.08
    0.1
    Absorbance
    Theoretical absorbance trough NOMAD
    CO
    2
    N
    2
    O
    2
    CO
    H
    2
    O
    O
    3
    CH
    4
    Solar
    19
    Results
    Schmidt, F. et al., 2020, JQRST, under review
    ).
    ationale, we propose to first correct for the continuum
    e space:
    T⇤
    (

    ) =
    T
    (

    )
    C
    (

    )
    (3)
    ctra into absorbance:
    X
    (

    ) = 1
    T⇤
    (

    )
    (4)
    inear mixture :
    X
    (

    )

    NS
    X
    i=1
    Si(

    )
    .Ai
    (5)
    spectra and Ai
    the spectral abundance. In this de-
    eaning of Si(

    )
    and Ai
    is lost but the apparent SNR
    endmembers (sources)
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  30. all SO from order 134
    • 2358 orbits
    • 365 985 spectra
    3015 3020 3025 3030 3035
    Wavenumber (cm-1)
    0.005
    0.01
    0.015
    0.02
    0.025
    Intensity
    psNMF : Endmember spectra
    No1 relevance =31.07%
    No2 relevance =19.58%
    No3 relevance =19.44%
    No4 relevance =15.87%
    No5 relevance =14.04%
    3015 3020 3025 3030 3035
    Wavenumber (cm-1)
    0
    0.02
    0.04
    0.06
    0.08
    0.1
    Absorbance
    Theoretical absorbance trough NOMAD
    CO
    2
    N
    2
    O
    2
    CO
    H
    2
    O
    O
    3
    CH
    4
    Solar
    19
    H2O
    New lines CO2
    Simulation
    • Detection of H2O
    • New lines CO2 magnetic
    dipole
    • No detection of CH4
    Results
    Schmidt, F. et al., 2020, JQRST, under review
    Trokhimovskiy, A.; et al., 2020, http://
    dx.doi.org/10.1051/0004-6361/202038134
    Villanueva, G, 2018,http://dx.doi.org/10.1016/j.jqsrt.2018.05.023
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  31. 2675 2680 2685 2690 2695
    Wavenumber (cm-1)
    0.002
    0.004
    0.006
    0.008
    0.01
    0.012
    0.014
    0.016
    0.018
    0.02
    0.022
    Intensity
    psNMF : Endmember spectra
    No1 relevance =29.87%
    No2 relevance =24.52%
    No3 relevance =23.23%
    No4 relevance =12.67%
    No5 relevance =9.72%
    2675 2680 2685 2690 2695
    Wavenumber (cm-1)
    0
    0.005
    0.01
    0.015
    0.02
    0.025
    0.03
    0.035
    Absorbance
    Theoretical absorbance trough NOMAD
    CO
    2
    N
    2
    O
    2
    CO
    H
    2
    O
    O
    3
    CH
    4
    Solar
    Figure 8: Results of the psNMF algorithm for the diffraction order 119 for NS = 5. The
    18
    Figure 8: Results of the psNMF algorithm for the diffraction order 119 for NS = 5. The
    18
    all SO from order 119
    • 821 orbits
    • 134 045 spectra
    Results
    Schmidt, F. et al., 2020, JQRST, under review
    • Detection of H2O and
    CO2
    • No detection of CH4
    H2O
    CO2
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  32. Conclusion
    • New fast unsupervised data
    quicklook
    • Validated on simulation datasets
    • Results on real NOMAD data
    SO from order 119, 134, 136
    • Detection of H2O, CO2
    • New lines CO2 magnetic
    dipole
    • No detection of CH4
    Trokhimovskiy, A.; et al., 2020, http://
    dx.doi.org/10.1051/0004-6361/202038134
    Schmidt, F. et al., 2020, JQSRT, under review
    [email protected]@universite-paris-saclay.fr
    F. Schmidt Blind Linear Unmixing
    samedi 26 septembre 20

    View Slide

  33. Plan
    • Different use of linear unmixing for detection
    • blind linear unmixing (what and where?)
    •supervised linear unmixing (where?
    How much?)
    F. Schmidt Plan
    samedi 26 septembre 20

    View Slide

  34. Instrument
    Surface
    …!
    Surface!
    Atmosphere!
    target!
    11th!
    6th!
    1st!
    …!
    "≈70º!
    "≈70º!
    "≈0º!
    "≈30º!
    "≈30º!
    $
    in!
    $
    out!
    Atmosphere
    nadir
    image
    incoming
    EPF
    outgoing
    EPF
    Targeted observations
    11 multi-angle images
    Murchie et al., JGR, 2007
    10 off-nadir images (180 m/pxl)
    eme±70°, constant inc
    1 nadir image (20 m/pxl) credit: http://crism.jhuapl.edu
    Hyperspectral image
    438+107 bands
    (0.36 to 3.92 μm)
    CRISM (Compact Reconnaissance Imaging Spectrometer for
    Mars in Mars Reconnaissance Orbiter spacecraft)
    samedi 26 septembre 20

    View Slide

  35. Spectral analysis
    1. Classification
    • Supervised: knowing laboratory spectra
    • GOAL: Where are the reference spectra ?
    2. Radiative transfer inversion
    • Quantitative estimation of surface properties
    samedi 26 septembre 20

    View Slide

  36. Mathematical problem
    • Estimation of abundances, under constraints
    d length, (iv) the emergence direction is always t
    Thus, based on this model and using the geogra
    ture assumption, the radiance factor at location (
    at wavelenght λ satisfies the following observatio
    L(x, y, λ) =
    ρa
    (λ) + Φ(λ)
    P
    p=1
    αp
    (x, y) ρp
    (λ) cos [θ(x
    where Φ(λ) is the spectral atmospheric tran
    θ(x, y) the angle between the solar direction and
    face normal (solar incidence angle), P the n
    endmembers in the region of coordinates (x, y),
    us, based on this model and using the geographic mix-
    re assumption, the radiance factor at location (x, y) and
    wavelenght λ satisfies the following observation model:
    L(x, y, λ) =
    ρa
    (λ) + Φ(λ)
    P
    p=1
    αp
    (x, y) ρp
    (λ) cos [θ(x, y)] (1)
    ere Φ(λ) is the spectral atmospheric transmission,
    x, y) the angle between the solar direction and the sur-
    e normal (solar incidence angle), P the number of
    dmembers in the region of coordinates (x, y), ρp
    (λ) the
    ectrum of the p-th endmember, αp
    (x, y) its weight in
    e mixture and ρa
    (λ) the radiation that did not arrive
    ectly from the area under view. This mixture model can
    o be written as:
    min ||↵p.⇢p L||
    min ||↵p.⇢p L||, ↵p > 0,
    P
    ↵p
    = 1
    •Property : linearly
    dependent
    spectra in the
    database give one
    single solution !
    1 1.5 2 2.5
    0.1
    0.2
    0.3
    0.4
    0.5
    0.6
    0.7
    0.8
    0.9
    Wavelength (µm)
    Reflectance (%)
    Endmembers
    Clay Montmorillonite Bentonite
    Flat 1
    Flat 0.0001
    Positivity,
    sum-to-one
    samedi 26 septembre 20

    View Slide

  37. 1 1.5 2 2.5
    0.4
    0.45
    0.5
    0.55
    0.6
    Wavelength (µm)
    Reflectance (%)
    Absorption band thickness
    Clay Montmorillonite Bentonite (0%)
    Clay Montmorillonite Bentonite (10%)
    Clay Montmorillonite Bentonite (20%)
    Clay Montmorillonite Bentonite (30%)
    Clay Montmorillonite Bentonite (40%)
    Clay Montmorillonite Bentonite (50%)
    Clay Montmorillonite Bentonite (60%)
    Clay Montmorillonite Bentonite (70%)
    Clay Montmorillonite Bentonite (80%)
    Clay Montmorillonite Bentonite (90%)
    Clay Montmorillonite Bentonite (100%)
    Absorption band depth
    samedi 26 septembre 20

    View Slide

  38. 1 1.5 2 2.5
    0.35
    0.4
    0.45
    0.5
    0.55
    0.6
    0.65
    Wavelength (µm)
    Reflectance (%)
    Continuum shift
    Flat 1 (0%)
    Flat 1 (1%)
    Flat 1 (2%)
    Flat 1 (3%)
    Flat 1 (4%)
    Flat 1 (5%)
    Flat 1 (6%)
    Flat 1 (7%)
    Flat 1 (8%)
    Flat 1 (9%)
    Flat 1 (10%)
    Continuum shift
    samedi 26 septembre 20

    View Slide

  39. To model the effect of :
    • aerosols
    • continuum
    • surface roughness of the regolith
    • grain size, shape roughness
    • mixture (opaque, feldspar, ...)
    Adding
    other
    spectra ?
    Remove discrepancies between
    observed spectra and database
    samedi 26 septembre 20

    View Slide

  40. 1 1.5 2 2.5
    0.35
    0.4
    0.45
    0.5
    0.55
    0.6
    Wavelength (µm)
    Reflectance (%)
    Maximum shift
    sin 1/2 (0%)
    sin 1/2 (1%)
    sin 1/2 (2%)
    sin 1/2 (3%)
    sin 1/2 (4%)
    sin 1/2 (5%)
    sin 1/2 (6%)
    sin 1/2 (7%)
    sin 1/2 (8%)
    sin 1/2 (9%)
    sin 1/2 (10%)
    Maximum shift ?
    samedi 26 septembre 20

    View Slide

  41. Spectral database
    • Selection of spectra for Mars
    Schmidt, F.; Legendre, M. & Le Mouëlic, S. Minerals detection for hyperspectral images using adapted
    linear unmixing: LinMin Icarus, 2014, 237, 61-74, http://dx.doi.org/10.1016/j.icarus.2014.03.044
    fitting. We also showed that the knowledge of the noise covariance matrix, that
    can be estimated from dark current or using other techniques is important to
    asses the detection limits.
    IPLS is shown to be the best numerical algorithm to solve the MINESINE
    problem. Its fast GPU implementation is particularly relevant for the treatment
    of hyperspectral images. In the future, this methodology should be applied in
    various planetary cases in order to study the surface geology. Also a signifi-
    cant improvement of the mineral detection may be addressed by using spectral
    database adapted to the context.
    Acknowledgement
    We acknowledge support from the “Institut National des Sciences de l’Univers”
    (INSU), the "Centre National de la Recherche Scientifique" (CNRS) and "Cen-
    tre National d’Etude Spatiale" (CNES) and through the "Programme National
    de Planétologie".
    Appendix
    Name of the 32 spectra:
    1 Inosilicate (Hypersthene OPX PYX02.h >250u) 12 Sulfate; Gypsum 23 Carbonate; Siderite
    2 Inosilicate (Diopside CPX CRISM) 13 Sulfate; Jarosite 24 Phyllosilicate (Chlorite)
    3 Olivine Fayalite CRISM 14 Sulfate; Kieserite 25 Muscovite GDS116 Tanzania
    4 Olivine Forsterite CRISM 15 Epsomite USGS GDS149 26 Alunite GDS83 Na63
    5 Phyllosilicate (Clay Montmorillonite Bentonite) 16 Oxide; Goethite 27 Atmospheric Transmission
    6 Phyllosilicate (Clay Illite Smectite) 17 Oxide; Hematite 28 H2O grain 1
    7 Phyllosilicate (Serpentine Chrysotile Clinochry.) 18 Oxide; Magnetite 29 H2O grain 100
    8 Phyllosilicate (Serpentine Lizardite) 19 Ferrihydrite USGS GDS75 Sy F6 30 H2O grain 1000
    9 Phyllosilicate (Clay Illite) 20 Maghemite USGS GDS81 Sy (M-3) 31 CO2 grain 100
    10 Phyllosilicate (Clay Kaolinite) 21 Carbonate; Calcite 32 CO2 grain 10 000
    11 Phyllosilicate (Nontronite) 22 Carbonate; Dolomite
    Name of the 12 additional spectra:
    33 Flat 1 37 cos 1/4 41 cos 1/2
    34 Flat 0.0001 38 sin 1/4 42 sin 1/2
    35 Slope Increasing 39 -cos 1/4 43 -cos 1/2
    36 Slope Decreasing 40 -sin 1/4 44 -sin 1/2
    References
    1 1.5 2 2.5
    0
    0.1
    0.2
    0.3
    0.4
    0.5
    0.6
    0.7
    0.8
    0.9
    1
    Endmember
    Wavelength (microns)
    Reflectance
    Figure 2: 32 Reference spectra of minerals, ice and atmospheric gas representing major cla
    of contributions of surface spectra. See Appendix for the names.
    • 90% is a flat component at reflectance 0.35 in agreement with OMEG
    studies Vincendon in order to reproduce the low level and flatness of act
    Martian spectra. The flatness is also representative of the Moon or ot
    planetary surface.
    • 10% of a random mixture of two over 32 reference spectra with rand
    uniform mixing coe cients, noted as A0
    . For each endmember i, th
    is ~30 spectral mixture with non-null mixing coe cient (noted Apositi
    and ~970 with null mixing coe cient (Afalse).
    We add synthetic noise, simulating the noise level of a typical OMEG
    observation after gas correction (covariance matrix from dark current noise
    ORB41_1, which has an wavelength-average standard deviation of 1.3 10
    Other noise statistics can be used, such MNF shift di⌧erence from ENVI so
    ware, but we point the fact that OMEGA dark current is archived, estimat
    the minimum noise statistics (excluding spike and other non-linear e⌧ects).
    On a desktop with Dual Core at 2.53 Ghz with 4Go RAM memory,
    typical computation time to solve the non-normalized problem of eq. 2
    M = 104000 spectra with (1000 mixture of spectra and 104 values of AOT,
    next section), N = 110, N = 44 are : 3.3 h for BI-ICE, 4.5 min for IPLS,
    h for Simplex Projection. The estimated mixing coe cients are noted AIP
    5
    samedi 26 septembre 20

    View Slide

  42. Algorithm
    • FCLS is not doing the job
    • Primal-dual interior-point
    • GPU implementation
    Chouzenoux, E.; Legendre, M.; Moussaoui, S. & Idier, J. Fast Constrained Least Squares
    Spectral Unmixing Using Primal-Dual Interior-Point Optimization Selected Topics in
    Applied Earth Observations and Remote Sensing, IEEE Journal of, 2014, 7, 59-69, http://
    dx.doi.org/10.1109/JSTARS.2013.2266732
    ted
    olution
    e chan-
    els. We
    hannels
    be dis-
    ocused
    ta cube
    Mars in
    st were
    p with
    73 µm
    es from
    ter cal-
    express
    Thus, based on this model and using t
    ture assumption, the radiance factor a
    at wavelenght λ satisfies the following
    L(x, y, λ) =
    ρa
    (λ) + Φ(λ)
    P
    p=1
    αp
    (x, y) ρp

    where Φ(λ) is the spectral atmosp
    θ(x, y) the angle between the solar dir
    face normal (solar incidence angle),
    endmembers in the region of coordina
    spectrum of the p-th endmember, αp
    ted
    chan-
    ls. We
    annels
    be dis-
    ocused
    a cube
    ars in
    t were
    p with
    73 µm
    s from
    er cal-
    xpress
    etween
    nd the
    n pho-
    ugh its
    at wavelenght λ satisfies the following observation m
    L(x, y, λ) =
    ρa
    (λ) + Φ(λ)
    P
    p=1
    αp
    (x, y) ρp
    (λ) cos [θ(x, y)
    where Φ(λ) is the spectral atmospheric transm
    θ(x, y) the angle between the solar direction and th
    face normal (solar incidence angle), P the num
    endmembers in the region of coordinates (x, y), ρp
    (
    spectrum of the p-th endmember, αp
    (x, y) its wei
    the mixture and ρa
    (λ) the radiation that did not
    directly from the area under view. This mixture mod
    also be written as:
    L(x, y, λ) =
    P
    α′
    p
    (x, y) · ρ′
    p
    (λ) + E(x, y, λ)
    min ||↵p.⇢p L||
    min ||↵p.⇢p L||, ↵p > 0,
    P
    ↵p
    = 1
    Legendre, M.; Capriotti, L.; Schmidt, F.; Moussaoui, S. & Schmidt, A.
    GPU implementation issues for fast unmixing of hyperspectral images
    EGU General Assembly Conference Abstracts, 2013, 15, 11686,
    Schmidt, F.; Legendre, M. & Le Mouëlic, S. Minerals detection for hyperspectral images using adapted
    linear unmixing: LinMin Icarus, 2014, 237, 61-74, http://dx.doi.org/10.1016/j.icarus.2014.03.044
    Heinz, D. & I-Chang, C., TGRS, 2001, 39, 529-545,
    samedi 26 septembre 20

    View Slide

  43. Synthetic test 1
    • Synthetic spectra of:
    • 90% Flat at 0.35 (average Mars)
    • 10% random mixture of one/two components
    • Radiative transfer :
    • DISORT : non-linear
    • Martian aerosols from AOT=0 to AOT=100
    • Adding instrumental noise from dark current
    Lin, Z.; et al.,Improved discrete ordinate solutions in the presence of an anisotropically reflecting
    lower boundary: Upgrades of the DISORT computational tool Journal of Quantitative Spectroscopy
    and Radiative Transfer, 2015, 157, 119 - 134, http://dx.doi.org/http://dx.doi.org/10.1016/j.jqsrt.
    2015.02.014
    Wolff, M. et al., Wavelength dependence of dust aerosol single scattering albedo as observed by the Compact Reconnaissance Imaging Spectrometer J. Geophys.
    Res., 2009, 114, E00D04-, http://dx.doi.org/10.1029/2009JE003350
    Schmidt, F.; Legendre, M. & Le Mouëlic, S. Minerals detection for hyperspectral images using adapted
    linear unmixing: LinMin Icarus, 2014, 237, 61-74, http://dx.doi.org/10.1016/j.icarus.2014.03.044
    samedi 26 septembre 20

    View Slide

  44. 1 1.5 2 2.5
    0.34
    0.36
    0.38
    0.4
    0.42
    0.44
    0.46
    0.48
    Detection of Gypsum (10 %)
    Wavelength
    Reflectance (offset for clarity)
    aot=0 (ab = 9.6 %)
    aot=1 (ab = 5.1 %)
    aot=5 (ab = 1.1 %)
    aot=100 (ab = 0 %)
    1
    0.32
    0.34
    0.36
    0.38
    0.4
    0.42
    0.44
    0.46
    Reflectance (offset for clarity)
    fa
    fa
    fa
    fa
    fa
    Figure 4: Examples of results using the IPLS algorithm, with renormalizat
    More
    aerosols
    0.7%
    37%
    100%
    0.0%
    Atmospheric
    transmittance
    samedi 26 septembre 20

    View Slide

  45. Synthetic test 2
    • Pure mineral spectra :
    • Grain size factor using Shkuratov theory
    from x1/1000 to x1000
    • Adding instrumental noise from dark
    current
    Schmidt, F.; Legendre, M. & Le Mouëlic, S. Minerals detection for hyperspectral images using adapted
    linear unmixing: LinMin Icarus, 2014, 237, 61-74, http://dx.doi.org/10.1016/j.icarus.2014.03.044
    Shkuratov, Y.; Starukhina, L.; Hoffmann, H. & Arnold, G. A Model of
    Spectral Albedo of Particulate Surfaces: Implications for Optical
    Properties of the Moon Icarus, 1999, 137, 235-246, http://
    www.sciencedirect.com/science/article/B6WGF-45GMFKB-5T/
    2/2b056567d27e74edbaf976c01f89d10f
    samedi 26 septembre 20

    View Slide

  46. 1 1.5 2 2.5
    0.32
    0.34
    0.36
    0.38
    0.4
    0.42
    0.44
    0.46
    Detection of Gypsum (10 %)
    Wavelength
    Reflectance (offset for clarity)
    fact=1000 (ab = 0 %)
    fact=10 (ab = 5.9 %)
    fact=1 (ab = 9.7 %)
    fact=1/10 (ab = 5.1 %)
    fact=1/1000 (ab = 0.5 %)
    Larger
    grain
    size
    Smaller
    grain
    size
    samedi 26 septembre 20

    View Slide

  47. Results on Ophir, Valles
    Marineris
    • Indurated rocks
    made of :
    • Orthopyroxene
    • Olivine
    • Kieserite
    • ultramafic rock,
    altered by aqueous
    processes
    COMPOSITIONAL+MAPPING
    Preliminary results
    samedi 26 septembre 20

    View Slide

  48. More sophisticated estimates through Mixed Integer Programming (MIP)
    y =
    S
    a ( + errors )
    S
    ; high number of endmembers
    Spectral variability [Zare and Ho, 2014, Meyer et al., 2016]
    )
    add more constraints to the problem
    Models based on binary variables encoding the presence of each member in the data
    bn
    = 1 , an
    6= 0
    ; reformulated as 0  an
     bn
    : Mixed Integer Programs
    Sparsity: most abundances are zero [Iordache et al., 2011]
    Structuration of the dictionary into groups [Meyer et al., 2016, Drumetz et al., 2019]
    Minimum values on the nonzero coe cients [never seen . . . ]
    3 / 11
    .1cm
    Collaboration with S. Bourguignon and R. Ben Mhenni
    samedi 26 septembre 20

    View Slide

  49. Exact
    `0
    -norm sparsity
    Only a small number of elementary spectra are used for representing the mixture
    0
    0
    0
    0
    0
    0
    0.2
    0.3
    0.5
    s
    1
    s
    2
    s
    3
    s
    4
    s
    5
    s
    6
    s
    7
    s
    8
    s
    9
    ਚ S a

    ! a must be sparse:
    kak0
    = Card(n|an
    6= 0)  K
    Standard sparse methods perform poorly (`1
    -norm, greedy algorithms)
    ;
    Exact
    `0
    -norm constraint:
    min
    a2[0,1]N , b2{0,1}N
    1
    2
    ky
    S
    ak2 s.t.
    8
    >
    <
    >
    :
    0
     a  b
    P
    N
    n=1
    bn
     K
    P
    N
    n=1
    an
    = 1
    4 / 11
    .1cm
    Collaboration with S. Bourguignon and R. Ben Mhenni
    samedi 26 septembre 20

    View Slide

  50. To sum up: reformulations as Mixed-integer Programs (MIPs)
    `0
    -norm sparsity
    min
    a2[0,1]N , b2{0,1}N
    1
    2
    ky
    S
    ak2 s.t.
    8
    >
    <
    >
    :
    0
     a  b
    P
    N
    n=1
    bn
     K
    P
    N
    n=1
    an
    = 1
    Group Exclusivity
    min
    a2[0,1]N , b2{0,1}N
    1
    2
    ky
    S
    ak2 s.t.
    8
    >
    <
    >
    :
    0
     a  b
    P
    i2Gj
    bn
     1
    P
    N
    n=1
    an
    = 1
    Significant Abundances
    min
    a2[0,1]N , b2{0,1}N
    1
    2
    ky
    S
    ak2 s.t.
    (
    ⌧b  a  b
    P
    N
    n=1
    an
    = 1
    ! E cient resolution via numerical MIP solvers (ex. CPLEX)
    ! These constraints can be mixed
    7 / 1
    .1cm
    Collaboration with S. Bourguignon and R. Ben Mhenni
    SA
    GE
    samedi 26 septembre 20

    View Slide

  51. Ben Mhenni, R.; Bourguignon, S.; Ninin, J. & Schmidt, F. Spectral Unmixing with Sparsity and Structuring Constraints Workshop on Hyperspectral Image and Signal Processing: Evolution in
    Remote Sensing (WHISPERS), 2018 9th Workshop, held 23-26 September 2018, Amsterdam, The Netherlands, 2018, 1-4
    FCLS GE `0
    Backward
    Abundance
    100 200 300 400
    0
    0.1
    0.2
    0.3
    0.4
    100 200 300 400
    0
    0.1
    0.2
    0.3
    0.4
    100 200 300 400
    0
    0.1
    0.2
    0.3
    0.4
    100 200 300 400
    0
    0.1
    0.2
    0.3
    0.4
    Reflectance
    1 1.5 2 2.5
    0
    0.1
    0.2
    1 1.5 2 2.5
    0
    0.1
    0.2
    1 1.5 2 2.5
    0
    0.1
    0.2
    1 1.5 2 2.5
    0
    0.1
    0.2
    (µm) (µm) (µm) (µm)
    Fig. 1. Example of result, K
    =
    4 spectra, SNR = 40 dB. Top: estimated (black squares) and true (red stars) abundances. Bottom:
    true (solid red line) and estimated (dashed black line) endmembers weighted by their abundances, and noise (blue solid line).
    k
    b
    a
    ˚
    a
    k
    2
    1 2 3 4 5 6 7
    0
    0.1
    0.2
    0.3
    0.4
    1 2 3 4 5 6 7
    0
    0.1
    0.2
    0.3
    0.4
    Fig. 2. Error on estimated abundances for FCLS (+), GE (4),
    `0
    (o), SA+GE (

    ), and backward [12] (x).
    lowed for MIP resolution (1000 s).
    SNR Method K=3 K=4 K=5 K=6 K=7
    55 dB GE 1.3 2.2 3.3 5.9 9.8
    `0
    1.5 3 11.8 124(3) 468(9)
    `0
    +GE 1.6 3.1 7.1 82(1) 320(7)
    SA+GE 2.5 11.5 29 76 264(5)
    40 dB GE 1.3 2.2 2.1 3.2 2.9
    `0
    9.2 107(1) 574(13) 963(27) 1000(30)
    `0
    +GE 8.1 75 515(11) 929(25) 982(28)
    SA+GE 282(6) 426(8) 746(19) 936(26) 997(29)
    Table 1. Computing times (s) for optimization of MIP prob-
    lems, averaged over 30 instances. In parentheses: number of
    Abu
    100 200 300 400
    0
    0.1
    0.2
    100 200 300 400
    0
    0.1
    0.2
    100 200
    0
    0.1
    0.2
    Reflectance
    1 1.5 2 2.5
    0
    0.1
    0.2
    1 1.5 2 2.5
    0
    0.1
    0.2
    1 1.5
    0
    0.1
    0.2
    (µm) (µm) (µm
    Fig. 1. Example of result, K
    =
    4 spectra, SNR = 40 dB. Top: estimated (black s
    true (solid red line) and estimated (dashed black line) endmembers weighted b
    k
    b
    a
    ˚
    a
    k
    2
    1 2 3 4 5 6 7
    0
    0.1
    0.2
    0.3
    0.4
    1 2 3 4 5 6 7
    0
    0.1
    0.2
    0.3
    0.4
    Fig. 2. Error on estimated abundances for FCLS (+), GE (4),
    `0
    (o), SA+GE (

    ), and backward [12] (x).
    lowed for MIP resolution (1000 s).
    Table 1 gives average computing times obtained on
    SNR M
    55 dB
    `
    S
    40 dB
    `
    S
    Table 1. Com
    lems, averag
    instances for
    •SA : significant abundance
    •GE : group exclusivity
    •Backward : J. B. Greer, “Sparse
    demixing of hyperspectral images,” IEEE
    Transactions on Image Processing, vol. 21,
    no. 1, pp. 219–228, Jan 2012.
    Average of 50
    realizations
    samedi 26 septembre 20

    View Slide

  52. Plan
    • Linear mixing is useful to handle complex
    dataset !
    • for characterizing unexpected spectra
    • to pick/detect out of a large dataset
    (Interior-Point VS Mixed Integer Program)
    • Next generation of algorithm (non-linear?!) ?
    New discoveries in planetary science...
    F. Schmidt Conclusion
    samedi 26 septembre 20

    View Slide

  53. So many fruitful
    collaboration...
    • Jennifer Fernando, François Andrieu, Ines Belgacem,
    Guillaume Cruz-Mermy, Chiara Marmo...
    • Matthieu Kowalski, Nicolas Gac, Alina Meresescu...
    • Saïd Moussaoui, Sébastien Bourguignon, Ramzi Ben
    Mhenni, Maxime Legendre...
    • Christian Jutten, Jocelyn Chanussot, Sylvain Douté...
    F. Schmidt Conclusion
    samedi 26 septembre 20

    View Slide