python 上的牛顿-高斯方法进行逼近

Question

我正在尝试在 python 上实施牛顿-高斯方法。这是我的完整代码：

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

def gauss_newton(X, Y, max_iter=1000, eps=1e-6): 
    P0 = [1, 1, 1]
    J = np.zeros([len(X), len(P0)])

    for i in range(max_iter):
        j1 = 1
        j2 = P0[0]
        j3 = P0[2]*X
        J[:,0] = j1
        J[:,1] = j2
        J[:,2] = j3

        r = Y - (P0[0] + P0[1]*X + P0[2]*X**2)
        t1 = np.linalg.inv(np.dot(J.T, J))
        t2 = np.dot(t1, J.T)
        t3 = np.dot(t2, r)

        P1 = P0 - t3
        t4 = abs(P1-P0)
        if max(t4) <= eps:
            break
        P0 = P1

    return P0[0], P0[1], P0[2]

X = np.asarray([1, 2, 3, 4, 5, 6,])
Y = np.asarray([5, 7, 9, 11, 13, 11])
C1, C2, C3 = gauss_newton(X, Y)
pred = C1*X/(C2+X)
plt.figure(1, figsize=(6,4), dpi=120)
plt.scatter(x=X, y=Y, c='green', marker='o', label='Data')
plt.plot(X, pred, '--m', label='Model')
plt.legend()
plt.show()

不幸的是，行中存在“奇异矩阵”错误：

t1 = np.linalg.inv(np.dot(J.T, J))

我完全重复了视频 https://www.youtube.com/watch?v=weuKzGFVQV0&t=109s 中的所有要点，除了使用的模型 - [Y = C1 + C2X + C3X**2]，不是[Y = C1*X/C2+X]。因此导数 j1、j2 和 j3（因为此模型中有 3 个参数）的计算方式不同。参数r也是如此。以下是视频中模型的计算方式：

 j1 = -(X/(P0[1]+X))
 j2 = ((P0[0]*X)/(P0[1]+X)**2)
 r = Y - ((P0[0]*X)/(P0[1]+X))
 # also matrix P0 is initialized as P0 = [1,1]

你能解释一下我上面的代码可能有什么问题吗？

Answer 1

不是python问题而是导数问题。

有 Y = P0[0] + P0[1]*X + P0[2]*X**2，发现 J 是关于常数而不是变量的偏导数。

在这种情况下：

j1 = d(Y)/d(P[0]) = 1
j2 = d(Y)/d(P[1]) = X
j3 = d(Y)/d(P[2]) = X**2

奇异矩阵是由于有两行具有相同的值 (j1=1, j2=1)，因为两个向量都是线性相关的。

python 上的牛顿-高斯方法进行逼近

Newton-Gauss method on python for approximation

python

mathematical-optimization

nonlinear-optimization