为什么在Matlab和Octave中inv()和pinv()的输出不相等?

mym*_*314 14 precision matlab floating-accuracy matrix-inverse numerical-analysis

我注意到如果A是NxN矩阵并且它具有逆矩阵.但是inv()和pinv()函数输出是不同的. - 我的环境是Win7x64 SP1,Matlab R2012a,Cygwin Octave 3.6.4,FreeMat 4.2

看看Octave的例子:

A = rand(3,3)
A =
0.185987   0.192125   0.046346
0.140710   0.351007   0.236889
0.155899   0.107302   0.300623

pinv(A) == inv(A)
ans =
0 0 0
0 0 0
0 0 0
Run Code Online (Sandbox Code Playgroud)
  • 这一切都同ans通过上面运行在Matlab同一命令的结果.

  • 我计算inv(A)*AA*inv(A),结果是Octave和Matlab中的身份3x3矩阵.
  • 结果A*pinv(A)pinv(A)*A在Matlab和FreeMat身份3x3矩阵.
  • 的结果 A*pinv(A)是Octave中的身份3x3矩阵.
  • 结果pinv(A)*A不是在倍频身份3x3矩阵.

我不知道原因 inv(A) != pinv(A),我已经考虑了矩阵中元素的细节.这似乎是导致这个问题的浮动精度问题.

点点后10位数可能不同,如下所示:

  • 6.65858991579923298331777914427220821380615200000000inv(A)(1,1)反对的元素

  • 6.65858991579923209513935944414697587490081800000000 元素 pinv(A)(1,1)

got*_*oth 15

这个问题很老了,但无论如何我都会回答它,因为它似乎在一些谷歌搜索中几乎排在最前面.

我将在我的例子中使用魔法(N)函数,该函数返回N-by-N魔方.

我将创建一个3x3魔方M3,取伪逆PI_M3并乘以它们:

   prompt_$ M3 = magic(3) , PI_M3 = pinv(M3) , M3 * PI_M3
  M3 =

    8   1   6
    3   5   7
    4   9   2

  PI_M3 =

     0.147222  -0.144444   0.063889
    -0.061111   0.022222   0.105556
    -0.019444   0.188889  -0.102778

  ans =

     1.0000e+00  -1.2212e-14   6.3283e-15
     5.5511e-17   1.0000e+00  -2.2204e-16
    -5.9952e-15   1.2268e-14   1.0000e+00

正如您所看到的,答案是单位矩阵,以节省一些舍入误差.我将用4x4魔方重复操作:

   prompt_$ M4 = magic(4) , PI_M4 = pinv(M4) , M4 * PI_M4
 
  M4 =

     16    2    3   13
      5   11   10    8
      9    7    6   12
      4   14   15    1

  PI_M4 =

     0.1011029  -0.0738971  -0.0613971   0.0636029
    -0.0363971   0.0386029   0.0261029   0.0011029
     0.0136029  -0.0113971  -0.0238971   0.0511029
    -0.0488971   0.0761029   0.0886029  -0.0863971

  ans =

     0.950000  -0.150000   0.150000   0.050000
    -0.150000   0.550000   0.450000   0.150000
     0.150000   0.450000   0.550000  -0.150000
     0.050000   0.150000  -0.150000   0.950000

结果不是单位矩阵,这意味着4x4魔方不具有逆.我可以通过尝试Moore-Penrose伪逆的一个规则来验证这一点:

   prompt_$ M4 * PI_M4 * M4
  
ans =

   16.00000    2.00000    3.00000   13.00000
    5.00000   11.00000   10.00000    8.00000
    9.00000    7.00000    6.00000   12.00000
    4.00000   14.00000   15.00000    1.00000

满足规则A*B*A = A. 这表明pinv在可用时返回逆矩阵,而当逆矩阵不可用时返回伪逆.这就是为什么在某些情况下你会得到一个小的差异,只是一些舍入错误,在其他情况下你会得到更大的差异.为了显示它,我将获得两个魔术象限的倒数并从伪逆中减去它们:

   prompt_$ I_M3 = inv(M3) , I_M4 = inv(M4) , DIFF_M3 = PI_M3 - I_M3, DIFF_M4 = PI_M4 - I_M4
  I_M3 =

     0.147222  -0.144444   0.063889
    -0.061111   0.022222   0.105556
    -0.019444   0.188889  -0.102778

  warning: inverse: matrix singular to machine precision, rcond = 1.30614e-17
  I_M4 =

     9.3825e+13   2.8147e+14  -2.8147e+14  -9.3825e+13
     2.8147e+14   8.4442e+14  -8.4442e+14  -2.8147e+14
    -2.8147e+14  -8.4442e+14   8.4442e+14   2.8147e+14
    -9.3825e+13  -2.8147e+14   2.8147e+14   9.3825e+13

  DIFF_M3 =

     4.7184e-16  -1.0270e-15   5.5511e-16
    -9.9226e-16   2.0470e-15  -1.0825e-15
     5.2042e-16  -1.0270e-15   4.9960e-16

  DIFF_M4 =

    -9.3825e+13  -2.8147e+14   2.8147e+14   9.3825e+13
    -2.8147e+14  -8.4442e+14   8.4442e+14   2.8147e+14
     2.8147e+14   8.4442e+14  -8.4442e+14  -2.8147e+14
     9.3825e+13   2.8147e+14  -2.8147e+14  -9.3825e+13


Ste*_*fin 8

在我看来,你在底部回答了自己的问题.原因是浮点运算.算法inv()pinv()不完全相同,因为pinv()必须能够处理非方形矩阵.因此答案将不完全相同.

如果你看一下它的值pinv(A)*A,你会发现它非常接近单位矩阵.

我明白了:

ans =

   1.0000e+00   6.1062e-16  -3.0809e-15
  -5.8877e-15   1.0000e+00   6.3942e-15
   2.4425e-15  -3.0184e-16   1.0000e+00
Run Code Online (Sandbox Code Playgroud)

而不是比较矩阵==,使用< tolerance_limit

c = A*pinv(A);
d = pinv(A)*A;

(c-d) < 1e-10
Run Code Online (Sandbox Code Playgroud)

边注:

x = A^-1*b不应该解决x = inv(A)*b;,而是x = A \ b;看到链接Shai发布的解释.


Dan*_*iel 2

浮点运算有一定的精度,不能依赖等式。为了避免此类错误,您可以尝试使用 matlab 的符号工具箱。

非常简单的八度代码行来演示问题:

>>> (1/48)*48==(1/49)*49
ans = 0
>>> (1/48)*48-(1/49)*49
ans =  1.1102e-16
>>>
Run Code Online (Sandbox Code Playgroud)