在 python 中拟合自定义函数

Fit a custom function in python

我正在尝试使用以下函数拟合我的数据:

我使用的数据如下:

X1: 
0        1.0
1      101.0
2      201.0
3      301.0
4      401.0
5      501.0
6      601.0
7      701.0
8      801.0
9      901.0
10    1001.0
11    1101.0
12    1201.0
13    1301.0
14    1401.0
15    1501.0
16    1601.0
17    1701.0
18    1801.0
19    1901.0

Y1: 
0     0.121159
1     0.195525
2     0.167305
3     0.125499
4     0.094138
5     0.071610
6     0.053848
7     0.039890
8     0.031099
9     0.023976
10    0.018271
11    0.013807
12    0.010596
13    0.008033
14    0.006710
15    0.005222
16    0.004299
17    0.003376
18    0.002976
19    0.002659

我调用该函数的代码如下所示:

def logN(X1, mu, SD1):

return  A/X1 * np.exp(-0.5 * (np.log(X1/mu)**2/np.log(SD1)**2))
        params, pcov = curve_fit(logN, X1,Y1) print (params)


plt.plot(X1, Y1, "o") 
plt.plot(X1, logN(X1 ,params[0], params[1]))  
plt.show()

此函数的结果显示参数等于 1,我收到以下警告:

minpack.py:829: OptimizeWarning: Covariance of the parameters could not be estimated

类别=优化警告)

我想知道我是否正确调用了我的函数语法错误的函数。一些想法?

观察结果

您面临多项挑战:

  • 如您所说,您的问题是非线性回归(根据系数),可以使用非线性算法解决,例如 Levenberg Marquardt(在 scipy.optimize.curve_fit 中实现)
  • 您在优化过程中没有考虑 A 系数,但它在您的函数中明确说明(因此它采用的是您的 post 中未详细说明的全局值)并且此 A 系数与 sigma 相关,因为前者包含后者。
  • 您的某些数据不符合对数正态分布(x=1 处的点似乎很可疑)并且没有估计 y 不确定性。这可能会在执行参数优化时阻止正确收敛,然后算法无法计算协方差矩阵。

改进建议:

  • 可以将您的问题重写为涉及二阶多项式的经典 OLS。那么我们就不必依赖 NLLS 算法了。只需将 log-log transform 应用于您的关系以确认它是可以忍受的并获得参数转换公式。如果可用,总是优先选择 OLS 而不是 NLLS。

  • 删除或惩罚(加权)可疑点,最好使用 objective 标准。

  • 调整你的模型函数(这里不考虑)。

MCVE

根据您提供的数据:

import io
import numpy as np
from scipy import optimize
import pandas as pd
import matplotlib.pyplot as plt

data = io.StringIO("""id;x;y;sy
0;1.0;0.121159;1
1;101.0;0.195525;1
2;201.0;0.167305;1
3;301.0;0.125499;1
4;401.0;0.094138;1
5;501.0;0.071610;1
6;601.0;0.053848;1
7;701.0;0.039890;1
8;801.0;0.031099;1
9;901.0;0.023976;1
10;1001.0;0.018271;1
11;1101.0;0.013807;1
12;1201.0;0.010596;1
13;1301.0;0.008033;1
14;1401.0;0.006710;1
15;1501.0;0.005222;1
16;1601.0;0.004299;1
17;1701.0;0.003376;1
18;1801.0;0.002976;1
19;1901.0;0.002659;1
""")
df = pd.read_csv(data, sep=";", index_col="id")

将您的模型函数重写为:

def func(x, A, mu, sigma):
    return (A/x)*np.exp(-((np.log(x/mu)/np.log(sigma))**2)/2)

修改签名

然后我们可以通过向优化算法提供数据和足够智能的初始条件来天真地拟合函数:

popt, pcov = optimize.curve_fit(func, df.x, df.y, sigma=df.sy,
                                p0=(50, 100, 0.1), method="lm")

但是结果不是很理想(未加权):

并且由于可疑点而容易发生变化(用 w=100 惩罚 x=1):

因此 y 测量值的不确定性有助于调整拟合度。

无论如何,由于问题可以线性化,我们应该依赖于此 属性,也可以在 OLS 中使用权重。

线性化

如果您愿意,可以使用 scipy.optimize.least_squares 执行 OLS。我将使用非常方便的 sklearn 框架:

from sklearn.preprocessing import PolynomialFeatures
from sklearn.linear_model import LinearRegression
from sklearn.pipeline import make_pipeline

让我们删除第一个可疑点:

df = df.loc[1:,:]

然后,我们调整输入并执行对数-对数转换:

X = np.log(df.x.values).reshape(-1, 1)
y = np.log(df.y)

我们为二阶多项式创建 OLS 管道:

poly = PolynomialFeatures(2)
linreg = LinearRegression()
model = make_pipeline(poly, linreg)

最后我们将模型调整为数据:

model.fit(X, y)
model.score(X, y) # 0.9982242621455882

它导致:

这似乎是对二次方程的合理调整。然后就是将系数转换回您想要的数量。