Unimodality constraint

[Bio-data-tricks]

Hello,
I would like to use a Unimodality constraint to follow a kinetic process (time-ordered matrix ). Is it possible to perform this kind of study with your package ?

Thank you

Nicolas

[CCampJr]

Hi Nicolas,

Yes, you can do unimodality with a little bit of coding. In your case, does the mode need to follow an analytical shape or is it just enforcing that the global and local maximum are the same?

For enforcing an analytical modal shape, here's an example: https://pages.nist.gov/pyMCR/usage.html#single-gaussian-distribution-spectral-fitting-with-lmfit

For just enforcing unimodality with no particular shape requirements, I could code up something quickly that you could copy-and-paste until I can push it to the full package. (Unimodality used to be in an older version of the code, but I can't find it -- d'oh!).

Thanks!

[Bio-data-tricks]

Hi Charles, Thank you very much for your quick response ! To explain a little more, I study a process of adipocytes differentiation by Raman spectroscopy. So I have a well known initial state (protein spectrum) and I try to visualize the appearance of certain spectral components along the differentiation process (regarding of time). I'm searching to visualize a sequential apparition (or disparition) of signals, so i don't want to use the package limfit (why not if i want to study a specific vibration after this global study) I have also heard that when studying an ordered matrix it may be better to use evolving factor analysis (EFA) rather than simplisma for initialization, could you advise me on this subject? Regards Nicolas

…

________________________________ De : Charles <[email protected]> Envoyé : vendredi 5 février 2021 16:17 À : usnistgov/pyMCR <[email protected]> Cc : colapili <[email protected]>; Author <[email protected]> Objet : Re: [usnistgov/pyMCR] Unimodality constraint (#38) Hi Nicolas, Yes, you can do unimodality with a little bit of coding. In your case, does the mode need to follow an analytical shape or is it just enforcing that the global and local maximum are the same? For enforcing an analytical modal shape, here's an example: https://pages.nist.gov/pyMCR/usage.html#single-gaussian-distribution-spectral-fitting-with-lmfit For just enforcing unimodality with no particular shape requirements, I could code up something quickly that you could copy-and-paste until I can push it to the full package. (Unimodality used to be in an older version of the code, but I can't find it -- d'oh!). Thanks! — You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub<#38 (comment)>, or unsubscribe<https://github.com/notifications/unsubscribe-auth/AGENLS2JPJNSW2TQBKMY7GTS5QD2BANCNFSM4XENIP3Q>.

[Bio-data-tricks]

Hi Charles, I have the impression that optimized spectra I'm getting with MCR-ALS are a bit noisy. I don't feel that my data are so noisy. Which parameter could influence this result ? Thank you Nicolas

…

________________________________ De : Nicolas Goffin <[email protected]> Envoyé : samedi 6 février 2021 10:30 À : usnistgov/pyMCR <[email protected]> Objet : RE: [usnistgov/pyMCR] Unimodality constraint (#38) Hi Charles, Thank you very much for your quick response ! To explain a little more, I study a process of adipocytes differentiation by Raman spectroscopy. So I have a well known initial state (protein spectrum) and I try to visualize the appearance of certain spectral components along the differentiation process (regarding of time). I'm searching to visualize a sequential apparition (or disparition) of signals, so i don't want to use the package limfit (why not if i want to study a specific vibration after this global study) I have also heard that when studying an ordered matrix it may be better to use evolving factor analysis (EFA) rather than simplisma for initialization, could you advise me on this subject? Regards Nicolas

________________________________ De : Charles <[email protected]> Envoyé : vendredi 5 février 2021 16:17 À : usnistgov/pyMCR <[email protected]> Cc : colapili <[email protected]>; Author <[email protected]> Objet : Re: [usnistgov/pyMCR] Unimodality constraint (#38) Hi Nicolas, Yes, you can do unimodality with a little bit of coding. In your case, does the mode need to follow an analytical shape or is it just enforcing that the global and local maximum are the same? For enforcing an analytical modal shape, here's an example: https://pages.nist.gov/pyMCR/usage.html#single-gaussian-distribution-spectral-fitting-with-lmfit For just enforcing unimodality with no particular shape requirements, I could code up something quickly that you could copy-and-paste until I can push it to the full package. (Unimodality used to be in an older version of the code, but I can't find it -- d'oh!). Thanks! — You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub<#38 (comment)>, or unsubscribe<https://github.com/notifications/unsubscribe-auth/AGENLS2JPJNSW2TQBKMY7GTS5QD2BANCNFSM4XENIP3Q>.

[CCampJr]

@colapili

Hi Nicolas,

That's a great question and one that I don't think has a singular answer. I have seen it result from preprocessing steps, I've seen it when some of the data are lots of 0's (exactly 0), and it can depend on the type of regressor used.

The first thing I would suspect is that regular/ordinary least-squares (OLS) and non-negative least squares (NNLS) are most likely to cause a noise-like appearance. This is because they will tend to give some weight (regression coefficients) to each and every feature. I would try using Ridge regression (https://pages.nist.gov/pyMCR/usage.html#ridge-regression-from-scikit-learn) with a small weight (alpha) to start -- maybe 1e-5 and continue to lower it by order-of-magnitude until the noise returns. You may even have to raise the alpha to 1 or 10, but I would start small. If Ridge regression is going to converge to OLS, it's usually at 1e-8 to 1e-10 in my experience -- though it doesn't always give the same answer (a case of numerical vs theoretical "should").

Just a side-note: I have found the Ridge SVD solver to be more consistent if using a LOT of data -- maybe your data isn't that big.

Ridge(alpha=1e-5, solver='svd')

Related to that, I can't tell you whether to make the S or the C regressors (or both) Ridge -- just going to need to test it out.

I hope that helps!

Charlie

[CCampJr]

@colapili

Here is an example function that can be used to make a unimodal constraint (see if you can get it going using this example: https://pages.nist.gov/pyMCR/usage.html#single-gaussian-distribution-spectral-fitting-with-lmfit as a template)

import numpy as np

import matplotlib.pyplot as plt
%matplotlib inline

x = np.arange(1001)-500
y = np.exp(-x**2/20**2) + 1.2*np.exp(-(x-100)**2/200**2)
y = np.vstack((y, 2*y))

def unimodal(signal, axis=-1):
    """ Enforces unimodality assuming the peak is upwards"""
    
    assert axis == -1  # Too lazy to implement other currently
    
    if signal.ndim == 1:
        signal = signal[None,:]
        
    loc_max = np.argmax(signal, axis=-1)
    
    # Gradients
    dx_signal_l = np.gradient(signal, axis=axis)  # Left-side of peak
    dx_signal_r = 1*dx_signal_l  # Right-side of peak
    
    for num, lm in enumerate(loc_max):
        dx_signal_l[num,lm:] = 0
        dx_signal_r[num,:lm] = 0
    
        dx_signal_l[num, dx_signal_l[num,:]<0] = 0
        dx_signal_r[num, dx_signal_r[num,:]>0] = 0
    
    return np.cumsum(dx_signal_l + dx_signal_r, axis=-1)
    
    
plt.plot(x,y.T)
plt.plot(x,unimodal(y).T)