0

我有一个要最小化的超定系统，与线性最小二乘非常相似（但不完全相同）。它是这样的形式：

AX = Y

其中 A 需要为 2x3，因为 X 是 3xN 而 Y 是 2xN。

这就是问题所在：这里的未知数是变换矩阵 A，而不是通常最小二乘公式中的 X。（另外，我会提到在这种情况下 X 和 Y 是矩阵，它们通常是经典最小二乘法中的向量）

更复杂的是，我希望这是一个加权最小二乘法，其中对于 Y 中的每个 2x1 向量（Y 是 2xN，这些向量的堆栈），我有一个 2x1 权重向量，它按元素对 Y 中的向量进行加权（像 np.divide(Y[i], np.square(weights[i])))。

我首先尝试在没有权重的情况下解决问题，而我做到了这一点：

AX = Y

XT * AT = YT

X * XT * AT = X * YT

AT = (X * XT)^-1 * X * YT

A = Y * XT * (X * XT)^-1

这导致了这个 Python 函数：

def mat_solve(X,Y):
return Y@X.T@np.linalg.inv(X@X.T)

但是我不能明显地看到这真的最小化了问题的成本函数，所以我不确定它是否真的正确，即使没有权重。

我还尝试将成本函数写为总和，并取一个导数并将其设置为 0，这导致了这个 Python 函数：

def sum_solve(X,Y):
result1 = 0
result2 = 0
for i in range(len(X)):
    result1 += Y[:,i:i+1]@X.T[i:i+1,:]
    result2 += X[:,i:i+1]@X.T[i:i+1,:]
return result1@np.linalg.inv(result2)

这两个函数返回一个非常相似的矩阵，但并不完全相同。如果它们都是真正解决问题的方法，那应该是相同的。我要提到的是，被最小化的空间不是非线性的，而且我非常有信心只有一个最小值。

我想我会走一条更强力的路线，并尝试使用一个名为 JAX 的 autodiff 包，它有自己的 minimize() 函数，与 scipy 非常相似（但以梯度下降的方式），只需将成本函数写出来明确地说，即使有权重。权重只是变成了一堆单位矩阵，然后乘以它们的 Y[i]s，如您所见（我将成本函数称为 myLSQ，然后在其上调用了最小化（），因此它相对于梯度下降A):

import jax 
import jax.numpy as jnp
from jax.scipy.optimize import minimize

uncertainty_new = [np.identity(2)/(x**2) for x in uncertainty_list]
freq_list_new = [np.append(x,[1]) for x in freq_list]

def myLSQ(A,X,Y,W):
    Anew = A.reshape((2,3))
    temp = []
    result = 0
    for i in range(len(X)):
        temp = jnp.subtract(jnp.dot(Anew,X[i]), Y[i])
        result += jnp.dot(jnp.dot(temp,W[i]), temp.T)
    return result



A_guess = np.array([[1.,1.],[1.,1.],[1.,1.]])
%prun result = minimize(myLSQ, A_guess.flatten(),method = "BFGS",args = (freq_list_new, beam_peak_list, uncertainty_new,))

但是这个东西运行了 20 多分钟，然后返回一个绝对可怕的矩阵，它没有接近于最小化问题。

我为这篇庞大的文章道歉，但这是我的问题。

  1. 为什么我的两种（未加权）方法给出的结果略有不同，我如何将权重添加到它们中？

  2. 有一个更好的方法吗？如果是这样，它是什么？

  3. 我希望对拟合进行加权，但我不想在以后使用它时给 A 矩阵权重（一旦它正确）。我希望能够只做 AX = Y，在 Y 上没有权重。权重只是为了修改拟合。我该怎么做呢？

4

0 回答 0