用CURE_FIT估计不同大小数据集上的常用模型参数

2022-05-23 00:00:00 python numpy arrays curve-fitting

问题描述

我正在处理一个曲线拟合问题,其中我打算在不同大小的几个数据集上全局估计共享的模型参数。我从下面链接中的代码开始工作,其中线性回归y=a*x+b的公共a参数是在具有公共x向量的三个不同y向量上估计的。How to use curve_fit from scipy.optimize with a shared fit parameter across multiple datasets?

我设法使代码样本适应更一般的情况,使用三个不同的x向量,每个y数据向量对应一个x向量。但是,当我想进一步扩展它以处理大小不等的数据集时,我遇到了以下错误:&Quot;ValueError:使用序列设置数组元素。&Quot;.

请找到下面的代码示例。如有任何帮助,我们不胜感激!

干杯

import matplotlib.pyplot as plt
import numpy as np
from scipy.optimize import curve_fit

x = [[0, 1, 2, 3],
[0.2, 1.2, 2.2, 3.2],
[0.3, 1.3, 2.3]]

y = [[-0.80216234,  1.41125365,  1.42565202,  2.42567754],
[ 1.34166743,  1.29731851,  2.98374731,  3.32110875],
[ 1.71398203,  3.29737756,  3.81456949]]

x = np.array(x)
y = np.array(y)


def f(x, a, b):
    return a * x + b


def g(x, a, b_1, b_2, b_3):
     return np.concatenate((f(x[0], a, b_1), f(x[1], a, b_2), f(x[2], a, b_3)))

(a, *b), _ = curve_fit(g, x, y.ravel())

for x_i, y_i, b_i in zip(x, y, b):
    plt.plot(x_i, f(x_i, a, b_i), label=f"{a:.1f}x{b_i:+.1f}")
    plt.plot(x_i, y_i, linestyle="", marker="x",  color=plt.gca().lines[-1].get_color())
plt.legend()
plt.show()

有关具有多个大小相等的x向量的工作示例的代码,请参见下面的代码:

import matplotlib.pyplot as plt
import numpy as np
from scipy.optimize import curve_fit

x = [[0, 1, 2, 3],
 [0.2, 1.2, 2.2, 3.2],
 [0.3, 1.3, 2.3, 3.3]]

y = [[-0.80216234,  1.41125365,  1.42565202,  2.42567754],
 [ 1.34166743,  1.29731851,  2.98374731,  3.32110875],
 [ 1.71398203,  3.29737756,  3.81456949, 4.25]]

x = np.array(x)
y = np.array(y)


def f(x, a, b):
    return a * x + b


def g(x, a, b_1, b_2, b_3):
    return np.concatenate((f(x[0], a, b_1), f(x[1], a, b_2),     f(x[2], a, b_3)))

(a, *b), _ = curve_fit(g, x, y.ravel())

for x_i, y_i, b_i in zip(x, y, b):
   plt.plot(x_i, f(x_i, a, b_i), label=f"{a:.1f}x{b_i:+.1f}")
   plt.plot(x_i, y_i, linestyle="", marker="x",     color=plt.gca().lines[-1].get_color())

plt.legend()
plt.show()

解决方案

总的来说,我同意最小二乘拟合的含义相当复杂...以下是我的快速想法:

  • 您能确保从不同长度的数据集中估计得到的b参数同样有效吗?
  • 获得的b参数越多,其估计的不确定性就越大,因为您只优化组合的匹配性能,而不是每个单独的匹配性能
  • 我也不确定雅可比的数值计算在这种情况下会有多好…可能值得实现一个自定义jac函数,该函数以一种精确的方式计算雅可比
  • .我肯定还有更多我目前不知道的问题:D
不过,您当然可以欺骗scipy.optimize做您想做的事情...
但是,您必须深入一步,直接使用scipy.optimize.least_squares,而不是使用更高级别的scipy.optimize.curve_fit函数。

通过这种方式,您可以更改残差的计算方式,以接受不同长度的数据集。

.以下是它如何工作的快速实现:

import matplotlib.pyplot as plt
import numpy as np
from scipy.optimize import least_squares

x = [[0.0, 1.0, 2.0, 3.0, 4.0, 5.0],
     [0.2, 1.2, 2.2, 3.2],
     [0.3, 1.3, 2.3     ]]

y = [[-0.80216234, 1.41125365, 1.42565202, 2.42567754, 3, 4],
     [ 1.34166743, 1.29731851, 2.98374731, 3.32110875],
     [ 1.71398203, 3.29737756, 3.81456949            ]]


def f(x, a, b):
    return a * x + b


def fun(parameters):
    # separate a and b parameters
    a, *b = parameters

    # calculate function-results based on a shared a- and variable b- parameters
    res = (f(xi, a, bi) for (xi, bi) in zip(map(np.array, x), b))

    # calculate the residuals
    errs = []
    for i, j in zip(res, map(np.array, y)):
        errs += (i - j).tolist()
    return np.array(errs)


# set start-values
start_values = (1, 1, 2, 3)
# do the fit
a, *b = least_squares(fun, start_values).x

for x_i, y_i, b_i in zip(map(np.array, x), y, b):
    plt.plot(x_i, f(x_i, a, b_i), label=f"{a:.1f}x{b_i:+.1f}")
    plt.plot(x_i, y_i, linestyle="", marker="x", color=plt.gca().lines[-1].get_color())

plt.legend()
plt.show()

相关文章