創建改變在一定的坡度公式，可用在curve_fit

我目前正與擬合曲線下降到實際生產數據的功能。我有好運氣與創建一個雙曲線，並使用curve_fit從scipy.optimize。目前的功能我用：創建改變在一定的坡度公式，可用在curve_fit

def hyp_func(x,qi,b,di): 
    return qi*(1.0-b*di*x)**(-1.0/b)

我想現在要做的，是在下降一定幅度，轉變爲一個指數函數。我將如何去做這件事，仍然可以在curve_fit中使用（我認爲下面的作品）？我正在嘗試下面的代碼，是這樣做的嗎？或者，還有更好的方法？

def hyp_func2(x,qi,b,di): 
    dlim = -0.003 
    hy = qi*(1.0-b*di*x)**(-1.0/b) 
    hdy = di/(1.0-b*di*x) 
    ex = x[hdy>dlim] 
    qlim = qi*(dlim/di)**(1/b) 
    xlim = ((qi/qlim)**b-1)/(b*-di) 
    ey = qlim*np.exp(dlim*(ex-xlim)) 
    y = np.concatenate((hy[hdy<dlim],ey)) 
    return y

hy是雙曲型方程
hdy是HY衍生
ex是x的一部分後衍生命中DLIM
ey是指數方程

我仍在制定的公式，我沒有獲得連續的功能。

編輯：數據here，並更新方程

來源

2017-09-21 Ditchbuster

請問您可以添加一個鏈接到數據？ –

[這裏描述]（https://stackoverflow.com/q/29382903/832621）的方法可以給你一些見解 –

謝謝，我會看看，並瞭解更多關於'numpy.piecewise（）' – Ditchbuster

對不起，是壞消息，但如果我明白你正在嘗試做的，我覺得是很難有scipy.optimize.curve_fit，或任何從scipy.optimize其他方法做你想做的事情。

最合適的算法設計與連續變量，通常（和curve_fit肯定）工作通過在參數值非常小的變化，找到正確的方向和步長採取提高的結果開始。

但你正在尋找的是一個離散變量到另一個（「指數」），該算法一個函數形式（大致，「電力法」）之間的斷點通常不會做出足夠大的變化在di參數作出哪個值被用作斷點的差異，並且可以決定di並不影響飛度（模型中應用在其他方面di，所以你可能會得到幸運和di可能對有影響適合

來源

2017-09-22 04:55:11

我想我沒有解釋清楚，我試圖阻止雜草。函數應該是連續的，基本上，它應該基於'dlim'限制下降。一旦第一部分下降小於'dlim'，就變成一個不斷下降的指數方程。第二部分全部基於第一部分，因此對變量的更改應改變整體功能。 – Ditchbuster

假設qi>0的斜率實際上是正值，所以我沒有得到-0.003的選擇，況且我認爲de rivative是錯誤的。您可以準確計算lope達到臨界值的值。現在，根據我的經驗，你有兩個選擇。如果你自己定義一個分段函數，你通常會遇到使用numpy數組進行函數調用的麻煩。我通常使用scipy.optimize.leastsq和一個自定義殘差函數。第二個選項是兩個函數之間的連續轉換。根據定義，您可以使其儘可能銳利，因爲價值和坡度已經合適。這兩種方案如下所示

進口matplotlib。pyplot如PLT 進口numpy的爲NP

def hy(x,b,qi,di): 
    return qi*(1.0-b*di*x)**(-1.0/b) 


def abshy(x,b,qi,di):#same as hy but defined for all x 
    return qi*abs(1.0-b*di*x)**(-1.0/b) 


def dhy(x,b,qi,di):#derivative of hy 
    return qi*di*(1.0-b*di*x)**(-(b+1.0)/b) 

def get_x_from_slope(s,b,qi,di):#self explaining 
    return (1.0-(s/(qi*di))**(-b/(b+1.0)))/(b*di) 


def exh(x,xlim,qlim,dlim):#exponential part (actually no free parameters) 
    return qlim*np.exp(dlim*(x-xlim)) 

def trans(x,b,qi,di, s0):#piecewise function 
    x0=get_x_from_slope(s0,b,qi,di) 
    if x<x0: 
     out= hy(x,b,qi,di) 
    else: 
     H0=hy(x0,b,qi,di) 
     out=exh(x,x0,H0,s0/H0) 
    return out 


def no_if_trans(x,b,qi,di, s0,sharpness=10):#continuous transition between the two functions 
    x0=get_x_from_slope(s0,b,qi,di) 
    H0=hy(x0,b,qi,di) 
    weight=0.5*(1+np.tanh(sharpness*(x-x0))) 
    return weight*exh(x,x0,H0,s0/H0)+(1.0-weight)*abshy(x,b,qi,di) 


xList=np.linspace(0,5.5,90) 
hyList=np.fromiter((hy(x,2.2,1.2,.1) for x in xList) ,np.float) 
t1List=np.fromiter((trans(x,2.2,1.2,.1,3.59) for x in xList) ,np.float) 
nt1List=np.fromiter((no_if_trans(x,2.2,1.2,.1,3.59) for x in xList) ,np.float) 


fig1=plt.figure(1) 
ax=fig1.add_subplot(1,1,1) 
ax.plot(xList,hyList) 
ax.plot(xList,t1List,linestyle='--') 
ax.plot(xList,nt1List,linestyle=':') 

ax.set_ylim([1,10]) 
ax.set_yscale('log') 
plt.show()

有兩個解決方案几乎沒有差別，但您使用scipy擬合功能的選項略有不同。第二種解決方案應該很容易與curve_fit

來源

2017-09-26 11:17:34 mikuszefski

創建改變在一定的坡度公式，可用在curve_fit

回答

相關問題