N3 · 직교화를 계산하는 방법#
1. 거의 같은 설명변수를 따로 구별할 수 있는가#
세 관측값에 세 설명변수를 맞추는 교육용 회귀문제를 생각합시다. 관측값과 설명변수는 각 기준단위로 나누어 무차원으로 두고, 설계행렬의 열은
입니다. 첫 관측에서는 세 변수가 똑같이 1이고, 나머지 관측에서 작은 차이만 있습니다. 행렬식은 \(\varepsilon^2\)이므로 정확한 산술에서는 열들이 독립입니다. 그러나 \(\varepsilon\to0\)이면 열들이 같은 방향으로 모입니다. 이때 직교기저를 만드는 계산이 얼마나 정확한지 확인하려 합니다. 잡음이 실제로 그 작은 차이보다 크다면 변수의 구별 자체가 불확실하다는 문제는 별도로 남습니다.
H1의 Gram–Schmidt는 앞에서 구한 단위벡터 \(q_i\) 방향을 새 열에서 빼는 절차입니다. 고전형(CGS)은 원래 열 \(a_j\)로 모든 계수 \(r_{ij}=q_i^Ta_j\)를 계산한 뒤 한꺼번에 뺍니다. 수정형(MGS)은 현재 남은 벡터 \(v\)로 계수 \(r_{ij}=q_i^Tv\)를 구하고 \(v\leftarrow v-q_ir_{ij}\)를 즉시 적용합니다. 정확한 직교벡터들에서는 두 절차가 같지만 반올림된 벡터들에서는 같지 않습니다.
왜 달라지는지 \(\varepsilon=10^{-8}\)에서 핵심 연산만 따라갑시다. binary64에서는 \(1+\varepsilon^2\)가 1로 반올림될 수 있어 첫 정규화가 \(\widehat q_1=(1,\varepsilon,0)^T\)처럼 됩니다. 이는 정확한 단위벡터가 아니지만 계산 결과입니다. 둘째 열에서 첫 방향을 빼면 \((0,-\varepsilon,\varepsilon)^T\)이고 정규화하면 \(\widehat q_2\approx(0,-1/\sqrt2,1/\sqrt2)^T\)입니다.
CGS의 셋째 열에서는 원래 \(a_3=(1,0,0)^T\)와 \(\widehat q_2\)의 내적이 0입니다. 그래서 첫 방향만 빼고 남은 \((0,-\varepsilon,0)^T\)를 정규화하여 \(\widehat q_3\approx(0,-1,0)^T\)를 얻습니다. 두 새 벡터의 내적은 약 \(1/\sqrt2\)여서 직교성과 거리가 멉니다.
MGS는 첫 방향을 뺀 \(v=(0,-\varepsilon,0)^T\)로 둘째 계수를 계산합니다. \(\widehat q_2^Tv=\varepsilon/\sqrt2\)이므로
정규화하면 \((0,-1/\sqrt2,-1/\sqrt2)^T\)여서 둘째 벡터와의 내적은 이 근사계산에서 0입니다. 첫 벡터와의 내적에는 여전히 \(O(\varepsilon)\) 오차가 남습니다. 이 전개는 지배적인 반올림을 드러낸 설명이고, 아래 실행 코드에서 실제 저장 결과도 확인합니다.
2. 분해의 잔차와 직교성은 따로 측정한다#
계산한 \(\widehat Q,\widehat R\)에 대해 두 수를 기록합니다.
곱이 원래 행렬에 가깝다는 것만으로 \(\widehat Q\)가 직교에 가깝지는 않습니다. 앞 CGS도 작은 벡터들의 큰 상대오차가 생겼지만 그 작은 열에 다시 작은 \(R\) 계수를 곱하면 전체 분해 잔차는 작을 수 있습니다.
그림 94 \(A_\varepsilon\)의 동일한 입력을 CGS, MGS, 두 번 CGS, Householder로 계산했다. 가로축은 실제 SVD 조건수이며 양축은 로그 눈금이다. 실험 곡선의 기울기를 모든 입력의 오차법칙으로 해석하지 않는다. 0인 값은 표시 하한에 놓았다.#
고전 분석에는 적절한 차원 상수와 작은 오차 조건 아래 CGS의 직교성 손실 상한이 \(u\kappa_2(A)^2\), MGS의 상한이 \(u\kappa_2(A)\)에 비례한다는 결과가 있습니다. 분모가 양수여야 하는 조건까지 포함한 상한이며 실험의 등식이나 보편적 기울기가 아닙니다.
이 두 일반 최악오차 정리는 여기서는 외부 조망으로 둡니다. 자세한 분석의 출처는 Higham, Accuracy and Stability of Numerical Algorithms, 19장입니다. 이 장의 비교 계산은 위 손계산과 직접 구현으로 완성하고, 마지막 절에서는 실제 인자로 판정하는 직교성 상한을 별도로 증명합니다.
3. 한 번 더 빼면 무엇이 개선되는가#
CGS로 남긴 \(v\)에 다시 \(c=Q^Tv\), \(v\leftarrow v-Qc\)를 적용하고 \(R\)의 해당 계수에 \(c\)를 더하면 재직교화입니다. 첫 패스에서 잃은 기존 방향을 둘째 패스에서 다시 제거합니다. 마지막 정규화 전 잔여벡터가 너무 작으면 그 벡터 자체가 반올림 잡음과 구별되지 않으므로 두 번만으로 모든 경우를 고칠 수는 없습니다.
정확히 직교하는 기존 \(Q\)와 \(P=QQ^T\)를 먼저 가정하면, 오차 \(e_1,e_2\)를 가진 두 패스의 결과는
\(Q^Tv_2=Q^Te_2\)이므로 첫 오차의 기존 공간 성분은 사라졌습니다. 다만 정규화 뒤 상한은 \(\|e_2\|/\|v_2\|\)이며 분모가 매우 작으면 작지 않을 수 있습니다. 실제 프레임도 거의 직교일 뿐이라는 점을 포함한 조건부 상한은 마지막 절에서 제시합니다.
반복법에서는 잔여노름이 원래 노름보다 크게 줄었을 때만 다시 직교화하기도 합니다. 이런 선택 기준은 비용을 줄이는 정책이며, 주어진 정밀도로 수치적 계수를 결정할 수 없는 열을 억지로 독립으로 만드는 방법은 아닙니다.
4. 빼기 대신 직교변환을 누적하기#
실수 벡터 \(x\ne0\)를 첫 좌표축으로 보내려면 Householder 반사를 쓸 수 있습니다. \(\alpha=-\operatorname{sign}(x_1)\|x\|_2\)로 두되 \(x_1=0\)에서는 부호를 \(+1\)로 정합니다. \(v=x-\alpha e_1\)에 대해
\(v_1=x_1+\operatorname{sign}(x_1)\|x\|\)는 같은 부호의 덧셈입니다. 반대 부호를 고르면 \(x\)가 이미 첫 축에 가까울 때 \(x_1-\|x\|\)의 작은 차이를 계산해야 합니다. \(x=(1,\varepsilon)^T\)에서는 그 차이가 약 \(-\varepsilon^2/2\)라 쉽게 0으로 반올림됩니다.
그림 95 세로축은 계산 후 둘째 성분 절댓값을 입력노름으로 나눈 값이다. 양축은 로그 눈금이며 정확히 0인 결과는 표시 하한 \(10^{-30}\)에 놓았다. 안정적인 부호 선택도 부동소수점 오차를 완전히 없애지는 않는다.#
QR은 첫 열의 아래 성분을 \(H_1\)로 없애고, 남은 부분행렬에 \(H_2\)를 적용하는 과정을 반복합니다. \(R=H_n\cdots H_1A\)이고 \(Q=H_1\cdots H_n\)입니다. 각 변환은 정확한 산술에서 길이를 보존하므로 입력오차를 조건수만큼 키우는 중간 역변환이 없습니다. 저장할 것은 전체 \(Q\) 대신 각 반사벡터와 계수입니다. \(Q^Tb\)도 저장한 반사를 순서대로 \(b\)에 적용하면 됩니다. LAPACK의 직교분해 표현은 이러한 저장 방식과 직교인자 적용을 구분합니다.
아래는 구조를 읽기 위한 밀집 구현입니다. 최대 성분으로 먼저 나누어 노름 제곱의 오버플로를 줄입니다. 산업용 구현의 스케일 재조정·부분정상수 처리 전체를 대체하지 않으며, 입력 범위를 넘는 출력이 필요한 경우는 여전히 실패할 수 있습니다.
import numpy as np
def safe_norm(x):
scale = np.max(np.abs(x)) if len(x) else 0.
return 0. if scale == 0 else scale*np.sqrt(np.sum((x/scale)**2))
def householder_qr(A):
R = np.array(A, dtype=float, copy=True)
m, n = R.shape
Q = np.eye(m) # 이 예에서만 직교성 검사 목적으로 형성
reflectors = []
for k in range(min(m, n)):
x = R[k:, k].copy()
scale = np.max(np.abs(x))
if scale == 0:
reflectors.append((k, None)); continue
v = x/scale
alpha = -np.copysign(safe_norm(v), v[0])
v[0] -= alpha
v /= safe_norm(v)
R[k:, k:] -= 2*np.outer(v, v@R[k:, k:])
R[k+1:, k] = 0.
Q[:, k:] -= 2*np.outer(Q[:, k:]@v, v)
reflectors.append((k, v.copy()))
return Q, R, reflectors
eps = 1e-8
A = np.array([[1., 1., 1.], [eps, 0., 0.], [0., eps, 0.]])
Q, R, reflectors = householder_qr(A)
assert np.linalg.norm(Q.T@Q-np.eye(3)) < 3e-15
assert np.linalg.norm(A-Q@R)/np.linalg.norm(A) < 3e-15
print("Householder 직교성 및 분해 잔차 확인")
Householder 직교성 및 분해 잔차 확인
5. 최소제곱 우변도 같은 변환에 포함하기#
\(m\ge n\)이고 \(A\)가 완전 열계수이면 \(Q^TA=\begin{pmatrix}R\\0\end{pmatrix}\), \(Q^Tb=\begin{pmatrix}c\\d\end{pmatrix}\)에서
따라서 \(Rx=c\)를 풀고 남은 \(\|d\|\)를 기록하면 됩니다. 이 계산은 \(A\)와 \(b\)에 같은 반사를 적용합니다. 직교성을 크게 잃은 \(\widehat Q\)를 만든 뒤 \(\widehat Q^Tb\)가 정확한 직교좌표인 것처럼 취급하는 것과는 다릅니다.
MGS에서도 확대행렬 \([A\ b]\)의 마지막 열을 같은 순서로 갱신하여 계수와 잔여벡터를 얻을 수 있습니다. Björck–Paige의 관계는 특정 구현의 MGS를 영행을 붙인 확대행렬의 Householder 과정과 연결하여 분석합니다. 임의의 두 라이브러리가 비트 단위로 같은 \(R\)을 내거나, 임의로 만든 \(\widehat Q^Tb\)가 후진안정하다는 주장은 아닙니다. 이 관계의 일반 부동소수점 정리는 조망으로 두고, 이 장에서 사용하는 최소제곱 식은 위 정확한 직교변환으로 증명했습니다. 자세한 안정성과 잔차 조건수는 N4에서 다룹니다.
6. 반사 여러 개를 블록으로 묶기, 두 행만 돌리기#
반사 \(H_i=I-\tau_i v_iv_i^T\)를 각각 적용하면 행렬·벡터곱을 반복합니다. 곱을 \(H_1\cdots H_k=I-VTV^T\)로 묶으면 큰 행렬곱을 사용할 수 있습니다. \(V=[v_1,\ldots,v_k]\)이고 \(T\)는 상삼각입니다. 이 compact WY 표현은 같은 수학적 변환을 다른 연산 순서로 적용하여 메모리 이동과 고속 행렬곱을 활용합니다. 마지막 절에 \(T\)의 재귀식을 유도합니다.
Givens는 두 좌표만 섞습니다. \(r=\operatorname{hypot}(a,b)\), \(c=a/r\), \(s=b/r\)이면
\(a=b=0\)이면 항등변환으로 둡니다. 조밀 정사각 QR에서는 회전이 \(n(n-1)/2\)개 필요하고 각 남은 열에 곱셈 4회와 덧셈 2회가 필요하여 주항 약 \(2n^3\)입니다. Householder의 약 \(4n^3/3\)보다 큽니다. 그러나 상 Hessenberg 행렬은 대각 바로 아래의 \(n-1\)개 성분만 지우면 되므로 \(O(n^2)\)입니다. 이미 계산한 삼각 \(R\)에 관측 행 하나를 붙일 때도 \(n\)번 회전으로 다시 삼각화하여 \(O(n^2)\)에 갱신합니다. 희소성 보존은 회전 자체가 두 행만 건드린다는 뜻이며 임의의 희소행렬에서 fill-in이 없다는 보장은 아닙니다.
7. 열을 바꾸면 수치적 계수가 드러나는가#
열피벗 QR은 현재 남은 열 중 잔여 2-노름이 가장 큰 것을 먼저 고릅니다. \(A\Pi=QR\)이며 \(\Pi\)는 설명변수 순열입니다. 정확한 피벗 선택에서는 \(|r_{11}|\ge|r_{22}|\ge\cdots\)입니다. 작은 후반 대각은 새 설명변수가 앞 변수들의 선형결합에 가깝다는 신호입니다. 하지만 \(|r_{kk}|\)를 \(\sigma_k(A)\)와 동일시하면 안 됩니다. LAPACK의 열피벗 QR 설명도 QR 분해와 수치적 계수 판정을 별도 단계로 다룹니다.
Kahan 계열을 계산해 봅시다. \(0<s<1\), \(c=\sqrt{1-s^2}\)이고 \(K=D U\), \(D=\operatorname{diag}(1,s,\ldots,s^{n-1})\), \(U\)는 대각 1, 대각 위 모든 성분 \(-c\)인 상삼각행렬입니다. \(k\)번째 단계의 각 남은 열의 잔여 제곱노름은 \(s^{2(k-1)}\)로 같습니다. 동률에서 앞 열을 유지하면 열피벗 QR은 순서를 바꾸지 않고 마지막 대각은 \(s^{n-1}\)입니다.
그런데 \(Kx=e_n\)을 뒤에서 풀면 \(x_n=s^{-(n-1)}\), \(x_i=c\sum_{j>i}x_j\)여서 \(x_1=c(1+c)^{n-2}s^{-(n-1)}\)입니다. 따라서
\(|r_{nn}|\)가 최소 특이값보다 지수 인자만큼 클 수 있습니다. 유한정밀도의 동률 처리는 순서를 바꿀 수 있으므로 다음 그림은 정확한 동률 규칙의 QR 대각과 계산한 SVD를 비교합니다. 모든 라이브러리의 기본 피벗 순서가 같다고 주장하지 않습니다.
그림 96 \(s=0.5\)로 고정했다. 세로축은 로그 눈금이고 점선은 증명된 비의 하한 \(c(1+c)^{n-2}\)이다. 이 예는 작은 QR 대각을 유용한 진단으로 쓰는 것과 특잇값을 정확히 대신하는 것을 구별한다.#
완전 공선성도 기준변수 선택 문제를 남깁니다. 절편과 두 집단 더미 \([\mathbf1,d_1,d_2]\)에는 \(\mathbf1=d_1+d_2\)라는 관계가 있습니다. 어느 한 열을 빼도 같은 적합공간을 만들 수 있지만 계수의 의미는 다릅니다. 소프트웨어가 열을 제외할 때에는 열 순서, 단위 스케일, 허용오차를 기록해야 재현할 수 있습니다. 수치적 계수는 선택한 노름과 오차 허용범위에 의존합니다.
8. 직접 비교할 수 있는 직교화 코드#
def gram_schmidt(A, mode="mgs"):
A = np.asarray(A, dtype=float)
m, n = A.shape
Q = np.zeros((m, n)); R = np.zeros((n, n))
for j in range(n):
v = A[:, j].copy()
if mode in ("cgs", "cgs2"):
R[:j, j] = Q[:, :j].T@v
v -= Q[:, :j]@R[:j, j]
if mode == "cgs2":
correction = Q[:, :j].T@v
v -= Q[:, :j]@correction
R[:j, j] += correction
elif mode == "mgs":
for i in range(j):
R[i, j] = Q[:, i]@v
v -= Q[:, i]*R[i, j]
else:
raise ValueError("cgs, cgs2, mgs 중 하나를 선택하세요")
R[j, j] = safe_norm(v)
if R[j, j] == 0:
raise ValueError("선택한 정밀도에서 영 잔여열입니다")
Q[:, j] = v/R[j, j]
return Q, R
for method in ("cgs", "mgs", "cgs2"):
q, r = gram_schmidt(A, method)
orth = np.linalg.norm(q.T@q-np.eye(3), 2)
residual = np.linalg.norm(A-q@r)/np.linalg.norm(A)
print(method, "직교성", orth, "분해 잔차", residual)
assert residual < 1e-14
cgs 직교성 0.7071067811865477 분해 잔차 5.90054702095257e-26
mgs 직교성 1.0000000000311164e-08 분해 잔차 8.003752726706583e-26
cgs2 직교성 3.352983893688406e-16 분해 잔차 5.080387220921977e-25
9. 연습과 전체 풀이#
1. \(x=(3,4)^T\)의 안전한 Householder 벡터와 반사행렬을 구하세요.
풀이. \(\|x\|=5\), \(\alpha=-5\), \(v=(8,4)^T\)이고 \(v^Tv=80\)입니다. 따라서 \(H=I-vv^T/40=\begin{pmatrix}-3/5&-4/5\\-4/5&3/5\end{pmatrix}\)입니다. \(Hx=(-5,0)^T\)이고 두 행의 제곱길이는 1, 내적은 0입니다.
2. \(x=(1,\varepsilon)\)에서 위험한 부호의 첫 반사벡터 성분을 유리화하세요.
풀이. \(1-\sqrt{1+\varepsilon^2}=-\varepsilon^2/(1+\sqrt{1+\varepsilon^2})\)입니다. 정확한 값의 크기는 약 \(\varepsilon^2/2\)인데 직접 뺄셈은 1 근처 수의 저장오차를 이 작은 값과 비교하게 됩니다. 안전한 선택은 \(1+\sqrt{1+\varepsilon^2}\ge2\)입니다.
3. \(G\)가 직교이고 \(\widehat y=Gx+e\)라면 이를 입력오차로 쓰세요.
풀이. \(\widehat y=G(x+G^Te)\)입니다. \(\Delta x=G^Te\)는 \(\|\Delta x\|_2=\|e\|_2\)입니다. 이 등식이 직교변환에서 국소 전진오차를 같은 크기의 후진오차로 바꾸는 근거입니다.
4. \(A=\begin{pmatrix}1&1\\0&0\\0&0\end{pmatrix}\)에서 QR의 대각과 열 선택의 의미를 설명하세요.
풀이. 첫 열 노름은 1이고 둘째 열은 첫 열과 같아 잔여가 0입니다. \(R\)은 첫 행 \((1,1)\), 나머지 0으로 둘 수 있습니다. 두 열 중 어느 쪽을 먼저 골라도 계수는 1입니다. 변수 이름만으로 어느 하나가 더 중요한지 결론낼 수 없습니다.
5. \(R=\begin{pmatrix}3&1\\0&2\end{pmatrix}\)에 행 \((4,0)\)을 추가할 때 첫 Givens를 구하세요.
풀이. 첫 행과 추가 행의 첫 성분 \((3,4)\)에 \(c=3/5\), \(s=4/5\)를 씁니다. 두 행은 \((5,3/5)\)와 \((0,-4/5)\)가 됩니다. 이제 둘째 행의 2와 추가 행의 \(-4/5\)를 회전하면 마지막 행이 0인 상삼각형을 얻습니다. 두 번째 피벗은 \(\sqrt{4+16/25}=2\sqrt{29}/5\)입니다.
6. Kahan 예의 \(n=3\)에서 \(K^{-1}e_3\)를 구하세요.
풀이. 마지막 행에서 \(x_3=s^{-2}\), 둘째에서 \(x_2=c s^{-2}\), 첫째에서 \(x_1=c(x_2+x_3)=c(1+c)s^{-2}\)입니다. 따라서 \(\|K^{-1}\|_2\ge c(1+c)s^{-2}\)이고 \(|r_{33}|/\sigma_{\min}\ge c(1+c)\)입니다. 일반 귀납식의 첫 비자명한 사례입니다.
10. 지금까지의 내용을 수학의 언어로 정리해 봅시다#
이제 사영과 반사를 알고리즘의 식으로 정리합니다. 실수행렬을 사용하며 복소수에서는 전치 대신 수반과 적절한 위상 선택이 필요합니다. 조건수와 무관한 상한에도 차원·연산오차·정상 범위 가정은 남습니다.
정리 1. 정확한 Gram–Schmidt와 재직교화#
완전 열계수 \(A\)에서 CGS와 MGS는 양의 \(R\) 대각을 택하면 같은 축소 QR을 만듭니다. 거의 직교인 \(Q\)에 대해 \(\|Q^TQ-I\|_2\le\delta\), \(P=QQ^T\)라 하고 두 패스가 \(v_1=(I-P)a+e_1\), \(v_2=(I-P)v_1+e_2\)이면
증명. 기존 열들이 직교한다고 가정하면 MGS의 \(i\)번째 계수는 \(q_i^T(a_j-\sum_{k<i}q_kr_{kj})=q_i^Ta_j\)입니다. 따라서 두 절차는 같은 계수와 잔여를 만듭니다. 이 잔여는 기존 열공간과 직교하며 완전 열계수 때문에 0이 아닙니다. 양의 길이로 나누면 같은 \(q_j\)를 얻어 귀납이 끝납니다.
둘째에는 \(E=Q^TQ-I\)라 두면 \(Q^T(I-P)=-EQ^T\)입니다. 그러므로 \(Q^Tv_1=-EQ^Ta+Q^Te_1\), \(Q^Tv_2=E^2Q^Ta-EQ^Te_1+Q^Te_2\)입니다. 유도노름과 삼각부등식을 적용하면 결론입니다. 최종 단위벡터의 기존 열과의 내적 상한은 이를 \(\|v_2\|\)로 나눈 값입니다. 따라서 \(\|v_2\|\)의 하한 없이 무조건 \(O(u)\)라고 말할 수 없습니다. ∎
보조정리 2. 실제 분해로 판정하는 직교성#
가역 \(\widehat R\)에 대해 \(A+F=\widehat Q\widehat R\), \(\widehat R^T\widehat R=A^TA+G\)이면
증명. \(\widehat Q=(A+F)\widehat R^{-1}\)를 대입하고 \(I=\widehat R^{-T}(\widehat R^T\widehat R)\widehat R^{-1}\)를 빼면 가운데 행렬은 \(A^TF+F^TA+F^TF-G\)입니다. 각 곱의 노름을 제한하여 결론을 얻습니다. 작은 분해 잔차 \(F\)만으로는 충분하지 않고, \(G\)와 역삼각인자의 크기도 필요하다는 것을 보입니다. 이 상한 자체는 특정 CGS·MGS의 정밀한 최악오차 정리를 대신하지 않습니다. ∎
정리 3. Householder와 QR#
\(v\ne0\)에 대해 \(H=I-2vv^T/(v^Tv)\)는 대칭·직교이며 \(H^2=I\)입니다. 4절의 안전한 \(v\)는 \(Hx=\alpha e_1\)을 만족합니다. 따라서 모든 실수 \(m\times n\) 행렬은 직교 \(Q\)와 상사다리꼴 \(R\)의 \(A=QR\)을 갖습니다.
증명. \(P=vv^T/(v^Tv)\)라 하면 \(P^T=P\), \(P^2=P\)입니다. \((I-2P)^2=I-4P+4P^2=I\)이고 대칭이므로 직교입니다. \(v=x-\alpha e_1\), \(\alpha^2=x^Tx\)에서 \(v^Tx=\|x\|^2-\alpha x_1\)이고 \(v^Tv=2(\|x\|^2-\alpha x_1)\)입니다. \(x\ne0\)과 안전한 부호에서 이 값은 양수입니다. 따라서 \(2v(v^Tx)/(v^Tv)=v\)이고 \(Hx=x-v=\alpha e_1\)입니다.
첫 열에 적용한 뒤 둘째 열의 둘째 행 이하에 같은 구성을 적용합니다. 이전 열은 그 영역에서 0이므로 바뀌지 않습니다. 영벡터를 만나면 항등변환을 쓰면 됩니다. 최대 \(\min(m,n)\)회 후 상사다리꼴이고, 반사들의 역순 곱의 전치가 직교 \(Q\)를 줍니다. ∎
정리 4. 국소 반사 오차에서 전체 QR 후진오차까지#
각 계산 단계에 정확한 직교행렬 \(H_j\)가 있어 \(\widehat A_j=H_j\widehat A_{j-1}+E_j\), \(\|E_j\|_F\le\eta\|\widehat A_{j-1}\|_F\)이고 \(\widehat A_0=A\)라 합시다. 마지막 결과가 상사다리꼴 \(\widehat R\)이면 정확한 직교 \(Q_*\)가 존재하여
증명. 직교성으로 \(\|\widehat A_j\|_F\le(1+\eta)^j\|A\|_F\)입니다. 재귀식을 전개하면 \(\widehat R=H_n\cdots H_1A+\sum_{j=1}^n H_n\cdots H_{j+1}E_j\)입니다. \(Q_*=H_1^T\cdots H_n^T\)를 곱하면 \(\Delta A\)는 직교행렬을 곱한 각 \(E_j\)의 합입니다. 따라서 노름은 \(\sum_j\eta(1+\eta)^{j-1}\|A\|_F=((1+\eta)^n-1)\|A\|_F\) 이하입니다. ∎
국소 가정은 어떻게 확인하는가. 저장된 비영 반사벡터 \(w\)를 수학적으로 정확히 정규화한 \(z=w/\|w\|\)로 두면 \(H=I-2zz^T\)는 정확히 직교합니다. 구현은 \(B-\tau w(w^TB)\)를 계산하고 \(\tau\)는 \(2/(w^Tw)\)의 근사입니다. 상대오차 \(|\widehat\tau/\tau-1|\le a\)와 N1의 내적·곱·뺄셈 상한을 결합하면, 정상 범위에서 내적 길이 \(m\)에 대해
라는 느슨한 상한을 쓸 수 있습니다. 실제로 내적 오차는 열마다 \(\gamma_m\|w\|\|b\|\) 이하이고, \(\tau\|w\|^2=2\)이므로 이 오차를 바깥 곱으로 옮긴 크기는 \(2(1+a)\gamma_m\|b\|\) 이하입니다. 스칼라·성분 곱의 두 반올림과 마지막 뺄셈에는 \(\gamma_3(\|b\|+2(1+a)(1+\gamma_m)\|b\|)\)를 더하면 됩니다. \(\gamma_{m+3}\le1/10\)에서 이 합은 위 상한 이하입니다. \(\tau\) 자체를 내적과 나눗셈으로 구하면 곱 보조정리에 의해 \(a\le\gamma_{m+2}\)로 둘 수 있습니다.
반사벡터 생성 오차도 확인해야 합니다. 정확한 \(x\)를 최대 절댓값으로 스케일한 후 노름을 계산하고 안전한 부호를 쓰면, 제곱·합·제곱근·덧셈·나눗셈의 상대오차를 N1로 묶을 수 있습니다. \(g=\gamma_{m+5}\le1/100\)이라 놓으면 계산한 스케일 벡터와 정확한 스케일 벡터의 차이는 상대 \(u\), 그 노름의 차이는 상대 \(3g\) 이하입니다. 여기서 제곱합의 상대오차는 \(\gamma_{m+2}\) 이하이고 \(|\sqrt{1+t}-1|\le|t|\) (\(|t|\le1/2\))를 썼습니다. 안전한 \(v=x-\alpha e_1\)에는 \(\|v\|\ge\|x\|\)이므로 성분 덧셈까지 포함한 상대 벡터오차는 \(5g\) 이하입니다. 정규화 부등식
은 두 항으로 분리하고 역삼각부등식을 쓰면 얻습니다. 따라서 저장벡터를 정확히 정규화한 방향은 이상적 방향에서 \(10g\) 이하이고, \(\|zz^T-tt^T\|_2\le2\|z-t\|\)에서 정확한 반사 둘의 차이는 \(40g\) 이하입니다. 이 때문에 계산 뒤 첫 열의 작은 아래 성분을 0으로 저장하는 조작도 \(O(g)\|B\|_F\)의 추가 오차입니다. 앞 적용 오차를 두 번 더하는 여유를 포함하여 \(\eta=100g\)를 쓰면 국소 단계 모형을 만족합니다. 이는 상수를 최적화한 결과가 아니라 조건수에 의존하지 않는 명시적 상한입니다. 실제 저장한 \(\widehat Q\)와 정확한 \(Q_*\)는 구별하며, 전자의 직교성도 반사를 형성하는 별도 반올림의 영향을 받습니다.
정리 5. compact WY와 열피벗 대각의 한쪽 상한#
\(H_1\cdots H_k=I-V_kT_kV_k^T\)일 때 새 반사를 붙인 표현은
정확한 열피벗 QR에서는 \(|r_{kk}|\ge\sigma_k(A)/\sqrt{n-k+1}\)입니다.
증명. \((I-VTV^T)(I-\tau vv^T)\)를 전개하면 \(I-VTV^T-\tau vv^T+\tau VTV^Tvv^T\)입니다. 제시한 블록 곱을 전개하면 같은 식이므로 귀납적으로 WY 표현이 성립합니다.
QR의 첫 \(k-1\)행만 남긴 근사는 계수가 최대 \(k-1\)이고 잔차는 뒤 블록 \(R_{22}\)입니다. 피벗 선택에서 이 블록의 모든 열노름은 첫 열의 노름 \(|r_{kk}|\) 이하입니다. 따라서 \(\|R_{22}\|_2\le\|R_{22}\|_F\le\sqrt{n-k+1}|r_{kk}|\)입니다. H5의 최적 저계수 근사 정리에서 어떤 계수 \(k-1\) 근사도 스펙트럼노름 오차가 \(\sigma_k(A)\) 이상이므로 결론입니다. 이는 \(r_{kk}\)가 특이값보다 지나치게 작아지는 것을 막는 한쪽 상한이며, Kahan 예처럼 지나치게 큰 경우를 막지 않습니다. ∎
별표: QR로 무작위 직교행렬 만들기#
확률을 배운 뒤 읽을 절입니다. 독립 표준정규 성분의 정사각행렬 \(G\)는 확률 1로 가역입니다. 특이행렬은 비영 다항식 \(\det G\)의 영점집합이고 이 집합의 Lebesgue 측도가 0이라는 사실을 외부 전제로 사용합니다. \(G=QR\)에서 \(R\) 대각을 양수로 맞추면 \(Q\)의 분포는 모든 고정 직교행렬 \(U\)에 대한 왼쪽 곱에 불변입니다.
정규밀도는 \(\exp(-\|G\|_F^2/2)\)에 비례하고 \(G\mapsto UG\)는 길이와 부피를 보존하므로 \(UG\)와 \(G\)의 분포가 같습니다. 양의 대각 QR의 유일성으로 \(UG=(UQ)R\)의 직교인자는 \(UQ\)여서 불변성이 따라옵니다. 이러한 정규화된 불변분포를 Haar 분포라 부릅니다. 일반 QR 루틴이 임의로 고른 대각 부호는 이 정규화를 대신하지 않으므로, \(d_i=\operatorname{sign}(r_{ii})\)로 \(Q\)의 열과 \(R\)의 행에 같은 부호를 곱해야 합니다.
직교화는 정확한 산술의 독립성과 컴퓨터에서의 구별 가능성이 다르다는 점을 보여 줍니다. 다음 장에서는 QR과 SVD를 최소제곱에 적용하고 작은 특이방향을 어떻게 처리할지 살펴봅니다. N4로 이어 읽기.