在 Python 中使用 Levenberg-Marquardt 算法优化方程组

Question

我希望使用 scipy 中的 scipy.optimize.leastsq() 方法来优化三个参数 a,b,c。我所拥有的是这两个方程式。

1*a+2*b+3*c = x1
4*a+5*b+6*c = x2

分析这组方程是不确定的，但在数值上我试图找到 a,b,c 以最小化测量给定结果的误差 [2,2]:

1*a+2*b+3*c - 2 = 0
4*a+5*b+6*c - 2 = 0

所以我写了一些代码。

def function(a,b,c,t):
    return np.array([1*a+2*b+3*c+t[1],4*a+5*b+6*c+t[1]])

a0 = 1
b0 = 1
c0 = 1
measdata = np.array([2,2])
t = [1,2]

def residual(x0,measdata,t):
    return measdata - function(x0[0],x0[1],x0[2],t)

erg = optimize.leastsq(func=residual,x0=(a0,b0,c0),args=(measdata,t))

它总是导致：

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-296-ab0fc90a2253> in <module>()
     14     return result - function(x0[0],x0[1],x0[2],t)
     15 
---> 16 erg = optimize.leastsq(func = residual, x0 = (a0,b0,c0) , args=(result,t), maxfev=10000)
     17 
     18 function(erg[0][0],erg[0][1])

    //anaconda/lib/python3.5/site-packages/scipy/optimize/minpack.py in leastsq(func, x0, args, Dfun, full_output, col_deriv, ftol, xtol, gtol, maxfev, epsfcn, factor, diag)
        378     m = shape[0]
        379     if n > m:
    --> 380         raise TypeError('Improper input: N=%s must not exceed M=%s' % (n, m))
        381     if epsfcn is None:
        382         epsfcn = finfo(dtype).eps
    TypeError: Improper input: N=3 must not exceed M=2

如何找到最小值？我知道这只是当地的最低要求，但我会很高兴。

Answer 1

错误告诉你你已经知道的，即系统未定，n 是参数的数量，m 是约束的数量。

如果您修改其中一个参数，使 n > m 为 False，代码将停止报错。例如，更改

def residual(x0,measdata,t):
    return measdata - function(x0[0],x0[1],x0[2],t)

erg = optimize.leastsq(func=residual,x0=(a0,b0,c0),args=(measdata,t))

至

def residual(x0,measdata,t):
    # we fix the value of `c` here
    return measdata - function(x0[0],x0[1],5,t)

# only two parameters for `x0`
erg = optimize.leastsq(func=residual,x0=(a0,b0),args=(measdata,t))

要回答如何做自己想做的问题，我不确定 scipy 是否可以做到。我发现这个 issue 说 scipy 无法处理不确定的系统：

Interesting, I assumed the MINPACK routines would handle also m < n, but apparently not. The reason why they don't is probably that for m < n the minimum is some manifold of points, which causes problems in the termination conditions.

So would be after all some interest in adding also a small-scale solver for underdetermined least-square problems.

即使那个 post 是 3 年前的，我仍然无法在文档中找到任何证据表明 scipy 可以做你想做的事。但是，我找到了一个声称的 SO 答案，但我还没有完全掌握数学知识以确保它是否适用于您的情况。由于我发现很难总结 post，我将只引用似乎最重要的部分。

In the case where A·x = b is underdetermined,
x1, res, rnk, s = np.linalg.lstsq(A, b)
will pick a solution x' that minimizes ||x||_L2 subject to ||A·x - b||_L2 = 0. This happens not to be the particular solution we are looking for, but we can linearly transform it to get what we want. In order to do that, we'll first compute the right null space of A, which characterizes the space of all possible solutions to A·x = b. We can get this using a rank-revealing QR decomposition

在 Python 中使用 Levenberg-Marquardt 算法优化方程组

Optimizing set of equations with Levenberg-Marquardt algorithm in Python

python

numpy

scipy

levenberg-marquardt