0Pricing
R Academy · 课时

使用 solve() 求解线性系统

求解 Ax = b 方程组并计算矩阵的逆

使用 solve() 求解线性系统 是 CoddyKit 上的免费 R Academy 课时。 这是第 2 节课,共 4 节。 你可以在下方免费阅读本课时的完整内容 — 然后在浏览器中使用内置代码编辑器和全天候 AI 导师进行实践。 这是 R Academy 学习路径的一部分,你的进度在网页和 CoddyKit 应用中同步。 R Academy 课程共包含 4 节课。

线性系统:Ax = b

线性方程组可以写成 Ax = b,其中 A 是系数矩阵,x 是未知向量,b 是等式右端项。通过解析方法求解 x,意味着计算 x = A⁻¹b。

# System of equations:
# 2x + y = 5
# x + 3y = 7

# Matrix form: A %*% x = b
A <- matrix(c(2, 1,
              1, 3), nrow = 2, byrow = TRUE)
b <- c(5, 7)

# What are A and b?
print(A)
print(b)
cat('We want to find x such that A %*% x = b')

solve(A, b):直接求解

solve(A, b) 求解 Ax = b 中的 x。它在内部使用 LU 分解,这比显式计算 A⁻¹ 再与 b 相乘具有更好的数值稳定性和效率。

A <- matrix(c(2, 1,
              1, 3), nrow = 2, byrow = TRUE)
b <- c(5, 7)

# Solve Ax = b
x <- solve(A, b)
print(x)  # x[1] = ?, x[2] = ?

# Verify: A %*% x should equal b
residual <- A %*% x - b
print(residual)  # Should be near zero

# Manual check:
# 2*(8/5) + (9/5) = 16/5 + 9/5 = 25/5 = 5 ✓
# 1*(8/5) + 3*(9/5) = 8/5 + 27/5 = 35/5 = 7 ✓

solve(A):矩阵逆

只使用一个参数调用 solve(A) 时,会返回 A 的逆矩阵 A⁻¹,满足 A %*% A⁻¹ = I。不要用它来求解 Ax=b——请直接使用 solve(A,b),这样更快且更稳定。

A <- matrix(c(4, 3,
              3, 2), nrow = 2, byrow = TRUE)

# Compute inverse
A_inv <- solve(A)
print(A_inv)

# Verify: A %*% A_inv = I
A %*% A_inv  # Should be identity matrix
round(A %*% A_inv, 10)

# Also A_inv %*% A = I
round(A_inv %*% A, 10)

# det(A) != 0 required for invertibility
det(A)  # -1 (nonzero, so invertible)

检查您的解

请始终通过计算 A %*% x - b 来验证解。由于浮点运算,残差不会恰好为零,但应接近机器精度(~1e-15)。使用 norm() 可得到单个残差的大小。

A <- matrix(c(3, -1, 2,
              1,  4, 0,
             -2,  1, 5), nrow = 3, byrow = TRUE)
b <- c(1, 2, 3)

# Solve
x <- solve(A, b)
cat('Solution x:\n'); print(x)

# Residual check
residual <- A %*% x - b
cat('Residual vector:\n'); print(residual)

# Residual norm (should be near 0)
resid_norm <- sqrt(sum(residual^2))
cat('Residual norm:', resid_norm, '\n')
# Expected: something like 2e-16

条件数:kappa()

条件数 A 衡量解对 b 中扰动的敏感程度。条件数较大意味着 b 中的微小误差会导致 x 中出现较大误差——该系统是病态的。

# Well-conditioned matrix
A_good <- matrix(c(2, 1, 1, 3), nrow = 2, byrow = TRUE)
kappa(A_good)  # Small -> good

# Ill-conditioned (nearly singular) matrix
A_bad <- matrix(c(1.000, 1.001,
                  1.001, 1.002), nrow = 2, byrow = TRUE)
kappa(A_bad)  # Very large -> bad!

# Rule of thumb: kappa > 1/machine_epsilon is trouble
.Machine$double.eps  # ~2.2e-16

# For A_bad: you lose about log10(kappa) digits of precision
cat('Digits lost:', log10(kappa(A_bad)), '\n')

使用 backsolve() 求解上三角系统

backsolve(R, b) 使用回代法求解 Rx = b,其中 R 是上三角矩阵。对于三角系统,它比通用的 solve() 快得多——复杂度为 O(n²),而不是 O(n³)。

# Upper triangular system: Rx = b
# 2x + 3y + z = 14
#      5y + 2z = 13
#           4z = 8
R <- matrix(c(2, 3, 1,
              0, 5, 2,
              0, 0, 4), nrow = 3, byrow = TRUE)
b <- c(14, 13, 8)

# Solve using back-substitution
x <- backsolve(R, b)
print(x)  # z=2, y=(13-4)/5=9/5, x=(14-3*9/5-2)/2

# Verify
all.equal(as.vector(R %*% x), b)  # TRUE

# Compare with general solve
x_general <- solve(R, b)
all.equal(x, x_general)  # TRUE (same result)

使用 forwardsolve() 求解下三角系统

forwardsolve(L, b) 使用前代法求解 Lx = b,其中 L 是下三角矩阵。它与 backsolve() 互为补充,两者共同构成 LU 分解求解器的基础。

# Lower triangular system: Lx = b
# 3x            = 6
# 2x + 4y       = 10
# x  + 2y + 5z  = 16
L <- matrix(c(3, 0, 0,
              2, 4, 0,
              1, 2, 5), nrow = 3, byrow = TRUE)
b <- c(6, 10, 16)

# Solve using forward-substitution
x <- forwardsolve(L, b)
print(x)  # x=2, y=(10-4)/4=1.5, z=(16-2-3)/5=2.2

# Verify
all.equal(as.vector(L %*% x), b)  # TRUE

# Use case: solving L*U*x = b
# forwardsolve(L, b) -> y, then backsolve(U, y) -> x

多个右端项

当 B 是矩阵时,solve(A, B) 会同时求解 AX = B,其中 B 的每一列对应一个解。这比对每一列分别调用 solve(A, b) 更高效。

A <- matrix(c(2, 1,
              1, 3), nrow = 2, byrow = TRUE)

# Solve for two right-hand sides simultaneously
B <- matrix(c(5, 7,   # first system
              3, 1),  # second system
            nrow = 2, byrow = TRUE)

# X[:,1] solves Ax = B[:,1]
# X[:,2] solves Ax = B[:,2]
X <- solve(A, B)
print(X)

# Verify both solutions
A %*% X  # Should equal B
all.equal(A %*% X, B)  # TRUE

检测奇异矩阵

对奇异矩阵调用 solve(A) 会抛出错误。求解前请检查 det(A) 或 rcond(A)(条件数的倒数)。使用 tryCatch() 可编写更健壮的代码。

# Singular matrix (rows are linearly dependent)
S <- matrix(c(1, 2,
              2, 4), nrow = 2, byrow = TRUE)
det(S)     # 0 -> singular
kappa(S)   # Inf

# Safe solve with tryCatch
safe_solve <- function(A, b) {
  tryCatch(
    solve(A, b),
    error = function(e) {
      cat('Matrix is singular or nearly so:\n')
      cat(e$message, '\n')
      return(NULL)
    }
  )
}

result <- safe_solve(S, c(1, 2))
print(result)  # NULL

使用 solve() 求解最小二乘问题

对于超定系统(方程多于未知数),不存在精确解。最小二乘解会最小化 ||Ax - b||²。它求解正规方程:A'Ax = A'b。

# Overdetermined: 4 equations, 2 unknowns (y = a + b*x)
set.seed(1)
x_vals <- c(1, 2, 3, 4)
y_vals <- c(2.1, 4.0, 5.9, 8.2)  # approx y = 0 + 2x
A <- cbind(1, x_vals)  # Design matrix (4x2)
b <- y_vals

# Normal equations: (A'A) beta = A'b
AtA <- crossprod(A)
Atb <- crossprod(A, b)
beta_ols <- solve(AtA, Atb)
print(beta_ols)  # Intercept ~0.1, slope ~2.0

# Residual sum of squares
y_hat <- A %*% beta_ols
rss <- sum((b - y_hat)^2)
cat('RSS:', rss)

使用 qr.solve() 提高稳定性

对于病态或超定系统,qr.solve(A, b) 的数值稳定性优于 solve()。它使用 QR 分解而不是 LU 分解。lm() 在内部使用此方法。

# For overdetermined system, qr.solve is preferred
x_vals <- c(1, 2, 3, 4, 5)
y_vals <- c(1.9, 4.1, 6.0, 7.8, 10.1)
A <- cbind(1, x_vals)
b <- y_vals

# qr.solve handles overdetermined systems directly
beta_qr <- qr.solve(A, b)
print(beta_qr)  # intercept, slope

# Equivalent to:
beta_lm <- coef(lm(y_vals ~ x_vals))
all.equal(beta_qr, beta_lm, check.names = FALSE)  # TRUE

# For well-determined square systems, solve() is fine
# For overdetermined or ill-conditioned: use qr.solve()

快速检查

请测试您对 R 中线性系统求解的理解。

回顾:求解线性系统

要点:solve(A, b) 直接求解 Ax=b(推荐)。solve(A) 计算 A⁻¹(不要用它来求解方程)。backsolve(R,b) 和 forwardsolve(L,b) 适合快速求解三角系统。kappa(A) 衡量条件性——值较大表示解不稳定。对于超定系统,请使用 qr.solve()。请始终使用 A %*% x - b 验证结果。

# Summary of solve() functions:
A <- matrix(c(3, 1, 1, 2), nrow = 2)
b <- c(9, 8)

# Solve Ax = b
x <- solve(A, b)
print(x)  # c(2, 3)

# Check conditioning
kappa(A)  # Small -> well-conditioned

# Verify
max(abs(A %*% x - b))  # Near zero

# For triangular systems:
R <- matrix(c(2, 3, 0, 4), nrow = 2, byrow = TRUE)
backsolve(R, c(8, 4))  # x=c(1, 1)

# For least squares (overdetermined):
# qr.solve(design_matrix, y)

常见问题解答

「使用 solve() 求解线性系统」课时是免费的吗?

是的 — 「使用 solve() 求解线性系统」的完整文本可在网页上免费阅读。要进行交互式练习(内置代码编辑器和全天候 AI 导师)并解锁 R Academy 课程的其余内容,请升级到 CoddyKit PRO。 R Academy 课程共包含 4 节课。

「使用 solve() 求解线性系统」这节课中我会学到什么?

求解 Ax = b 方程组并计算矩阵的逆 你通过在浏览器中直接运行的动手代码来练习 R Academy,全天候 AI 导师会在你学习这节课的过程中回答你的问题。

学习 R Academy 需要有经验吗?

无需任何先前经验。CoddyKit 上的 R Academy 课程适合初学者到高级学习者,你可以从这里开始或从头开始,按照自己的节奏学习。 这是第 2 节课,共 4 节。

「使用 solve() 求解线性系统」课时需要多长时间?

大多数 CoddyKit 课程大约需要 5–10 分钟。每节课都很精短且互动,所以你能稳步进步,并在网页和应用中从离开的地方继续。

我能在这节 R Academy 课中编写并运行代码吗?

能。每节 R Academy 课都包含内置代码编辑器,你可以在浏览器中直接编写并运行真实代码,并获得即时 AI 反馈 — 无需本地设置。

此课程中的所有课时

  1. 矩阵乘法与行列式
  2. 使用 solve() 求解线性系统
  3. 特征值与特征向量
  4. SVD、QR 与 Cholesky 分解
← 返回 R Academy