为什么我在 R studio 桌面和 R studio 云上得到相同算法的不同答案？

Question

我在我的 R studio 桌面（M1 pro MacBook Pro）上运行这个来自 class 的简单代码，答案与我在 R studio 云上得到的答案不同。在 MacBook 上，我得到 239.3093 与云上的 237.7821。我错过了什么吗？

A = matrix(NA,nrow = 50,ncol = 50)
for (i in 1:50) {
  for (j in 1:50) {
    A[i,j] = sin(i) + cos(j)
  }
}

pivot.above = function(A, r, c) {
  A[r, ] = A[r, ]/A[r,c]
  for (i in (r-1):1) {
    A[i, ] = A[i, ]-A[i,c]*A[r, ]
  }
   A
}

for (i in 0:48) {
  j = 50-i
  A = pivot.above(A,j,j)
}

A[1, ] = A[1, ]/A[1,1]

print(sum(A))

Answer 1

这个矩阵有一个非常大的条件数。来自 Wikipedia:

If the condition number is very large, then the matrix is said to be ill-conditioned. Practically, such a matrix is almost singular, and the computation of its inverse, or solution of a linear system of equations is prone to large numerical errors.

kappa(A)
[1] 1.352794e+19

这意味着微小的数值差异，例如编译器、系统线性代数库、操作系统、芯片等的差异，都会导致non-negligible结果的差异 — 如您所见。

正如@DaveArmstrong 评论的那样，OP 在这里所做的并不是完全维基百科文章中提到的任何操作（线性系统的求逆或求解）；我们通过高斯消元将矩阵简化为下三角形式。然而，它们是密切相关的（一旦我们将矩阵简化为三角形形式，我们就完成了求解线性系统的大部分工作，因为我们现在可以使用简单的 back-substitution）。我不知道为什么减少矩阵的第一列的总和应该是一个特别敏感的量，但我并不奇怪条件数与这种敏感性有关。在 these notes 中，Luke Tierney（R-core 成员 ...）描述了为什么简单的高斯消去法不稳定，以及为什么更复杂的 部分主元 方法更稳定.

这是一个非常合理的问题，但也是您在进行数值计算时必须学习的东西 about/get。问问你的导师实际上是一件好事。

Answer 2

抱歉，这不是一个真正的答案，但也许它是一个开始的地方。在回答 Ben 对我的问题的回答时，我了解到该算法不仅仅是将数字相加，但据我所知，它只是矢量化标量算术（加法、减法、乘法、除法）——没有矩阵运算继续。输入和输出以矩阵形式组织的事实是偶然的。考虑一下我认为在矩形数据框而不是方矩阵中组织的数据的一组等效操作。

B <- expand.grid(
  row = 1:50, 
  col = 1:50)

B <- B %>% 
  mutate(A = sin(row) + cos(col))

我重写了 pivot.above 函数来处理这个 long-form 数据：

pivot.above.long <- function(A, rowind, colind){
  A$A[which(A$row == rowind)]  <- A$A[which(A$row == rowind)]/A$A[which(A$row == rowind & A$col == colind)]
  for(i in (rowind - 1):1){
    A$A[which(A$row == i)] = A$A[which(A$row == i)]-A$A[which(A$row == i & A$col == colind)]*A$A[which(A$row == rowind)]
  }
  A
}

然后，我将其应用于 long-form 数据：

for (i in 0:48) {
  j = 50-i
  B = pivot.above.long(B,j,j)
}

B$A[which(B$row == 1)] <- B$A[which(B$row == 1)]/B$A[which(B$row == 1 & B$col == 1)]

而且，我得到的结果与我从矩阵运算得到的结果完全相同 - 原始问题中返回的两个值之一。

print(sum(B$A))
# [1] 239.3093

为什么我在 R studio 桌面和 R studio 云上得到相同算法的不同答案？

Why am I getting different answer for same algorithm on R studio desktop vs R studio cloud?

r

linear-algebra

apple-m1