r - 用 R 求解欠定线性系统

标签 r math linear-algebra numerical-methods

R 可以解决欠定线性系统:

A = matrix((1:12)^2,3,4,T)
B = 1:3
qr(A)$rank  # 3
qr.solve(A, B)  # solutions will have one zero, not necessarily the same one
# 0.1875 -0.5000  0.3125  0.0000
solve(qr(A, LAPACK = TRUE), B)
# 0.08333333 -0.18750000  0.00000000  0.10416667

(它给出了无数解决方案中的一个解决方案)。

但是,如果排名(此处为 2)低于行数(此处为 3),则它将不起作用:

A = matrix(c((1:8)^2,0,0,0,0),3,4,T)
B = c(1,2,0)
A
#      [,1] [,2] [,3] [,4]
# [1,]    1    4    9   16
# [2,]   25   36   49   64
# [3,]    0    0    0    0

qr.solve(A, B)  # Error in qr.solve(A, B) : singular matrix
solve(qr(A, LAPACK = TRUE), B)  # Error in qr.coef(a, b) : error code 3 

但是这个系统确实有解决方案!

我知道一般的解决方案是使用 SVD 或 A 的广义/伪逆(参见 this question 及其答案),但是:

是否有一种方法可以使用 solveqr.solve 自动将系统 AX=B 简化为仅等级(A )行,其中 qr.solve(C, D) 可以开箱即用?

示例:

C = matrix(c((1:8)^2),2,4,T)
D = c(1,2)
qr.solve(C, D)
# -0.437500  0.359375  0.000000  0.000000

最佳答案

qr.coefqr 似乎可以完成这项工作:

(A <- matrix(c((1:8)^2, 0, 0, 0, 0), nrow = 3, ncol = 4, byrow = TRUE))
#     [,1] [,2] [,3] [,4]
# [1,]    1    4    9   16
# [2,]   25   36   49   64
# [3,]    0    0    0    0
(B <- c(1, 2, 0))
# [1] 1 2 0
(X0 <- qr.coef(qr(A), B))
# [1] -0.437500  0.359375        NA        NA
X0[is.na(X0)] <- 0
X0
# [1] -0.437500  0.359375  0.000000  0.000000
# Verification:
A %*% X0
#      [,1]
# [1,]    1
# [2,]    2
# [3,]    0

第二个例子:

(A<-matrix(c(1, 2, 0, 0, 1, 2, 0, 0, 1, 2, 1, 0), nrow = 3, ncol = 4, byrow = TRUE))
#      [,1] [,2] [,3] [,4]
# [1,]    1    2    0    0
# [2,]    1    2    0    0
# [3,]    1    2    1    0
(B<-c(1, 1, 2))
# [1] 1 1 2
qr.solve(A, B)
# Error in qr.solve(A, B) : singular matrix 'a' in solve
(X0 <- qr.coef(qr(A), B))
# [1]  1 NA  1 NA
X0[is.na(X0)] <- 0
X0
# [1] 1 0 1 0
A %*% X0
#      [,1]
# [1,]    1
# [2,]    1
# [3,]    2

关于r - 用 R 求解欠定线性系统,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/53009467/

相关文章:

r - 从 4 个变量创建单列

c++ - Rcpp 移动平均 - 边界错误导致 fatal error

python - 厄密矩阵的 logm 函数返回非厄密矩阵

python - 当 x = y 时,Numpy 和 R 在线性回归中给出非零截距

matlab - 稀疏矩阵的线性代数库

r - 大量因果条件会在 QCA 分析中产生错误吗?

regex - 将字符分组后的字符替换为空

python - 我如何将矢量投影到由 Python 中的正交矢量定义的平面上?

c++ - 如何根据球击中 Racket 的位置使球以特定角度从 Racket 上反弹?

JAVA动圆与非动圆弹性碰撞