I am writing a program in Python that will fit Gaussian and Lorentzian shapes to some given resonance data. I originally began using scipy.optimize.leastsq
but changed to using optimize.curve_fit
after having difficulties in retrieving the errors in the optimized parameters from the covariance matrix.
I have defined a function to fit a sum of Gaussian and Lorentzian:
def mix(x,*p):
ng = numg
p1 = p[:3*ng]
p2 = p[3*ng:]
a = sumarray(gaussian(x,p1),lorentzian(x,p2))
return a
where p
is an array of the initial guesses at the fit parameters. Here is the instance where it is called using curve_fit
:
leastsq,covar = opt.curve_fit(mix,energy,intensity,inputtot)
At the moment numg
(the number of Gaussian shapes) is a global variable. Is there's any way that it can be incorporated into curve_fit
as an extra argument instead, as can be done with leastsq
?
The Python SciPy Optimize Curve Fit function is widely used to obtain the best-fit parameters. The curve_fit() function is an optimization function that is used to find the optimized parameter set for a stated function that perfectly fits the provided data set.
What does popt and pcov mean? popt- An array of optimal values for the parameters which minimizes the sum of squares of residuals. pcov-2d array which contains the estimated covariance of popt. The diagonals provide the variance of the parameter estimate.
The SciPy API provides a 'curve_fit' function in its optimization library to fit the data with a given function. This method applies non-linear least squares to fit the data and extract the optimal parameters out of it.
The returned parameter covariance matrix pcov is based on scaling sigma by a constant factor. This constant is set by demanding that the reduced chisq for the optimal parameters popt when using the scaled sigma equals unity. In other words, sigma is scaled to match the sample variance of the residuals after the fit.
The great thing about python is that you can define functions that return other functions, try currying:
def make_mix(numg):
def mix(x, *p):
ng = numg
p1 = p[:3*ng]
p2 = p[3*ng:]
a = sumarray(gaussian(x,p1),lorentzian(x,p2))
return a
return mix
and then
leastsq, covar = opt.curve_fit(make_mix(numg),energy,intensity,inputtot)
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With