使用 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 反馈 — 无需本地设置。
此课程中的所有课时
- 矩阵乘法与行列式
- 使用 solve() 求解线性系统
- 特征值与特征向量
- SVD、QR 与 Cholesky 分解