matrix - MAXIMA - 使用最小二乘估计解决超定系统

Question

我正在拼命尝试使用 Maxima 和 lsquares_estimate 函数来解决超定系统。

系统由以下矩阵方程描述: M.D = I

这是M:

M:matrix([19.707,0.06700000000000017,0,0,0,0,0,0,19.707-Phi01_8,0,0,0],
[0,0,-0.2690000000000019,0.7359999999999971,0,0,0,0,0,19.371-Phi01_8,0,0],
[-0.067234,-24.562234,0,0,0,0,0,0,-Phi12_8-0.067234,0,0,0],
[0,0,-1.753,3.741,0,0,0,0,0,22.742-Phi12_8,0,0],
[-0.26894,1.75226,0,0,0,0,0,0,-Phi23_8-0.26894,0,0,0],   
[0,0,92.99119999999999,14.967,0,0,0,0,0,90.97-Phi23_8,0,0],
[-0.73571,3.74149,0,0,0,0,0,0,-Phi34_8-0.73571,0,0,0],
[0,0,-14.9658,-655.2429999999999,0,0,0,0,0,-Phi34_8-19.443,0,0],
[-1.8393,9.3537,0,0,0,0,0,0,-Phi45_8-1.8393,0,0,0],
[0,0,-37.415,361.892,0,0,0,0,0,-Phi45_8-48.608,0,0],
[-8.6658,5.6092,0,0,0,0,0,0,-Phi56_8-8.6658,0,0,0],
[0,0,-22.437,166.598,0,0,0,0,0,-Phi56_8-36.712,0,0],
[-5.1995,3.3655,0,0,0,0,0,0,-Phi67_8-5.1995,0,0,0],
[0,0,-13.462,99.953,0,0,0,0,0,-Phi67_8-22.027,0,0])

这里是D:

D:matrix([0.0499998750003125],
[0.0399993670651863],
[0.009999946119925522],
[0.001250003764581481],
[4.999947987060861*10^-4],
[1.000099826698366*10^-4],
[1.666560996027836*10^-4],
[0.00499998750003125],
[D8],[D9],[D10],[D11])

这是我:

I:matrix([1],[0],[-1],[0],[0],[1],[0],[-1],[0],[0],[0],[0],[0],[0])

该系统有 14 个方程和 9 个变量。 (D8、D9、Phi01、Phi12、Phi23、Phi34、Phi45、Phi56、Phi67) 我想获得这些变量的最小二乘估计。我正在使用“lsquares_estimates”函数。为了给函数提供正确的输入，我将矩阵方程合并为一个矩阵。这就是我调用 lsquares_estimates 函数的方式:

D0:D[1][1]$
D1:D[2][1]$
D2:D[3][1]$
D3:D[4][1]$
D4:D[5][1]$
D5:D[6][1]$
D6:D[7][1]$
D7:D[8][1]$
D8:D[9][1]$
D9:D[10][1]$
D10:D[11][1]$
D11:D[12][1]$

C:addcol(M,I);
resultat:lsquares_estimates(C,
[c0,c1,c2,c3,c4,c5,c6,c7,c8,c9,c10,c11,Ie],
Ie=c0*D0+c1*D1+c2*D2+c3*D3+c4*D4+c5*D5+c6*D6+c7*D7+c8*D8+c9*D9+c10*D10+c11*D11,
[Phi01_8,Phi12_8,Phi23_8,Phi34_8,Phi45_8,Phi56_8,Phi67_8,D8,D9]);

显然，Maxima 无法解析该系统...它持续运行很长时间，有时，这是我收到的错误:

Maxima encountered a Lisp error:
-: ((MTIMES SIMP) -0.14285714285714285
 ((MPLUS SIMP)
  ((MTIMES SIMP)
   ((MPLUS SIMP) 0.01197250577547404
    ((MTIMES SIMP) -1.0 ((MPLUS SIMP) 19.707 ((MTIMES SIMP) -1.0 |$Phi01_8|))))
   ((MPLUS SIMP) 19.707 ((MTIMES SIMP) -1.0 |$Phi01_8|)))
  ((MTIMES SIMP)
   ((MPLUS SIMP) -0.014164494697229801
    ((MTIMES SIMP) -1.0
     ((MPLUS SIMP) -0.067234 ((MTIMES SIMP) -1.0 |$Phi12_8|))))
   ((MPLUS SIMP) -0.067234 ((MTIMES SIMP) -1.0 |$Phi12_8|)))
  ((MTIMES SIMP)
   ((MPLUS SIMP) -0.056642324551059296
    ((MTIMES SIMP) -1.0
     ((MPLUS SIMP) -0.26894 ((MTIMES SIMP) -1.0 |$Phi23_8|))))
   ((MPLUS SIMP) -0.26894 ((MTIMES SIMP) -1.0 |$Phi23_8|)))
  ((MTIMES SIMP)
   ((MPLUS SIMP) -0.11287182384424398
    ((MTIMES SIMP) -1.0
     ((MPLUS SIMP) -0.73571 ((MTIMES SIMP) -1.0 |$Phi34_8|))))
   ((MPLUS SIMP) -0.73571 ((MTIMES SIMP) -1.0 |$Phi34_8|)))
  ((MTIMES SIMP)
   ((MPLUS SIMP) -0.28217730962955834
    ((MTIMES SIMP) -1.0
     ((MPLUS SIMP) -1.8393 ((MTIMES SIMP) -1.0 |$Phi45_8|))))
   ((MPLUS SIMP) -1.8393 ((MTIMES SIMP) -1.0 |$Phi45_8|)))
  ((MTIMES SIMP)
   ((MPLUS SIMP) 0.2089244670356651
    ((MTIMES SIMP) -1.0
     ((MPLUS SIMP) -8.6658 ((MTIMES SIMP) -1.0 |$Phi56_8|))))
   ((MPLUS SIMP) -8.6658 ((MTIMES SIMP) -1.0 |$Phi56_8|)))
  ((MTIMES SIMP)
   ((MPLUS SIMP) 0.12535648020624035
    ((MTIMES SIMP) -1.0
     ((MPLUS SIMP) -5.1995 ((MTIMES SIMP) -1.0 |$Phi67_8|))))
   ((MPLUS SIMP) -5.1995 ((MTIMES SIMP) -1.0 |$Phi67_8|))))) is not a number
Automatically continuing.
To enable the Lisp debugger set *debugger-hook* to nil.

我使用 lsquares_estimates 的方式是否错误？当我减少方程数量以匹配变量数量时，我可以获得结果，但我希望使用所有方程得到准确的结果。

任何帮助将非常感激...... 提前致谢。

Answer 1

lsquares_estimates 构造一个均方误差表达式，然后尝试精确求解它，然后近似求解(如果精确求解不能产生结果)。不确定哪里出了问题。

但无论如何，您可以直接跳转到lsquares_estimates_approximate。试试这个:

e1 : M.D - I;
my_mse : transpose(e1) . e1;
load (lsquares);
lsquares_estimates_approximate (my_mse, [Phi01_8, Phi12_8, Phi23_8, 
    Phi34_8, Phi45_8, Phi56_8, Phi67_8, D8, D9], tol=1e-4);

请注意，my_mse 是误差的平方范数。

通过这些输入，我得到:

[[Phi01_8 = 22.70059866292853,Phi12_8 = -4.329661712144546,
  Phi23_8 = -17.40507905561839,Phi34_8 = -34.32504757490423,
  Phi45_8 = -85.85849926978096,Phi56_8 = 54.09091751735507,
  Phi67_8 = 32.46719922407783,D8 = -0.003346473781793617,
  D9 = 2.1359057759474724E-5]]