Help with a simple Gaussian fitting

[Info Center]

Dear all,

I have an absorption line data which  I thought could be easily fitted with a simple  1D Gaussian function. For some reason, the lmfit does not seem to give any meaningful result. The initial values for the fit are being guessed from the data, but I feel that the guess that lmfit makes is incorrect. Here, the question that really puzzles me is when is it that lmfit can make a correct initial guess and in what circumstances it fails?   I also tried setting the parameters manually closer to expected values but that too seems to fail. But that again begs the question as to how close should the initial guess be ?     

I am attaching my python code as well as the figure and the data used in this example.I'll really appreciate it if someone could tell me what I am doing wrong.... 

#--------------------------------------------------------------

from lmfit.models import GaussianModel 

import numpy as np

data = np.loadtxt('test.dat')  # Load the test data 

xx = data[:, 0]                 # read the x-column

yy = data[:, 1]                 # read the y-columns

model = GaussianModel ()        # initiate the Gaussian Model

params = model.guess(yy, x=xx)  # guess the initial parameters

result = model.fit(yy, params, x=xx) # fit a 1D Gaussian to data  

result.plot_fit()                    # display the data and the fit

#----------------------------------------------------------------------------------------------------

print (params)

Parameters([('amplitude', <Parameter 'amplitude', value=27.216, bounds=[-inf:inf]>), ('center', <Parameter 'center', value=5141.315789473684, bounds=[-inf:inf]>), ('sigma', <Parameter 'sigma', value=15.0, bounds=[0.0:inf]>), ('fwhm', <Parameter 'fwhm', value=35.3223, bounds=[-inf:inf], expr='2.3548200*sigma'>), ('height', <Parameter 'height', value=0.7238409091200001, bounds=[-inf:inf], expr='0.3989423*amplitude/max(2.220446049250313e-16, sigma)'>)])

-------------------------------------------------------------------------------------------------------------------

Thanks 

With Best Wishes

Steve

[Matt Newville]

On Tue, Jun 15, 2021 at 10:45 AM Info Center <info...@gmail.com> wrote:

Dear all,

I have an absorption line data which  I thought could be easily fitted with a simple  1D Gaussian function. For some reason, the lmfit does not seem to give any meaningful result. The initial values for the fit are being guessed from the data, but I feel that the guess that lmfit makes is incorrect. Here, the question that really puzzles me is when is it that lmfit can make a correct initial guess and in what circumstances it fails?   I also tried setting the parameters manually closer to expected values but that too seems to fail. But that again begs the question as to how close should the initial guess be ?     

Generally speaking, the initial values should be reasonably close, but not necessarily perfect.  For a simple Gaussian, one might say that the guess for centroid should be within 3 or 4 sigmas.  Guessing sigma larger than the expected value is probably a good idea -- even up to 10x the expected sigma should be fine.  The value for amplitude should be within a factor of 1000 or so.  

For Gaussians with a center of 1000 and a sigma of 25, the real danger happens if the initial guess of the center is -500.   Basically, guessing the center to be the position of maximum intensity is probably suitable.

The problem is that the data you show is not Gaussian.  In lmfit, a Gaussian function goes to 0 far away from the peak.  Your data varies from 0.3 to 0.9.  I would guess that your intention is that there is an "offset" around 1 and a Gaussian that has a negative amplitude, dropping to 0.3 at around 5140.  It is also possible that there are 2 peaks, but I guess that is not how you are looking at the data.

If that guess is correct, what you would want to do is model these data as ConstantModel (or perhaps a LinearModel?) to model the background plus a GaussianModel, perhaps something like

from lmfit.models import GaussianModel, ConstantModel

import numpy as np

data = np.loadtxt('test.dat')  # Load the test data
xx = data[:, 0]                 # read the x-column
yy = data[:, 1]                 # read the y-columns

g_model = GaussianModel()
params = g_model.guess(yy, x=xx)

model = g_model + ConstantModel()
params.add('c', value=0.9)

# since the Gaussian has a negative amplitude, make it sure it starts out negative

# and maybe even place an upper bound of 0 (probably not necessary here...)

params['amplitude'].value = -10
params['amplitude'].max = 0

result = model.fit(yy, params, x=xx) #

result.plot_fit()                    
print(result.fit_report())

--Matt

I am attaching my python code as well as the figure and the data used in this example.I'll really appreciate it if someone could tell me what I am doing wrong.... 

#--------------------------------------------------------------

from lmfit.models import GaussianModel 

import numpy as np

data = np.loadtxt('test.dat')  # Load the test data 

xx = data[:, 0]                 # read the x-column

yy = data[:, 1]                 # read the y-columns

model = GaussianModel ()        # initiate the Gaussian Model

params = model.guess(yy, x=xx)  # guess the initial parameters

result = model.fit(yy, params, x=xx) # fit a 1D Gaussian to data  

result.plot_fit()                    # display the data and the fit

#----------------------------------------------------------------------------------------------------

print (params)

Parameters([('amplitude', <Parameter 'amplitude', value=27.216, bounds=[-inf:inf]>), ('center', <Parameter 'center', value=5141.315789473684, bounds=[-inf:inf]>), ('sigma', <Parameter 'sigma', value=15.0, bounds=[0.0:inf]>), ('fwhm', <Parameter 'fwhm', value=35.3223, bounds=[-inf:inf], expr='2.3548200*sigma'>), ('height', <Parameter 'height', value=0.7238409091200001, bounds=[-inf:inf], expr='0.3989423*amplitude/max(2.220446049250313e-16, sigma)'>)])

-------------------------------------------------------------------------------------------------------------------

Thanks 

With Best Wishes

Steve

--
You received this message because you are subscribed to the Google Groups "lmfit-py" group.
To unsubscribe from this group and stop receiving emails from it, send an email to lmfit-py+u...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/lmfit-py/CAKjCJTqiXgja38yLTXx5Fbm7mxHFW6ye_6VQJhz_LJDe6%2BEeAA%40mail.gmail.com.

[Info Center]

Dear Matt,

Thanks for explaining it so clearly and suggesting an excellent  way to solve the  problem. In fact, it never occured to me that the shallow tail in data was responsible for Gaussian fit not performing  as expected. The ConstantModel is indeed a great hack to deal with real data.  Once again, thank you so much for pitching in and educating !!

To view this discussion on the web visit https://groups.google.com/d/msgid/lmfit-py/CA%2B7ESbrZmwvhSHzWCtHir%3D_E0uOmuceorxUWXn7k_eu9Z4c-dw%40mail.gmail.com.

[Thomas Bersano]

As Matt points out, your data is not Gaussian. Even with the constant offset you may still have issues because spectral lines tend not to be Gaussian in shape. I suggest you consider another distribution .e.g Lorentzian for a proper fit to your data.

https://en.wikipedia.org/wiki/Spectral_line_shape

Best,

Thomas