python uncertainty package: error propagation in curve_fit

Question

I would like to use the uncertainty package in python to propagate errors in a Gaussian fitting.

I have experimental data, x and y, and I have uncertainty on both of them. The y is approximated to follow a Gaussian relationship with x, so I fit it with a Gaussian curve. I would like to know the uncertainties of the fitting parameters.

I used the following code:

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
import uncertainties as un

def gaus_fit(x,a,x0,sigma):
    return a*np.exp(-(x-x0)**2/(2*sigma**2*sigma/abs(sigma)))

def curve_fit_float(x_arr,y_arr,p0):
    popt_X,pcov_X = curve_fit(gaus_fit,x_arr,y_arr,p0)
    return popt_X[2]

arr_x=np.array([    0.        ,     0.        ,     0.        ,     0.        ,
           0.        ,     0.        ,     0.        ,     0.        ,
           0.        ,     0.        ,     0.        ,     0.        ,
           0.        ,     0.        ,     0.        ,     0.        ,
           0.        ,     0.        ,     0.        ,     0.        ,
           0.        ,     0.        ,     0.        ,     0.        ,
           0.        ,   263.46478873,   380.        ,   614.74647887,
         946.07042254,  1838.18309859,  3226.18309859,  5361.07042254,
        7910.49295775, 10306.29577465, 12224.07042254, 12463.97183099,
       11558.87323944,  9225.8028169 ,  6181.09859155,  3764.05633803,
        2075.78873239,   983.43661972,   479.70422535,   218.84507042,
         118.09859155,    73.53521127,    61.29577465,    51.85915493,
          18.63380282,    22.77464789,    18.76056338,    20.84507042,
           0.        ,     0.        ,     0.        ,     0.        ,
           0.        ,     0.        ,     0.        ,     0.        ,
           0.        ,     0.        ,     0.        ,     0.        ,
           0.        ,     0.        ,     0.        ,     0.        ,
           0.        ,     0.        ,     0.        ])

arr_x_coor=np.array([-18.64222222, -18.08666667, -17.53111111, -16.97555556,
       -16.42      , -15.86444444, -15.30888889, -14.75333333,
       -14.19777778, -13.64222222, -13.08666667, -12.53111111,
       -11.97555556, -11.42      , -10.86444444, -10.30888889,
        -9.75333333,  -9.19777778,  -8.64222222,  -8.08666667,
        -7.53111111,  -6.97555556,  -6.42      ,  -5.86444444,
        -5.30888889,  -4.75333333,  -4.19777778,  -3.64222222,
        -3.08666667,  -2.53111111,  -1.97555556,  -1.42      ,
        -0.86444444,  -0.30888889,   0.24666667,   0.80222222,
         1.35777778,   1.91333333,   2.46888889,   3.02444444,
         3.58      ,   4.13555556,   4.69111111,   5.24666667,
         5.80222222,   6.35777778,   6.91333333,   7.46888889,
         8.02444444,   8.58      ,   9.13555556,   9.69111111,
        10.24666667,  10.80222222,  11.35777778,  11.91333333,
        12.46888889,  13.02444444,  13.58      ,  14.13555556,
        14.69111111,  15.24666667,  15.80222222,  16.35777778,
        16.91333333,  17.46888889,  18.02444444,  18.58      ,
        19.13555556,  19.69111111,  20.24666667])

arr_x_coor_uncertainty=[0.]*71
arr_x_uncertainty=[0.]*71

for ii in np.arange(0,71):
    arr_x_coor_uncertainty[ii]=un.ufloat(arr_x_coor[ii],0.555/2)
    arr_x_uncertainty[ii]=un.ufloat(arr_x[ii],np.sqrt(arr_x[ii]))

wrapped_fit=un.wrap(curve_fit_float)

# fit curve with a Gaussian
mean_X=(arr_x_coor[np.size(arr_x_coor)-1]+arr_x_coor[0])/2
sigma_X=1.0
wrapped_fit(arr_x_coor_uncertainty,arr_x_uncertainty,[max(arr_x),mean_X,sigma_X])

I recrived the following error:


  File "/opt/anaconda3/lib/python3.8/site-packages/uncertainties/core.py", line 669, in f_with_affine_output
    return f(*args, **kwargs)

  File "/Users/huangze/PhD/academics/post_PSI_data_analysis/debug_odr_uncertainty_package.py", line 18, in curve_fit_float
    popt_X,pcov_X = curve_fit(gaus_fit,x_arr,y_arr,p0)

  File "/opt/anaconda3/lib/python3.8/site-packages/scipy/optimize/minpack.py", line 734, in curve_fit
    ydata = np.asarray_chkfinite(ydata, float)

  File "/opt/anaconda3/lib/python3.8/site-packages/numpy/lib/function_base.py", line 486, in asarray_chkfinite
    a = asarray(a, dtype=dtype, order=order)

  File "/opt/anaconda3/lib/python3.8/site-packages/numpy/core/_asarray.py", line 102, in asarray
    return array(a, dtype, copy=False, order=order)

  File "/opt/anaconda3/lib/python3.8/site-packages/uncertainties/core.py", line 2712, in raise_error
    raise TypeError("can't convert an affine function (%s)"

TypeError: can't convert an affine function (<class 'uncertainties.core.Variable'>) to float; use x.nominal_value

Where should I use the nominal_value?

Can I use the uncertainty package with curve_fit?

I really appreciate any help.

score 0 · Answer 1 · answered Jun 25 '22 at 17:14

You're getting this error because scipy is expecting floats, not ufloat. Remove the loop where you create a bunch of ufloats and pass your stats (note:look up uncertainties.unumpy). Typically, errors in regressions are considered by weighing the data points by their errors. Curve_fit has that option. The fit parameter errors are obtained by the square root of the diagonal of the covariance matrix pcov. In summary, you don't use uncertainties with curve_fit

python uncertainty package: error propagation in curve_fit

1 Answers1