I have a time series dataset pr11 (shape is (151,)) which looks like the graph below when plotted. Note the very small numbers. I want to find the average slope of the data by doing a least square fit to a straight line.
I've tried two different methods from another StackExchange page to get the answer. I tried using scipy.optimize.curve_fit as below...
len = np.arange(pr11.shape[0])
def f(x, A, B):
return A*x + B
A,B = curve_fit(f,pr11,len)[0]
However, this gives me a slope (A) of 1.0, which I know isn't right, so something must be off here. The "fitted" data just ends up looking the exact same as my original data. I also tried scipy.stats...
slope, intercept, r_value, p_value, std_err = stats.linregress(len,pr11)
My slope this time is a number on the order of e-08. The issue with that is when I use the equation for a line slope*x + intercept, that number multiplies my time series data to a very low value (order e-15). Therefore when I plot the fitted line, the line is horizontal and doesn't fit my data at all.
How can I get a a fitted line for this data?
A package that I like to use for fitting is lmfit
. Once you install it, you can do:
from lmfit import minimize, Parameters, Parameter, report_fit
import numpy as np
# create data to be fitted
x = np.arange(150)/100.
data = 2e-6*x-5e-7 + np.random.normal(size=len(x), scale=5e-7)
# define objective function: returns the array to be minimized
def fcn2min(params, x, data):
""" model decaying sine wave, subtract data"""
slope = params['slope'].value
offset = params['offset'].value
model = slope * x + offset
return model - data
# create a set of Parameters
params = Parameters()
params.add('slope', value= 1., min=0)
params.add('offset', value= 0.)
# do fit, here with leastsq model
result = minimize(fcn2min, params, args=(x, data))
# calculate final result
final = data + result.residual
# write error report
report_fit(result.params)
# [[Variables]]
# slope: 2.1354e-06 +/- 9.33e-08 (4.37%) (init= 1)
# offset: -6.0680e-07 +/- 8.02e-08 (13.22%) (init= 0)
# [[Correlations]] (unreported correlations are < 0.100)
# C(slope, offset) = -0.865
# plot results
import matplotlib.pyplot as plt
plt.plot(x, data, 'k+')
plt.plot(x, final, 'r')
plt.show()