선형 시스템(Linear Systems)

Ax = b를 푸는 일은 지금도 신경망(neural network)을 움직이는 가장 오래된 수학 문제입니다.

유형: Build 언어: Python 선수 지식: Phase 1, Lesson 01(Linear Algebra Intuition), 02(Vectors & Matrices), 03(Matrix Transformations) 예상 시간: 약 120분

학습 목표

• 부분 피벗팅(partial pivoting)과 후진 대입(back substitution)을 사용하는 가우스 소거법(Gaussian elimination)으로 Ax = b를 풉니다.
• LU, QR, 촐레스키 분해(Cholesky decomposition)로 행렬(matrix)을 분해(factorization)하고 각 방법이 언제 적절한지 설명합니다.
• 최소제곱(least squares)의 정규방정식(normal equation)을 유도하고 선형 회귀(linear regression), 능형 회귀(ridge regression)와 연결합니다.
• 조건수(condition number)로 악조건 시스템(ill-conditioned system)을 진단하고 정규화(regularization)로 안정화합니다.

문제

선형 회귀(linear regression)를 학습할 때마다 선형 시스템(linear system)을 풉니다. 최소제곱 적합(least-squares fit)을 계산할 때마다 선형 시스템을 풉니다. 신경망 계층(neural network layer)이 y = Wx + b를 계산할 때는 선형 시스템의 한쪽을 평가(evaluate)하고 있는 것입니다. 정규화(regularization)를 추가하면 시스템을 수정하는 것이고, 가우시안 프로세스(Gaussian process)를 사용하면 행렬을 분해(factor)합니다. 마할라노비스 거리(Mahalanobis distance)를 위해 공분산 행렬(covariance matrix)을 역행렬로 변환(invert)할 때도 선형 시스템을 풉니다.

방정식(equation) Ax = b는 어디에나 등장합니다. A는 이미 알려진 계수 행렬(known coefficient matrix)이고, b는 이미 알려진 출력 벡터(known output vector)이며, x는 찾고 싶은 미지수 벡터(unknown vector)입니다. 선형 회귀에서는 A가 데이터 행렬(data matrix), b가 목표 벡터(target vector), x가 가중치 벡터(weight vector)입니다. 전체 모델(model)은 "Ax가 b에 최대한 가까워지는 x를 찾는다"로 줄어듭니다.

이 lesson에서는 이 방정식을 푸는 주요 방법(method)을 처음부터 구현합니다. 어떤 방법이 빠르고 어떤 방법이 안정적인지, 어떤 방법은 정방 시스템(square system)에서만 작동하고 어떤 방법은 과결정 시스템(overdetermined system)까지 다루는지, 행렬의 조건수가 왜 답을 신뢰할 수 있는지를 결정하는지 이해합니다.

개념

Ax = b의 기하학적 의미

선형 방정식 시스템(linear equation system)에는 기하학적 해석(geometric interpretation)이 있습니다. 각 방정식은 초평면(hyperplane)을 정의하며, 해(solution)는 모든 초평면이 만나는 점(point) 또는 점의 집합입니다.

2x + y = 5          2차원 평면 위의 두 직선입니다.
x - y  = 1          x=2, y=1에서 교차합니다.

graph LR
    A["2x + y = 5"] --- S["해: (2, 1)"]
    B["x - y = 1"] --- S

세 가지 경우가 일어날 수 있습니다.

graph TD
    subgraph "유일한 해"
        A1["직선들이 한 점에서 교차"]
    end
    subgraph "해 없음"
        A2["직선들이 평행 -- 교차하지 않음"]
    end
    subgraph "무한한 해"
        A3["직선들이 일치 -- 모든 점이 해"]
    end

행렬 형태(matrix form)에서 "유일한 해"는 A가 역행렬을 갖는(invertible) 행렬이라는 뜻이고, "해 없음"은 시스템이 모순(inconsistent)이라는 뜻이며, "무한한 해"는 A가 영공간(null space)을 갖는다는 뜻입니다. 대부분의 ML 문제는 "정확한 해가 없는" 범주(category)에 속합니다. 방정식(데이터 점)이 미지수(파라미터)보다 많기 때문입니다. 여기서 최소제곱(least squares)이 등장합니다.

열 관점과 행 관점(Column Picture vs Row Picture)

Ax = b를 읽는 방법은 두 가지입니다.

행 관점(Row picture). A의 각 행(row)은 하나의 방정식을 정의합니다. 각 방정식은 초평면이며, 해는 모든 초평면이 만나는 지점입니다.

열 관점(Column picture). A의 각 열(column)은 벡터입니다. 이때 질문은 "A의 열들을 어떤 선형 결합(linear combination)으로 더하면 b가 되는가?"가 됩니다.

A = | 2  1 |    b = | 5 |
    | 1 -1 |        | 1 |

행 관점: 2x + y = 5와 x - y = 1을 동시에 풉니다.

열 관점: 다음을 만족하는 x1, x2를 찾습니다.
  x1 * [2, 1] + x2 * [1, -1] = [5, 1]
  2 * [2, 1] + 1 * [1, -1] = [4+1, 2-1] = [5, 1]   확인.

열 관점이 더 본질적(fundamental)입니다. b가 A의 열공간(column space) 안에 있으면 시스템에는 해가 있습니다. 그렇지 않다면 열공간 안에서 b에 가장 가까운 점을 찾습니다. 이 가장 가까운 점이 바로 최소제곱 해(least-squares solution)입니다.

가우스 소거법(Gaussian Elimination)

가우스 소거법은 Ax = b를 상삼각 시스템(upper triangular system) Ux = c로 변환하고 후진 대입(back substitution)으로 푸는 가장 직접적인 방법입니다.

알고리즘:

1. 각 열 k(피벗 열)에 대해:
   a. 행 k 이하에서 열 k의 절댓값이 가장 큰 항목을 찾습니다. (부분 피벗팅)
   b. 그 행을 행 k와 맞바꿉니다.
   c. k 아래의 각 행 i에 대해:
      - 승수 m = A[i][k] / A[k][k]를 계산합니다.
      - 행 i에서 m 배의 행 k를 뺍니다.
2. 후진 대입: 마지막 방정식부터 위로 올라가며 풉니다.

예시:

원래 형태:
| 2  1  1 | 8 |       R2 = R2 - (2)R1     | 2  1   1 |  8 |
| 4  3  3 |20 |  -->  R3 = R3 - (1)R1 --> | 0  1   1 |  4 |
| 2  3  1 |12 |                            | 0  2   0 |  4 |

                       R3 = R3 - (2)R2     | 2  1   1 |  8 |
                                       --> | 0  1   1 |  4 |
                                           | 0  0  -2 | -4 |

후진 대입:
  -2 * x3 = -4    -->  x3 = 2
  x2 + 2  = 4     -->  x2 = 2
  2*x1 + 2 + 2 = 8 --> x1 = 2

가우스 소거법의 비용은 O(n^3)입니다. 1000x1000 시스템이면 약 10억 회의 부동소수점 연산(floating-point operation)입니다. 빠르긴 하지만, 같은 A로 여러 시스템을 풀어야 한다면 더 나은 방법이 있습니다.

부분 피벗팅(Partial Pivoting)이 중요한 이유

피벗팅 없이 가우스 소거법을 수행하면 실패하거나 엉뚱한 값을 만들 수 있습니다. 피벗 원소(pivot element)가 0이면 0으로 나누게 되고, 피벗이 작으면 반올림 오차(rounding error)를 증폭시킵니다.

나쁜 피벗:                       부분 피벗팅 적용:
| 0.001  1 | 1.001 |            먼저 행을 맞바꿉니다.
| 1      1 | 2     |            | 1      1 | 2     |
                                 | 0.001  1 | 1.001 |
m = 1/0.001 = 1000              m = 0.001/1 = 0.001
R2 = R2 - 1000*R1               R2 = R2 - 0.001*R1
| 0.001  1     | 1.001   |      | 1      1     | 2     |
| 0     -999   | -999.0  |      | 0      0.999 | 0.999 |

x2 = 1.000 (정답)               x2 = 1.000 (정답)
x1 = (1.001 - 1)/0.001          x1 = (2 - 1)/1 = 1.000 (정답)
   = 0.001/0.001 = 1.000        승수가 작으므로 수치적으로 안정적입니다.

정밀도가 제한된 부동소수점 연산에서는 피벗팅을 하지 않은 버전이 유효 숫자(significant digit)를 잃을 수 있습니다. 부분 피벗팅은 오차 증폭을 줄이기 위해 그 시점에서 사용할 수 있는 가장 큰 피벗을 항상 선택합니다.

LU 분해(LU Decomposition)

LU 분해는 A를 하삼각 행렬(lower triangular matrix) L과 상삼각 행렬(upper triangular matrix) U로 분해합니다. 즉 A = LU입니다. L은 가우스 소거법에서 나온 승수를 저장하고, U는 소거(elimination)의 결과입니다.

A = L @ U

| 2  1  1 |   | 1  0  0 |   | 2  1   1 |
| 4  3  3 | = | 2  1  0 | @ | 0  1   1 |
| 2  3  1 |   | 1  2  1 |   | 0  0  -2 |

왜 그냥 소거하지 않고 분해할까요? L과 U가 있으면 새로운 b에 대한 Ax = b를 O(n^2)만에 풀 수 있기 때문입니다.

Ax = b
LUx = b
y = Ux라고 두면:
  Ly = b    (전진 대입, O(n^2))
  Ux = y    (후진 대입, O(n^2))

O(n^3) 비용은 분해할 때 한 번만 지불합니다. 이후의 풀이(solve)는 O(n^2)입니다. 같은 A에 대해 b 벡터가 1000개 있다면 LU 분해가 전체 작업량을 크게 줄여 줍니다.

부분 피벗팅을 사용하면 PA = LU가 됩니다. P는 행 교환을 기록하는 순열 행렬(permutation matrix)입니다.

QR 분해(QR Decomposition)

QR 분해는 A를 직교 행렬(orthogonal matrix) Q와 상삼각 행렬 R로 분해합니다. 즉 A = QR입니다.

직교 행렬은 Q^T Q = I를 만족합니다. 그 열들은 정규 직교(orthonormal) 벡터이며, Q를 곱해도 길이와 각도가 보존됩니다.

A = Q @ R

Q는 정규 직교 열을 가집니다: Q^T Q = I
R은 상삼각 행렬입니다

Ax = b를 풀려면:
  QRx = b
  Rx = Q^T b    (역행렬 없이 Q^T만 곱하면 됨)
  후진 대입으로 x를 구함.

QR 분해는 최소제곱 문제를 풀 때 LU보다 수치적으로 안정적입니다. 그램-슈미트 과정(Gram-Schmidt process)은 Q를 한 열씩 차례로 만들어 갑니다.

A의 열이 a1, a2, ... 라고 할 때:

q1 = a1 / ||a1||

q2 = a2 - (a2 . q1) * q1        (q1 방향 사영을 뺌)
q2 = q2 / ||q2||                (정규화)

q3 = a3 - (a3 . q1) * q1 - (a3 . q2) * q2
q3 = q3 / ||q3||

R[i][j] = qi . aj    (i <= j 일 때)

각 단계는 이전 q 벡터 방향 성분을 모두 제거하여 새로운 직교 방향만 남깁니다.

촐레스키 분해(Cholesky Decomposition)

A가 대칭(symmetric, A = A^T)이고 양의 정부호(positive definite, 모든 고윳값이 양수)이면 A = L L^T로 분해할 수 있습니다. 이것이 촐레스키 분해입니다.

A = L @ L^T

| 4  2 |   | 2  0 |   | 2  1 |
| 2  5 | = | 1  2 | @ | 0  2 |

L[i][i] = sqrt(A[i][i] - sum(L[i][k]^2 for k < i))
L[i][j] = (A[i][j] - sum(L[i][k]*L[j][k] for k < j)) / L[j][j]    (i > j 일 때)

촐레스키 분해는 LU보다 두 배 정도 빠르고 저장 공간도 절반만 필요합니다. 대칭 양의 정부호 행렬에서만 작동하지만, ML에서는 이런 행렬이 자주 등장합니다.

• 공분산 행렬(covariance matrix)은 대칭 양의 준정부호(symmetric positive semi-definite)이며, 정규화를 더하면 양의 정부호가 됩니다.
• 가우시안 프로세스의 커널 행렬(kernel matrix)은 대칭 양의 정부호입니다.
• 볼록 함수(convex function)의 최솟값에서 헤시안(Hessian)은 대칭 양의 정부호입니다.
• A^T A는 항상 대칭 양의 준정부호입니다.

가우시안 프로세스에서는 커널 행렬 K를 촐레스키 분해로 분해하고 K alpha = y를 풀어 예측 평균(predictive mean)을 얻습니다. 촐레스키 인자(Cholesky factor)는 주변 가능도(marginal likelihood)의 로그 행렬식도 제공합니다. 즉 log det(K) = 2 * sum(log(diag(L)))입니다.

최소제곱(Least Squares): Ax = b에 정확한 해가 없을 때

A가 m x n이고 m > n이면 방정식이 미지수보다 많은 과결정 시스템(overdetermined system)입니다. 정확한 해는 없으며, 대신 제곱 오차(squared error)를 최소화합니다.

||Ax - b||^2 를 최소화한다

이는 잔차 제곱의 합과 같습니다:
  sum((A[i,:] @ x - b[i])^2 for i in range(m))

이 최소값을 만드는 점은 정규방정식(normal equation)을 만족합니다.

A^T A x = A^T b

유도 과정: ||Ax - b||^2 = (Ax - b)^T (Ax - b) = x^T A^T A x - 2 x^T A^T b + b^T b로 전개합니다. x에 대한 그래디언트를 구해 0으로 두면 2 A^T A x - 2 A^T b = 0이 됩니다.

원래 시스템 (과결정, 방정식 4개, 미지수 2개):
| 1  1 |         | 3 |
| 1  2 | x     = | 5 |       모든 방정식을 만족하는 정확한 x는 없습니다.
| 1  3 |         | 6 |
| 1  4 |         | 8 |

정규방정식:
A^T A = | 4  10 |    A^T b = | 22 |
        | 10 30 |            | 63 |

해: x = [1.5, 1.7]

이것이 곧 선형 회귀입니다. x[0]은 절편(intercept), x[1]은 기울기(slope)입니다.

정규방정식(Normal Equation) = 선형 회귀(Linear Regression)

이 연결은 정확하게 일치합니다. 선형 회귀에서 데이터 행렬 X는 샘플(sample)마다 한 행, 특성(feature)마다 한 열을 갖고, 목표 벡터 y는 샘플마다 하나의 항목을 갖습니다. 가중치 벡터 w는 다음을 만족합니다.

X^T X w = X^T y
w = (X^T X)^(-1) X^T y

이것이 선형 회귀의 닫힌 형태 해(closed-form solution)입니다. sklearn.linear_model.LinearRegression.fit() 호출은 이 식을 직접 계산하거나 QR/SVD를 거치는 같은 결과의 방법을 사용합니다.

행렬에 정규화 항 lambda * I를 더하면 능형 회귀(ridge regression)가 됩니다.

(X^T X + lambda * I) w = X^T y
w = (X^T X + lambda * I)^(-1) X^T y

정규화는 행렬의 조건을 더 좋게 만들어 역행렬 계산을 정확하게 만들고, 가중치를 0 쪽으로 줄여 과적합(overfitting)을 막습니다. lambda > 0이면 X^T X + lambda * I는 항상 대칭 양의 정부호이므로 촐레스키 분해로 풀 수 있습니다.

의사역행렬(Pseudoinverse, Moore-Penrose)

의사역행렬 A+는 행렬 역연산을 정방이 아니거나 특이한(non-square, singular) 행렬로 일반화합니다. 임의의 행렬 A에 대해 다음처럼 씁니다.

x = A+ b

여기서 A+ = V Sigma+ U^T    (SVD로 계산)

Sigma+는 0이 아닌 특잇값(singular value)의 역수를 취하고 전치(transpose)해 만듭니다. A = U Sigma V^T이면 A+ = V Sigma+ U^T입니다.

A = U Sigma V^T        (SVD)

Sigma = | 5  0 |       Sigma+ = | 1/5  0  0 |
        | 0  2 |                | 0  1/2  0 |
        | 0  0 |

A+ = V Sigma+ U^T

의사역행렬은 최소 노름 최소제곱 해(minimum-norm least-squares solution)를 제공합니다. 시스템 상태에 따라 다음처럼 작동합니다.

• 유일한 해가 있는 경우: A+ b가 그 해를 그대로 돌려줍니다.
• 해가 없는 경우: A+ b가 최소제곱 해를 돌려줍니다.
• 해가 무한히 많은 경우: A+ b가 그중에서 ||x||가 가장 작은 해를 돌려줍니다.

NumPy의 np.linalg.lstsq와 np.linalg.pinv는 모두 내부적으로 SVD를 사용합니다.

조건수(Condition Number)

조건수는 입력의 작은 변화에 해가 얼마나 민감한지를 측정합니다. 행렬 A의 조건수는 다음과 같습니다.

kappa(A) = ||A|| * ||A^(-1)|| = sigma_max / sigma_min

여기서 sigma_max와 sigma_min은 각각 가장 큰 특잇값과 가장 작은 특잇값입니다.

좋은 조건 (kappa ~ 1):               악조건 (kappa ~ 10^15):
b의 작은 변화 -->                    b의 작은 변화 -->
x의 작은 변화                        x의 큰 변화

| 2  0 |   kappa = 2/1 = 2          | 1   1          |   kappa ~ 10^15
| 0  1 |   안전하게 풀 수 있음      | 1   1+10^(-15) |   해는 사실상 쓰레기 값

경험칙:

• kappa < 100: 안전하며 해가 정확합니다.
• kappa ~ 10^k: 부동소수점 연산에서 약 k 자리의 정밀도를 잃습니다.
• kappa ~ 10^16 (float64 기준): 해는 의미가 없습니다. 행렬이 사실상 특이(singular)합니다.

ML에서는 특성(feature)들이 거의 공선(collinear)일 때 악조건이 발생합니다. 정규화(lambda * I를 더하는 작업)는 조건수를 sigma_max / sigma_min에서 (sigma_max + lambda) / (sigma_min + lambda)로 개선합니다.

반복법(Iterative Method): 켤레 그래디언트(Conjugate Gradient)

미지수가 수백만 개에 이르는 매우 큰 희소(sparse) 시스템에서는 LU나 촐레스키 분해 같은 직접 풀이법(direct method)이 너무 비쌉니다. 반복법은 추정값을 여러 번의 반복(iteration)에 걸쳐 개선하면서 해를 근사합니다.

켤레 그래디언트(Conjugate Gradient; CG)는 A가 대칭 양의 정부호일 때 Ax = b를 풉니다. 정확한 산술 환경에서는 최대 n 회의 반복 안에 정확한 해를 찾지만, 고윳값이 한 곳에 모여 있으면 보통 그보다 훨씬 빨리 수렴합니다.

알고리즘 개요:
  x0 = 초기 추정값 (보통 0)
  r0 = b - A x0           (잔차)
  p0 = r0                 (탐색 방향)

  k = 0, 1, 2, ... 에 대해:
    alpha = (rk . rk) / (pk . A pk)
    x_{k+1} = xk + alpha * pk
    r_{k+1} = rk - alpha * A pk
    beta = (r_{k+1} . r_{k+1}) / (rk . rk)
    p_{k+1} = r_{k+1} + beta * pk
    if ||r_{k+1}|| < 허용 오차: 멈춤

CG가 쓰이는 곳:

• 대규모 최적화(Newton-CG 방법)
• 편미분방정식(PDE) 이산화 풀이
• 커널 행렬이 너무 커서 분해하기 어려운 커널 방법(kernel method)
• 다른 반복 해법(iterative solver)을 위한 사전 조정(preconditioning)

수렴 속도는 조건수에 따라 달라집니다. 조건이 좋은 시스템일수록 더 빨리 수렴합니다. 정규화가 도움이 되는 또 하나의 이유입니다.

전체 그림: 어떤 방법을 언제 쓸까

ML과의 연결

이 lesson의 모든 방법은 실제 운영 환경의 ML에 등장합니다.

선형 회귀(Linear regression). 닫힌 형태 해는 정규방정식 X^T X w = X^T y를 풉니다. n이 작으면 촐레스키 분해, 수치적 안정성이 중요하면 QR, 계수 결핍(rank-deficient) 가능성이 있으면 SVD를 사용합니다.

능형 회귀(Ridge regression). X^T X에 lambda * I를 더합니다. 정규화된 시스템 (X^T X + lambda * I) w = X^T y는 lambda > 0일 때 대칭 양의 정부호이므로 촐레스키 분해로 풀 수 있습니다.

가우시안 프로세스(Gaussian processes). 예측 평균은 커널 행렬 K에 대해 K alpha = y를 푸는 문제입니다. K에 대한 촐레스키 분해가 표준 접근입니다. 로그 주변 가능도(log marginal likelihood)는 log det(K) = 2 sum(log(diag(L)))를 사용합니다.

신경망 초기화(Neural network initialization). 직교 초기화(orthogonal initialization)는 QR 분해를 사용해 열들이 정규 직교인 가중치 행렬을 만듭니다. 이는 깊은 신경망에서 신호 붕괴(signal collapse)를 막는 데 도움이 됩니다.

사전 조정(Preconditioning). 대규모 최적화기는 켤레 그래디언트 해법의 사전 조정자로 불완전 촐레스키(incomplete Cholesky)나 불완전 LU(incomplete LU)를 사용합니다.

특성 공학(Feature engineering). X^T X의 조건수는 특성들이 공선인지 알려줍니다. kappa가 크면 특성을 제거하거나 정규화를 추가합니다.

만들어 보기

Step 1: 부분 피벗팅을 사용하는 가우스 소거법

import numpy as np

def gaussian_elimination(A, b):
    n = len(b)
    Ab = np.hstack([A.astype(float), b.reshape(-1, 1).astype(float)])

    for k in range(n):
        max_row = k + np.argmax(np.abs(Ab[k:, k]))
        Ab[[k, max_row]] = Ab[[max_row, k]]

        if abs(Ab[k, k]) < 1e-12:
            raise ValueError(f"행렬이 피벗 {k}에서 특이하거나 거의 특이합니다")

        for i in range(k + 1, n):
            m = Ab[i, k] / Ab[k, k]
            Ab[i, k:] -= m * Ab[k, k:]

    x = np.zeros(n)
    for i in range(n - 1, -1, -1):
        x[i] = (Ab[i, -1] - Ab[i, i+1:n] @ x[i+1:n]) / Ab[i, i]

    return x

Step 2: LU 분해

def lu_decompose(A):
    n = A.shape[0]
    L = np.eye(n)
    U = A.astype(float).copy()
    P = np.eye(n)

    for k in range(n):
        max_row = k + np.argmax(np.abs(U[k:, k]))
        if max_row != k:
            U[[k, max_row]] = U[[max_row, k]]
            P[[k, max_row]] = P[[max_row, k]]
            if k > 0:
                L[[k, max_row], :k] = L[[max_row, k], :k]

        for i in range(k + 1, n):
            L[i, k] = U[i, k] / U[k, k]
            U[i, k:] -= L[i, k] * U[k, k:]

    return P, L, U

Step 3: 촐레스키 분해

def cholesky(A):
    n = A.shape[0]
    L = np.zeros_like(A, dtype=float)

    for i in range(n):
        for j in range(i + 1):
            s = A[i, j] - L[i, :j] @ L[j, :j]
            if i == j:
                if s <= 0:
                    raise ValueError("행렬이 양의 정부호가 아닙니다")
                L[i, j] = np.sqrt(s)
            else:
                L[i, j] = s / L[j, j]

    return L

Step 4: 정규방정식으로 최소제곱 풀기

def least_squares_normal(A, b):
    AtA = A.T @ A
    Atb = A.T @ b
    return gaussian_elimination(AtA, Atb)

def ridge_regression(A, b, lam):
    n = A.shape[1]
    AtA = A.T @ A + lam * np.eye(n)
    Atb = A.T @ b
    L = cholesky(AtA)
    y = np.zeros(n)
    for i in range(n):
        y[i] = (Atb[i] - L[i, :i] @ y[:i]) / L[i, i]
    x = np.zeros(n)
    for i in range(n - 1, -1, -1):
        x[i] = (y[i] - L.T[i, i+1:] @ x[i+1:]) / L.T[i, i]
    return x

Step 5: 조건수

def condition_number(A):
    U, S, Vt = np.linalg.svd(A)
    return S[0] / S[-1]

사용하기

실제 데이터에서 선형 회귀와 능형 회귀를 함께 구성하면 다음과 같습니다.

np.random.seed(42)
X_raw = np.random.randn(100, 3)
w_true = np.array([2.0, -1.0, 0.5])
y = X_raw @ w_true + np.random.randn(100) * 0.1

X = np.column_stack([np.ones(100), X_raw])

w_ols = least_squares_normal(X, y)
print(f"OLS 가중치 (직접 구현): {w_ols}")

w_np = np.linalg.lstsq(X, y, rcond=None)[0]
print(f"OLS 가중치 (numpy):    {w_np}")
print(f"최대 차이: {np.max(np.abs(w_ols - w_np)):.2e}")

w_ridge = ridge_regression(X, y, lam=1.0)
print(f"능형 회귀 가중치 (직접 구현): {w_ridge}")

산출물 만들기

이 lesson의 최종 산출물은 다음과 같습니다.

• code/linear_systems.py: 가우스 소거법, LU 분해, 촐레스키 분해, 최소제곱, 능형 회귀를 처음부터 직접 구현한 코드를 포함합니다.
• outputs/prompt-linear-solver.md: 행렬의 성질에 따라 적절한 선형 시스템 해법을 추천해 주는 prompt입니다.

연습문제

1. (쉬움) [[1,2,3],[4,5,6],[7,8,10]] x = [6, 15, 27] 시스템을 직접 구현한 가우스 소거법, LU 풀이기, np.linalg.solve로 각각 풉니다. 세 답이 부동소수점 허용 오차 안에서 일치하는지 확인합니다.

2. (중간) 50x5 크기의 무작위 행렬 X와 목표 벡터 y = X @ w_true + noise를 생성합니다. 정규방정식, QR(np.linalg.qr), SVD(np.linalg.svd), np.linalg.lstsq로 w를 구합니다. 네 가지 해를 비교하고, X^T X의 조건수를 측정하여 어떤 방법을 더 신뢰할 수 있는지 설명합니다.

3. (중간) 두 열을 거의 동일하게 만들어 거의 특이한 행렬을 구성합니다. 예: 열 2 = 열 1 + 1e-10 * noise. 조건수를 계산한 뒤, 정규화 없이 그리고 0.01 * I를 더해 Ax = b를 각각 풉니다. 해와 잔차(residual)를 비교하고 정규화가 왜 도움이 되는지 설명합니다.

4. (어려움) 100x100 크기의 무작위 대칭 양의 정부호 행렬에 대해 켤레 그래디언트 알고리즘을 직접 구현합니다. 허용 오차 1e-8에 도달하기까지 걸린 반복 횟수를 세고, 이론적 최댓값인 n회 반복과 비교합니다.

5. (어려움) 크기 10, 50, 200, 500인 대칭 양의 정부호 행렬에서 촐레스키 풀이기, LU 풀이기, np.linalg.solve의 수행 시간을 측정합니다. 결과를 그래프로 그려 보고 촐레스키 분해가 LU보다 대략 2배 빠른지 확인합니다.

핵심 용어

더 읽을거리

• MIT 18.06: Linear Algebra (Gilbert Strang) — 선형 시스템과 행렬 분해를 다루는 대표 강의입니다.
• Numerical Linear Algebra (Trefethen & Bau) — 수치적 안정성, 조건수, 알고리즘이 실패하는 이유를 이해하는 데 쓰이는 표준 참고서입니다.
• Matrix Computations (Golub & Van Loan) — 모든 행렬 알고리즘을 망라한 백과사전 같은 참고 문헌입니다.
• 3Blue1Brown: Inverse Matrices — Ax = b를 푸는 일이 기하학적으로 무엇을 의미하는지 시각적으로 설명합니다.