python中带有两个变量的非线性曲线拟合

Question

我正在尝试定义一个适合表单输入 x 和 y 数据的函数：

def nlvh(x,y, xi, yi, H,C):

    return ((H-xi*C)/8.314)*((1/xi) - x) + (C/8.314)*np.log((1/x)/xi) + np.log(yi)

x 和 y 数据是相同长度的一维 numpy 数组。我想对数据进行切片，以便我可以选择 x 和 y 的前 5 个点，通过优化模型中的 C 和 H 来拟合它们，然后向前移动一个点并重复。我有一些代码可以对相同的数据进行线性拟合：

for i in np.arange(len(x)):
    xdata = x[i:i + window]
    ydata = y[i:i+window]
    a[i], b[i] = np.polyfit(xdata, ydata,1)
    xdata_avg[i] = np.mean(xdata)
    if i == (lenx - window):
        break

但是在上面定义的等式上做同样的事情似乎有点棘手。 x 和 y 显示为自变量和因变量，但也有参数 xo 和 yo，它们是每个窗口中 x 和 y 的第一个值。

我想要的最终结果是两个带有 H[i] 和 C[i] 的新数组，其中 i 指定每个后续窗口。有人对我如何开始有一些见解吗？

Answer 1

根据您对我之前的答案的评论（您建议您希望 xi 和 yi 成为每个“切片”x 和 y 数组中的初始值），我正在添加另一个答案.这个答案引入了函数nlvh 的变化，并完全实现了你想要的。就像我之前的回答一样，我们将使用 scipy.optimize 中的curve_fit。

在下面提到的代码中，我使用 python 中的globals() 函数来定义xi 和yi。对于每个切片的x 和y 数组，xi 和yi 存储各自切片数组的第一个值。这是修改后的代码：

from __future__ import division #For decimal division.
import numpy as np
from scipy.optimize import curve_fit

def nlvh(x, H, C):
    return ((H-xi*C)/8.314)*((1/xi) - x) + (C/8.314)*np.log((1/x)/xi) + np.log(yi)

xdata = np.arange(1,21) #Choose an array for x.

#Choose an array for y.
ydata = np.array([-0.1404996,  -0.04353953,  0.35002257,  0.12939468, -0.34259184, -0.2906065,
     -0.37508709, -0.41583238, -0.511851,   -0.39465581, -0.32631751, -0.34403938,
     -0.592997,   -0.34312689, -0.4838437,  -0.19311436, -0.20962735, -0.31134191, 
     -0.09487793, -0.55578775])


H_lst, C_lst = [], []
for i in range( len(xdata)-5 ):
    #Select 5 consecutive points of xdata (from index i to i+4).
    xnew = xdata[i: i+5]
    globals()['xi'] = xnew[0]

    #Select 5 consecutive points of ydata (from index i to i+4).
    ynew = ydata[i: i+5]
    globals()['yi'] = ynew[0]  

    #Fit function nlvh to data using scipy.optimize.curve_fit
    popt, pcov = curve_fit(nlvh, xnew, ynew, maxfev=100000)

    #Optimal values for H from minimization of sum of the squared residuals.
    H_lst += [popt[0]] 

    #Optimal values for C from minimization of sum of the squared residuals.
    C_lst += [popt[1]] 

H_arr, C_arr = np.asarray(H_lst), np.asarray(C_lst) #Convert list to numpy arrays.

H_arr 和 C_arr 的输出现在将如下所示：

print H_arr 
>>>[1.0, 1.0, -23.041138662879327, -34.58915200575536, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0]

print C_arr
>>>[1.0, 1.0, -8.795855063863234, -9.271561975595562, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0]

以下是您针对上面选择的数据（xdata、ydata）获得的图表。

Answer 2

您可以使用 scipy.optimize 中的curve_fit。它将使用非线性最小二乘法将函数nlvh 的参数(H, C, xi, yi) 拟合到x 和y 的给定输入数据。

试试下面的代码。在下面提到的代码中，H_arr 和 C_arr 是 numpy 数组，当函数 nlvh 拟合到 xdata 和 @987654332 的 5 个连续点的窗口时，它们分别包含 H 和 C 的拟合参数@（xdata和ydata是我为x和y选择的数组。你可以在这里选择不同的数组。）

from __future__ import division #For decimal division.
import numpy as np
from scipy.optimize import curve_fit

def nlvh(x, H, C, xi, yi):
    return ((H-xi*C)/8.314)*((1/xi) - x) + (C/8.314)*np.log((1/x)/xi) + np.log(yi)

xdata = np.arange(1,21) #Choose an array for x

#Find an array yy for chosen values of parameters (H, C, xi, yi)
yy = nlvh(xdata, H=1.0, C=1.0, xi=1.0, yi=1.0)  

print yy
>>>[ 0. -0.08337108 -0.13214004 -0.16674217 -0.19358166 -0.21551112 -0.23405222 -0.25011325 -0.26428008 -0.27695274 -0.28841656 -0.2988822 -0.30850967 -0.3174233  -0.3257217  -0.33348433 -0.3407762 -0.34765116 -0.35415432 -0.36032382]

#Add noise to the initally chosen array yy.
y_noise = 0.2 * np.random.normal(size=xdata.size)
ydata = yy + y_noise    

print ydata

>>>[-0.1404996  -0.04353953  0.35002257  0.12939468 -0.34259184 -0.2906065 -0.37508709 -0.41583238 -0.511851   -0.39465581 -0.32631751 -0.34403938 -0.592997   -0.34312689 -0.4838437  -0.19311436 -0.20962735 -0.31134191-0.09487793 -0.55578775]

H_lst, C_lst = [], []
for i in range( len(xdata)-5 ):
    #Select 5 consecutive points of xdata (from index i to i+4).
    xnew = xdata[i: i+5]

    #Select 5 consecutive points of ydata (from index i to i+4).
    ynew = ydata[i: i+5]

    #Fit function nlvh to data using scipy.optimize.curve_fit
    popt, pcov = curve_fit(nlvh, xnew, ynew, maxfev=100000)

    #Optimal values for H from minimization of sum of the squared residuals.
    H_lst += [popt[0]] 

    #Optimal values for C from minimization of sum of the squared residuals.
    C_lst += [popt[1]]   

H_arr, C_arr = np.asarray(H_lst), np.asarray(C_lst) #Convert list to numpy arrays.

以下将是 H_arr 和 C_arr 对 xdata 和 ydata 的选定值的输出。

print H_arr

>>>[ -11.5317468   -18.44101926   20.30837781   31.47360697  -14.45018355 24.17226837   39.96761325   15.28776756 -113.15255865   15.71324201 51.56631241  159.38292301  -28.2429133   -60.97509922  -89.48216973]

print C_arr

>>>[0.70339652 0.34734507 0.2664654  0.2062776  0.30740565 0.19066498 0.1812445  0.30169133 0.11654544 0.21882872 0.11852967 0.09968506 0.2288574  0.128909   0.11658227]