usnistgov / pyMCR

pyMCR: Multivariate Curve Resolution for Python
https://pages.nist.gov/pyMCR
Other
80 stars 27 forks source link

Unimodality constraint #38

Open Bio-data-tricks opened 3 years ago

Bio-data-tricks commented 3 years ago

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 commented 3 years ago

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 commented 3 years ago

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 notifications@github.com Envoyé : vendredi 5 février 2021 16:17 À : usnistgov/pyMCR pyMCR@noreply.github.com Cc : colapili colapili@hotmail.fr; Author author@noreply.github.com 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 GitHubhttps://github.com/usnistgov/pyMCR/issues/38#issuecomment-774095557, or unsubscribehttps://github.com/notifications/unsubscribe-auth/AGENLS2JPJNSW2TQBKMY7GTS5QD2BANCNFSM4XENIP3Q.

Bio-data-tricks commented 3 years ago

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 colapili@hotmail.fr Envoyé : samedi 6 février 2021 10:30 À : usnistgov/pyMCR reply@reply.github.com 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 notifications@github.com Envoyé : vendredi 5 février 2021 16:17 À : usnistgov/pyMCR pyMCR@noreply.github.com Cc : colapili colapili@hotmail.fr; Author author@noreply.github.com 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 GitHubhttps://github.com/usnistgov/pyMCR/issues/38#issuecomment-774095557, or unsubscribehttps://github.com/notifications/unsubscribe-auth/AGENLS2JPJNSW2TQBKMY7GTS5QD2BANCNFSM4XENIP3Q.

CCampJr commented 3 years ago

@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 commented 3 years ago

@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)