<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Intro 
Step1:  LMFIT package
Step2:  Lets visualize how a quadratic curve fits to it
Step3:  Lets build a general quadratic model
Step4:  Questions ?
leastsq function from scipy
Step5:  LMFIT 
Step6:  Fit values are the same as before !
Step7:  Questions ?
Manipulating parameters 
Step8:  Another way of defining the parameters 
Step9:  Other manipulations 
Step10:  Challenge 
Step11:  stderr
Step12:  correl
Step13:  report_fit
Step14:  Choosing Different Fitting Methods 
Step15:  Challenge 
Step16:  Complete report 
Step17:  Using expressions
<ASSISTANT_TASK:>
Python Code:
from IPython.display import Image
Image(filename='minimize.png', width=500, height=500) 
Explanation: Intro :
End of explanation
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
x= np.array([0.,1.,2.,3.])
data = np.array([1.3,1.8,5.,10.7])
Explanation: LMFIT package:
https://lmfit.github.io/lmfit-py/index.html
Example :
End of explanation
plt.scatter(x,data)
xarray=np.arange(-1,4,0.1)
plt.plot(xarray, xarray**2,'r-')  # Not the best fit
plt.plot(xarray, xarray**2+1,'g-')
Explanation: Lets visualize how a quadratic curve fits to it:
End of explanation
def get_residual(vars,x, data):
    a= vars[0]
    b=vars[1]
    model =a* x**2 +b
    return data-model
vars=[1.,0.]
print get_residual(vars,x,data)
print sum(get_residual(vars,x,data))
vars=[1.,1.]
print sum(get_residual(vars,x,data))
vars=[2.,0.]
print sum(get_residual(vars,x,data))
Explanation: Lets build a general quadratic model:
End of explanation
from scipy.optimize import leastsq
vars = [0.,0.]
out = leastsq(get_residual, vars, args=(x, data))
print out
vars=[1.06734694,  0.96428571]
print sum(get_residual(vars,x,data)**2)
vars=[1.06734694,  0.96428571]
plt.scatter(x,data)
xarray=np.arange(-1,4,0.1)
plt.plot(xarray, xarray**2,'r-')
plt.plot(xarray, xarray**2+1,'g-')
fitted = vars[0]* xarray**2 +vars[1]
plt.plot(xarray, fitted,'b-')
Explanation: Questions ?
leastsq function from scipy:
End of explanation
from lmfit import minimize, Parameters
params = Parameters()
params.add('amp', value=0.)
params.add('offset', value=0.)
def get_residual(params,x, data):
    amp= params['amp'].value
    offset=params['offset'].value
    
    model =amp* x**2 +offset
    
    return data-model
out = minimize(get_residual, params, args=(x, data))
dir(out)
out.params 
dir(out.params)
out.params.values
Explanation: LMFIT :
Using Parameter objects instead of plain floats as variables. A parameter value:
can be varied in the fit
have a fixed value
have upper and/or lower bounds
constrained by an algebraic expression of other Parameter values
Ease of changing fitting algorithms. Once a fitting model is set up, one can change the fitting algorithm used to find the optimal solution without changing the objective function.
Improved estimation of confidence intervals. While scipy.optimize.leastsq() will automatically calculate uncertainties and correlations from the covariance matrix, the accuracy of these estimates are often questionable. To help address this, lmfit has functions to explicitly explore parameter space to determine confidence levels even for the most difficult cases.
Improved curve-fitting with the Model class. This extends the capabilities of scipy.optimize.curve_fit(), allowing you to turn a function that models for your data into a python class that helps you parametrize and fit data with that model.
Many pre-built models for common lineshapes are included and ready to use.
minimize & Parameters:
End of explanation
out.__dict__
out.params['amp'].__dict__ 
Explanation: Fit values are the same as before !
End of explanation
params['amp'].vary = False
out = minimize(get_residual, params, args=(x, data))
print out.params
print out.chisqr
params['amp'].value = 1.0673469387778385
out = minimize(get_residual, params, args=(x, data))
print out.chisqr
Explanation: Questions ?
Manipulating parameters :
parameter class gives a lot of flexibility in manipulating the model parameters !
End of explanation
def get_residual(params,x, data):
    #amp= params['amp'].value
    #offset=params['offset'].value
    #xoffset=params['xoffset'].value    
    
    parvals = params.valuesdict()
    amp = parvals['amp']
    offset = parvals['offset']
    
    model =amp* x**2 +offset
    
    return data-model
Explanation: Another way of defining the parameters :
End of explanation
params = Parameters()
params.add('amp', value=0.)
#params['amp'] = Parameter(value=..., min=...)
params.add('offset', value=0.)
params.add('xoffset', value=0.0, vary=False)
out = minimize(get_residual, params, args=(x, data))
print out.params
Image(filename='output.png', width=500, height=500) 
Explanation: Other manipulations :
End of explanation
params['offset'].min = -10.
params['offset'].max = 10.
out = minimize(get_residual, params, args=(x, data))
print out.params
Explanation: Challenge : Set parameter bound for 'amp' using "min" and "max"
End of explanation
print out.params['amp'].stderr
Explanation: stderr:
End of explanation
print out.params['amp'].correl
Explanation: correl:
End of explanation
from lmfit import minimize, Parameters, Parameter, report_fit
result = minimize(get_residual, params, args=(x, data))
help(report_fit)
# write error report
report_fit(result.params)
Explanation: report_fit: For a better report :
End of explanation
Image(filename='fitting.png', width=500, height=500) 
Explanation: Choosing Different Fitting Methods :
End of explanation
result2 = minimize(get_residual, params, args=(x, data), method='tnc')
report_fit(result2.params)
result3 = minimize(get_residual, params, args=(x, data), method='powell')
report_fit(result3.params)
Explanation: Challenge : Run with two other methods e.g., 'tnc' and 'powell' and compare the results:
End of explanation
print(report_fit(result3))
Explanation: Complete report :
End of explanation
params.add('amp2',     expr='(amp-offset)**2')
def get_residual(params,x, data):
    #amp= params['amp'].value
    #offset=params['offset'].value
    #xoffset=params['xoffset'].value    
    
    parvals = params.valuesdict()
    amp = parvals['amp']
    offset = parvals['offset']
    amp2 = parvals['amp2']
    
    model =amp* x**2 + amp2*x**4 +offset
    
    return data-model
result4 = minimize(get_residual, params, args=(x, data))
report_fit(result4.params)
Explanation: Using expressions:
End of explanation
<END_TASK>