pythonscipyleast-squaresnon-linear-regressionrobust

How to get a robust nonlinear regression fit using scipy.optimize.least_squares?


My specific issue is that I cannot seem to get my data to converted to floating points. I have data and simply want to fit a robust curve using my model equation:

y = a * e^(-b*z)

This cookbook is my reference: click

Below is my attempt. I am getting this:

TypeError: 'data type not understood'

which I believe is because my columns are strings, so I tried pd.Series.astype().

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import least_squares

for i in range(1):

    def model(z, a, b):
        y = a * np.exp(-b * z)
        return y
    data = pd.read_excel('{}.xlsx'.format(600+i), names = ['EdGnd','380','395','412','443','465','490','510','520','532','555','560','565','589','625','665','670','683','694','710','Temp','z','EdZTemp','Tilt','Roll','EdZVin'])
    data.dropna(axis = 0, how = 'any')
    data.astype('float')
    np.dtype(data)
    data.plot.scatter('z','380') 
    def fun(x, z, y):
        return x[0] * np.exp(-x[1] * z) - y
    x0 = np.ones(3)
    rbst1 = least_squares(fun, x0, loss='soft_l1', f_scale=0.1, args=('z', 'ed380'))
    y_robust = model('z', *rbst1.x)
    plt.plot('z', y_robust, label='robust lsq')
    plt.xlabel('$z$')
    plt.ylabel('$Ed$')
    plt.legend();

Solution

  • I think the problem is that you pass 'z' in args which is a string and can therefore not be used in the multiplication.

    Below is some code using curve_fit which uses least_squares but might be slightly easier to use:

    import matplotlib.pyplot as plt
    import numpy as np
    from scipy.optimize import curve_fit
    
    
    # your model definition
    def model(z, a, b):
        return a * np.exp(-b * z)
    
    # your input data
    x = np.array([20, 30, 40, 50, 60])
    y = np.array([5.4, 4.0, 3.0, 2.2, 1.6])
    
    # do the fit with some initial values
    popt, pcov = curve_fit(model, x, y, p0=(5, 0.1))
    
    # prepare some data for a plot
    xx = np.linspace(20, 60, 1000)
    yy = model(xx, *popt)
    
    plt.plot(x, y, 'o', xx, yy)
    plt.title('Exponential Fit')
    
    plt.show()
    

    This will plot

    enter image description here

    You could try to adapt this code for your needs.

    If you want to use f_scale you can use:

    popt, pcov = curve_fit(model, x, y, p0=(5, 0.1), method='trf', f_scale=0.1)
    

    See the documentation:

    kwargs

    Keyword arguments passed to leastsq for method='lm' or least_squares otherwise.

    If you have an unbound problem, by default method='lm' is used which uses leastsq which does not accept f_scale as a keyword. Therefore, we can use method='trf' which then uses least_squares which accepts f_scale.