Python曲线拟合
标签:
python |
在进行数据处理时经常会用到曲线的拟合。Matlab和Python都可以完成这项任务。在介绍Python曲线拟合之前先说一下使用Matlab。它对曲线的拟合可以使用cftool(curve fitting tool)曲线拟合工具箱。比如我们有了一组数据(x,y)在命令窗口输入cftool调用曲线拟合工具箱。在以下窗口中选择要拟合的数据。
http://s2/mw690/002kh4Qkzy76NIKbo1r01&690
在下图在红色方框中第一个下拉菜单可以选择使用拟合的方法。其中包括自定义方程、指数拟合、插值拟合、多项式拟合、高斯拟合等等。比如我选择了多项式拟合,在第一个下拉菜单中选择多项式的阶数,我们选择4。我们在窗口的右上角一般勾选对勾。这样拟合曲线会实时显现,如下图所示。
http://s12/mw690/002kh4Qkzy76NILtLM72b&690
在上图中椭圆标识的地方分为两部分。第一部分是关于拟合的方程以及拟合方程中的参数,第二部分是关于拟合的效果。其中SSE:sum of squares due to error(和方差).用来评估拟合值和真实值之间的偏差。
http://s5/mw690/002kh4Qkzy76NIRlldif4&690
这个值是越小越好。
RMSE:标准差(Root Mean Squared Error),越接近0越好。
对于Python来说可以借助Scipy库来进行拟合。安装Scipy必须要安装Numpy+MKL。下载地址如下:
http://www.lfd.uci.edu/~gohlke/pythonlibs/#numpy
在Scipy中的optimize中就包括了curve_fit拟合函数。我们查找它的帮助文档看看它的用法:
https://docs.scipy.org/doc/scipy-0.18.1/reference/generated/scipy.optimize.curve_fit.html
以一个例子来说明它的简单用法。假定我们有了一组数据(x,y):
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
data=np.loadtxt('median.txt')
x=data[1:,0]
y=data[1:,1]
plt.plot(x,y,'*')
plt.show()
我们在进行拟合之前要观察一下看看这个曲线适合用什么方程进行拟合,是多项式还是指数等等。
http://s13/mw690/002kh4Qkzy76NIV5X1a5c&690
我们假设这组数据的拟合公式是:
f=a*e^(bx)+c*e^(dx)
我们对这个公式进行定义:
def func(x,a,b,c,d):
return a*np.exp(b*x)+c*np.exp(d*x)
接着调用curve_fit函数:
popt,pcov=curve_fit(func,x,y)
popt是一个一维数组,表示得到的拟合方程的参数,在本例子中是a,b,c,d。
pcov是一个二维数组,是在popt参数下得到的协方差。
我们将得到的popt和pcov进行输出,得到结果如下:
http://s16/mw690/002kh4Qkzy76NIZNvyLef&690
我们发现协方差大的离谱,我们再看得到的拟合曲线图:
http://s4/mw690/002kh4Qkzy76NJ2gU4b23&690
与原图相去甚远。这是什么原因呢?这可能是因为curve_fit的算法不是很好,在默认情况下曲线拟合使用的a,b,c,d都有一个初始值,默认全部为1。导致拟合曲线不佳的情况可能是初始值与终值相差太大,所以我们需要对初始值进行优化。
关于如何选取初始值有较多的文献讨论,有兴趣的可以去下载查看,这里我们选择一种比较简单的方法,那是啥?那就是观察!嗯,是的,观察,这很玄!我们观察到要拟合的曲线并不是一直是快速增加的,所以我们猜测c和d可能是负的。对程序进行修改:
popt,pcov=curve_fit(func,x,y,[1,1,-1,-1])
得到的a,b,c,d参数和协方差以及拟合效果图,如下图:
http://s9/mw690/002kh4Qkzy76NJ3Adkcc8&690
看,效果已经很好了。将以上代码进行汇总如下:
#使用指定函数进行曲线拟合
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
data=np.loadtxt('median.txt')
x=data[1:,0]
y=data[1:,1]
def
func(x,a,b,c,d):
return a*np.exp(b*x)+c*np.exp(d*x)
popt,pcov=curve_fit(func,x,y,[1,1,-1,-1])
print('popt=',popt)
print('pcov=',pcov)
a=popt[0]
b=popt[1]
c=popt[2]
d=popt[3]
yvals=func(x,a,b,c,d)
plot1=plt.plot(x,y,'*')
plot2=plt.plot(x,yvals,'r')
plt.show()
我们比较一下使用Matlab的拟合工具箱得到的结果:
http://s7/mw690/002kh4Qkzy76NJay0DQd6&690
|
|
a |
b |
c |
d |
|
Python |
10.71923515 |
0.38822168 |
-9.56577783 |
-109.89613443 |
|
Matlab |
10.72 |
0.3882 |
-9.566 |
-109.9 |
|
差值 |
0.000765 |
2.17E-5 |
0.000222 |
0.003866 |
我们可以看到这两者拟合方程参数相差其实很小,Python给出的值小数点后更多一点而已。
其实曲线拟合背后的原理两者差不多,但是Matlab的拟合工具箱实在太好用了,点几下就能实时看结果,非常方便,而且还能自动生成代码。一般情况下我们会使用cftool查看最好的拟合方式,然后用代码写出,比如上边这个例子:
p=fit(x,y,'exp2')
plot(p,x,y)
没了,就两行代码。Python可以替代Matlab吗?能做到,但是比较麻烦。

加载中…