使用scipy.optimize拟合多个数据集，保持相同参数

Question

我正在尝试写一段代码，能够同时处理多个（超过100个）数据集，这些数据集有相同的参数。

举个简单的例子：

假设我有4个不同的数据集，它们大致遵循这个函数：y(x)=A*Sin(bx+c)+dx。它们的参数a、b和d都是一样的。唯一不同的参数是c（相位），不过对于每个数据集，我都知道c应该是什么。

现在我想同时拟合这4个函数，并得到一组拟合参数a、b和d。

我写了一段代码，现在是分别拟合每个函数，这样我得到了4组不同的参数，这并不是我想要的结果。

谢谢你们的帮助。

from scipy import optimize
import matplotlib
from matplotlib import pyplot as plt
from scipy import *

Mat_ydata=[]
Mat_angle=[]
Mat_xdata=[]

定义一个用于拟合的函数

def sin_func(x,a,b,c,d):
   return a*sin(b*x+c)+d*x

生成一些数据：4个正弦函数，它们的幅度、频率和线性偏移是相同的，但每个函数的相位不同，而这些相位我已经知道了。此外，每个函数的定义范围也不同。函数的x和y坐标存储在列表的列表Mat_ydata和Mat_xdata中。相位存储在列表Mat_phase中。

for phase in range(0,181,60):
    angle=float(phase)/180*pi
    num_points=phase+50
    xdata=linspace(num_points/10,num_points/5, num_points)
    ydata=2.11*sin(pi/2*xdata+angle)+3.16*xdata+2*((0.5-rand(num_points))*exp(2*rand(num_points)**2))
    Mat_angle.append(angle)
    Mat_ydata.append(ydata)
    Mat_xdata.append(xdata)

这里是拟合的部分。我现在是分别拟合每个函数，因此每个函数的相位、幅度和线性偏移的值都不同。我希望能够把它们一起拟合。最后，我还会绘制数据和拟合结果。

for i in range(0,len(Mat_angle)):
    fitfunc = lambda p, x: sin_func(x,p[0],p[1],Mat_angle[i],p[2])
    errfunc = lambda p, x, y: fitfunc(p, x) - y 
    p0 = [2,1.5,3] 
    p1, success = optimize.leastsq(errfunc, p0[:], args=(Mat_xdata[i], Mat_ydata[i]))
    print p1
    plt.plot(Mat_xdata[i],Mat_ydata[i],'o',Mat_xdata[i],fitfunc(p1,Mat_xdata[i]),'-')
plt.show()

Answer 1

我这几天其实已经搞明白了。我会提供一些代码，可能对有需要的人有帮助。我还发现，拟合正弦函数挺难的，所以我把我的假数据改成了洛伦兹函数。

首先是导入模块和生成假数据，这些数据会保存在一个列表的列表里：

from scipy import optimize
import matplotlib
from matplotlib import pyplot as plt
from scipy import *
import numpy as np


Mat_ydata=[]
Mat_angle=[]
Mat_xdata=[]

for c in range(0,100,20):
    num_points=c+100
    xdata=linspace(0,num_points/2, num_points)
    ydata=5.1/((xdata-c)**2+2.1**2)+0.05*((0.5rand(num_points))*exp(2*rand(num_points)**2))
    Mat_angle.append(c)
    Mat_ydata.append(ydata)
    Mat_xdata.append(xdata)

接下来是定义拟合函数和误差函数：

def lor_func(x,c,par):
    a,b,d=par
    return a/((x-c)**2+b**2)

def err (p,c,x,y):
    return lor_func(x,c,p)-y

def err_global(p,Mat_a,Mat_x,Mat_y):
    err0=[]
    for i in range(0, len(Mat_a)):
        errc=err(p,Mat_a[i],Mat_x[i],Mat_y[i])
        err1=np.concatenate((err0,errc))
        err0=err1
    return err0

然后进行实际的拟合，并显示结果：

p_global=[1,1,1]
p_best,success=optimize.leastsq(err_global, p_global,args=(Mat_angle,Mat_xdata,Mat_ydata),maxfev=40000)
toplot=[]
for i in range(0,len(Mat_angle)):
    toplot.append(lor_func(Mat_xdata[i],Mat_angle[i],p_best))
err_toplot=err_global(p_best,Mat_angle,Mat_xdata,Mat_ydata)
print p_best

for i in range(0,len(Mat_angle)):
    plt.plot(Mat_xdata[i],Mat_ydata[i],'o',Mat_xdata[i],toplot[i],'-')
plt.show()

最后，这里是显示的结果：