为什么在 R 中求解会给出奇怪的逆？答案

【问题标题】：Why does solve in R give weird inverse?为什么在 R 中求解会给出奇怪的逆？
【发布时间】：2014-09-22 07:33:01
【问题描述】：

M矩阵

     1000          793504755
793504755    629649829362825

这是我在 R 中的矩阵 M solve(M) 给出这个错误：

solve(M)
Error in solve.default(M) : 
  system is computationally singular: reciprocal condition number = 8.36282e-20

而 excel 中的 MMINV 为这个 2 x 2 矩阵提供所需的输出

R 是否使用某种分解算法进行逆运算？为什么它没有为我的矩阵提供输出？矩阵行列式也不为零。事实上它是一个很大的值。

【问题讨论】：

请提供 Excel 结果。我不相信它在合理的数字误差范围内是正确的。你为什么要计算这个矩阵的逆矩阵？也许你可以避免这一步？
@Roland，这确实很奇怪。您可以自己反转它并通过执行M <- matrix(c(1000, 793504755, 793504755, 629649829362825), 2, 2); M2 <- matrix(c(629649829362825, -793504755, -793504755, 1000), 2, 2);(1/(M[1,1]*M[2,2] - M[1,2]*M[1,2]))*M2 之类的操作来获得结果，这确实提供了所需的结果
你可以通过solve(m,tol=1e-40)得到an答案，但我不相信它！
@BenBolker，为什么不呢？它给出了正确的结果，就像你用手计算一样
@DavidArenburg 浮点精度是这里的一个问题。

标签： r matrix inverse matrix-inverse

【解决方案1】：

您遇到了麻烦，因为您的矩阵是病态的。您可以使用Matrix 包解决它，而无需担心中间步骤：

library(Matrix)
solve(Matrix(M))
2 x 2 Matrix of class "dsyMatrix"
              [,1]          [,2]
[1,]  1.899097e+04 -2.393303e-02
[2,] -2.393303e-02  3.016117e-08

但结果中存在一些数值不精确：

M %*% solve(Matrix(M))
2 x 2 Matrix of class "dgeMatrix"
            [,1] [,2]
[1,] 1.000000004    0
[2,] 0.001953125    1

这是因为倒数条件数小于双精度数。

rcond(M)
[1] 8.362817e-20
storage.mode(M)
[1] "double"
.Machine$double.eps
[1] 2.220446e-16

但是，gmp 包使用有理数提供了精确的精度，从而为您提供所需的结果：

library(gmp)
solve(as.bigz(M))
Big Rational ('bigq') 2 x 2 matrix:
     [,1]                  [,2]              
[1,] 932814562019/49118837 -5877813/245594185
[2,] -5877813/245594185    40/1326208599     
M %*% solve(matrix(as.bigz(c(M)),2))
Big Rational ('bigq') 2 x 2 matrix:
     [,1] [,2]
[1,] 1    0   
[2,] 0    1

【讨论】：

【解决方案2】：

这会遇到数字困难。您可以使用 David 建议的“手动”方法，也可以使用 Choleski 分解来获得一个合理的估计：

chol2inv(chol(M))
#              [,1]          [,2]
#[1,]  1.899097e+04 -2.393303e-02
#[2,] -2.393303e-02  3.016117e-08

chol2inv(chol(M)) %*% M
#             [,1]        [,2]
#[1,] 1.000000e+00 0.001953125
#[2,] 3.552714e-15 1.000000000

但是，我建议避免计算此矩阵的逆矩阵。

【讨论】：