我正在嘗試計(jì)算兩種形式的指數(shù)對(duì)某些x,y數(shù)據(jù)的最佳擬合(數(shù)據(jù)文件可以從here下載)

這是代碼：

from scipy.optimize import curve_fit
import numpy as np

# Get x,y data
data = np.loadtxt('data.txt', unpack=True)
xdata, ydata = data[0], data[1]

# Define first exponential function
def func(x, a, b, c):
    return a * np.exp(b * x)   c

# Get parameters estimate
popt, pcov = curve_fit(func, xdata, ydata)

print popt

# Define second exponential function (one more parameter)
def func2(x, a, b, c, d):
    return a * np.exp(b * x   c)   d

# Get parameters estimate
popt2, pcov2 = curve_fit(func2, xdata, ydata)

print popt2

對(duì)于popt,第一個(gè)指數(shù)給出與(PDF here)完全相同的值：

[  7.67760545e-15   1.52175476e 00   2.15705939e-02]

但是第二個(gè)給出了popt2明顯錯(cuò)誤的值：

[ -1.26136676e 02  -8.13233297e-01  -6.66772692e 01   3.63133641e-02]

對(duì)于相同的第二個(gè)函數(shù),這是值(PDF here)：

a = 6.2426224704624871E-15
b = 1.5217697532005228E 00
c = 2.0660424037614489E-01
d = 2.1570805929514186E-02

我嘗試將列表數(shù)組作為推薦在這里Strange result with python’s (scipy) curve fitting,但這沒有幫助.我在這做錯(cuò)了什么？

加1

我猜這個(gè)問題與缺少初始值有關(guān)我正在提供我的功能(如下所述：gaussian fit with scipy.optimize.curve_fit in python with wrong results)

如果我將估計(jì)從第一個(gè)指數(shù)提供給第二個(gè)指數(shù)(如此)(使新參數(shù)d最初為零)：

popt2, pcov2 = curve_fit(func2, xdata, ydata, p0 = [popt[0], popt[1], popt[2], 0]) 

與相比,我得到的結(jié)果更合理但仍然錯(cuò)誤：

[  1.22560853e-14   1.52176160e 00  -4.67859961e-01   2.15706930e-02]

所以現(xiàn)在問題改為：如何自動(dòng)為我的第二個(gè)函數(shù)提供更合理的參數(shù)？

解決方法:

請(qǐng)注意,zunzun和第一個(gè)模型的估算中a = 0.所以他們只是估計(jì)一個(gè)常數(shù).因此,在第一種情況下b和在第二種情況下b和c是無關(guān)緊要的并且未被識(shí)別.

Zunzun也使用差分進(jìn)化作為全局求解器,這是我最后一次看到它. Scipy現(xiàn)在將流水作為全局優(yōu)化器看起來相當(dāng)不錯(cuò),在可能的局部最小值的情況下值得一試.

我的“便宜”方式,因?yàn)閰?shù)在您的示例中沒有大范圍：嘗試隨機(jī)起始值

np.random.seed(1)
err_last = 20
best = None

for i in range(10):
    start = np.random.uniform(-10, 10, size=4)
    # Get parameters estimate
    try:
        popt2, pcov2 = curve_fit(func2, xdata, ydata, p0=start)
    except RuntimeError:
        continue
    err = ((ydata - func2(xdata, *popt2))**2).sum()
    if err < err_last:
        err_last = err
        print err
        best = popt2


za = 6.2426224704624871E-15
zb = 1.5217697532005228E 00
zc = 2.0660424037614489E-01
zd = 2.1570805929514186E-02

zz = np.array([za,zb,zc,zd])
print 'zz', zz
print 'cf', best

print 'zz', ((ydata - func2(xdata, *zz))**2).sum()
print 'cf', err_last

最后一部分打印(zz是zunzun,cf是curve_fit)

zz [  6.24262247e-15   1.52176975e 00   2.06604240e-01   2.15708059e-02]
cf [  1.24791299e-16   1.52176944e 00   4.11911831e 00   2.15708019e-02]
zz 9.52135153898
cf 9.52135153904

與Zunzun不同的參數(shù)為b和c,但殘差平方和相同.

加成

a * np.exp(b * x c)d = np.exp(b * x(c np.log(a)))d

要么

a * np.exp(b * x c)d =(a * np.exp(c))* np.exp(b * x)d

第二個(gè)功能與第一個(gè)功能沒有什么不同. a和c未單獨(dú)標(biāo)識(shí).因此,使用衍生信息的優(yōu)化器也會(huì)遇到問題,因?yàn)槿绻艺_地看到這個(gè),Jacobian在某些方向上是單數(shù)的.