열 피벗 QR

1. 개요

열 피벗 QR QR with Column Pivoting
형태	AΠ = QR (Π: 열 순열)
피벗 규칙	남은 열 중 (갱신된) 노름이 가장 큰 열을 앞으로
보장	|r11| ≥ |r22| ≥ … (대각 비증가)
제안	비징거·골럽, 1965
LAPACK	xGEQP3 (BLAS-3판), 최소제곱 xGELSY
약점	카한 행렬처럼 랭크 결핍을 못 드러내는 반례가 있음

SVD 는 비싸고, 그냥 QR 은 랭크를 모른다. 그 사이에서 “대개는 맞는” 값싼 타협.

열 피벗 QR(QR factorization with column pivoting, QRCP)은 하우스홀더 QR 분해를 하면서 매 단계 남은 열 중 노름이 가장 큰 열을 앞으로 끌어오는 열 교환을 함께 수행해 $A\Pi = QR$ 를 얻는 분해다. 그 결과 $R$ 의 대각이 크기순으로 내려가서, 대각이 어디서 뚝 떨어지는지 보면 행렬의 수치 랭크를 대충 읽을 수 있다. 비징거와 골럽이 1965년 랭크가 모자란 최소제곱 문제를 풀려고 제안했다.

QR 분해 본문이 이것을 한 절로 요약했다 — “SVD 보다 훨씬 싸고 대개는 맞다”. 이 문서는 그 “대개”의 내용을 연다. 피벗 규칙이 왜 대각을 정렬하는지, 열 노름을 싸게 갱신하는 요령과 그 요령이 실제로 LAPACK 을 틀리게 만든 사건, 카한 반례가 피벗을 어떻게 속이는지, 그리고 이 분해가 지금 의외로 많이 쓰이는 곳 — 열 고르기·센서 배치·축소차수모델의 보간점 선택 — 이다.

2. 알고리즘

$m \times n$ 행렬 $A$ 에 대해 $k = 1, 2, \dots$ 단계마다:

1. 아직 처리 안 된 열 $j \ge k$ 각각의 남은 부분 $A^{(k)}_{k:m,\,j}$ 의 노름 $\nu_j$ 를 본다.
2. 가장 큰 $\nu_j$ 를 가진 열을 $k$ 번째 열과 교환하고, 그 교환을 $\Pi$ 에 기록한다.
3. 하우스홀더 반사로 $k$ 열의 $k+1$ 행 아래를 0 으로 만든다. 이제 $|r_{kk}| = \max_j \nu_j$ 다.
4. 남은 열들의 노름을 갱신한다.

3 에서 바로 보이듯 $|r_{kk}|$ 는 그 시점 남은 열 노름의 최댓값이고, 다음 단계의 남은 열 노름은 하우스홀더가 한 행을 떼어 간 만큼 줄기만 하므로

$$|r_{11}| \ge |r_{22}| \ge \cdots \ge |r_{nn}|, \qquad |r_{kk}|^2 \ge \sum_{i=k}^{j} |r_{ij}|^2 \quad (j > k)$$

가 성립한다. 둘째 부등식이 ”$R$ 의 각 열에서 대각이 그 아래 남은 부분을 지배한다”는 뜻이다. 랭크가 $r$ 인 행렬이라면 정확한 산술에서는 $r$ 단계 뒤 남은 열 노름이 모두 0 이 되고 분해가 멈춘다.

2.1. 열 노름 갱신 — 싼 요령과 그 대가

4 를 매번 처음부터 계산하면 단계당 $O(mn)$ 이 더 든다. 대신 직교 변환이 노름을 보존한다는 걸 이용해

$$\nu_j^{\text{new}} = \sqrt{\nu_j^2 - r_{kj}^2} = \nu_j\sqrt{1 - (r_{kj}/\nu_j)^2}$$

로 $O(1)$ 에 갱신한다. 문제는 $\nu_j$ 가 $r_{kj}$ 와 거의 같을 때 생기는 상쇄다. 몇 단계 이어지면 갱신된 노름에 상대 정확도가 하나도 남지 않고, 엉뚱한 열이 피벗으로 뽑힌다. 그래서 갱신값이 원래 값보다 충분히 작아지면($\sqrt{\varepsilon}$ 기준 근처) 그 열 노름을 다시 직접 계산하는 안전장치가 들어간다.

이 안전장치의 문턱이 원래 LAPACK 에서 부정확했고, 드르마치와 부야노비치가 피벗이 실제로 틀어져 랭크 판정이 어긋나는 행렬들을 보고했다(LAPACK Working Note 176). 수정판이 LAPACK 3.1 의 xGEQP3·xLAQPS 에 들어갔다. “휴리스틱이 틀린다”와 “휴리스틱을 구현한 코드가 틀린다”는 다른 사고인데, 랭크 판정에서는 둘 다 실제로 일어났다는 점이 이 분해의 교훈이다.¹

xGEQP3 은 퀸타나-오르티·선·비숍(1998)의 BLAS-3 판으로, 이전 xGEQPF 를 대체했다. 다만 피벗을 고르려면 매 단계 노름이 갱신돼야 해서 연산의 절반가량은 여전히 행렬-벡터 연산으로 남는다. 같은 크기에서 피벗 없는 xGEQRF 보다 눈에 띄게 느린 이유다.

3. 랭크 드러내기와 그 한계

대각이 정렬되면 ”$|r_{kk}|$ 가 허용오차보다 작아지는 첫 $k$” 를 수치 랭크로 잡고 싶어진다. 실제로 대부분 잘 된다.

검산 — 진짜 랭크 20. $200\times100$ 의 무작위 랭크-20 행렬에 크기 $10^{-10}$ 의 잡음을 더해 열 피벗 QR 을 했다.

	19번째	20번째	21번째	22번째
$\lvert r_{kk}\rvert$	34.9	23.9	$3.2\times10^{-9}$	$2.8\times10^{-9}$
$\sigma_k$	82.0	69.0	$2.2\times10^{-9}$	$2.1\times10^{-9}$

20 과 21 사이의 10 자릿수 절벽이 그대로 보인다. 개별 값은 특이값과 다르지만 간격은 정확히 같은 곳에 있다.

3.1. 카한 행렬 — 피벗이 아무것도 안 하는 경우

카한 행렬은 $s^2 + c^2 = 1$ 에 대해

$$K_n = \mathrm{diag}(1, s, s^2, \dots, s^{n-1}) \begin{pmatrix} 1 & -c & -c & \cdots & -c\\ & 1 & -c & \cdots & -c\\ & & \ddots & & \vdots\\ & & & & 1 \end{pmatrix}$$

이다. 이미 상삼각이라 $R = K_n$ 이고, 계산해 보면 모든 열의 노름이 정확히 1 이다. 열 피벗은 동점에서 순서를 바꿀 이유가 없다(실무에선 동점이 반올림으로 깨지므로 열 $j$ 에 $1 - 10^{-12} j$ 정도를 곱해 원래 순서를 고정한다). 그러면 $R$ 의 대각은 $s^{k-1}$ 로 완만하게 줄어 아무 절벽도 없다.

그런데 $n = 90$, $c = 0.285$ 에서 실측하면:

• 열 피벗 QR 의 순열 = 항등(아무 교환도 없음)
• 가장 작은 대각 $|r_{nn}| = s^{89} \approx 2.3\times10^{-2}$
• 실제 최소 특이값 $\sigma_{90} \approx 8.8\times10^{-12}$, 그 위 $\sigma_{89} \approx 2.7\times10^{-2}$

행렬은 수치적으로 랭크 89 인데(간격 9 자릿수), $R$ 의 대각은 그걸 $10^9$ 배 과대평가해서 풀랭크라고 말한다. 일반적으로 열 피벗 QR 이 보장하는 것은 $|r_{kk}|$ 와 $\sigma_k$ 의 비가 최악 $2^k$ 꼴의 지수적 한계 안에 있다는 것뿐이고, 카한 행렬이 그 지수를 실현한다.

이 간극을 메우려는 것이 강한 RRQR(rank-revealing QR, 구·아이젠슈타트 1996)이다. 기본 열 피벗 뒤에 열을 몇 번 더 맞바꿔 $\sigma_k(R_{11})$ 과 $\sigma_{k+1}$ 이 다항식 인자 안에서 일치하도록 보장한다. 이론적으로 확실하지만 실무 라이브러리의 기본값은 여전히 보통 열 피벗이다 — 카한 같은 행렬을 실전에서 만날 확률이 낮고, 진짜 확실해야 하면 특이값 분해로 가면 되기 때문이다.

4. 랭크 결핍 최소제곱 — 기본 해와 최소 노름 해

$A$ 의 수치 랭크가 $r < n$ 이면 $A\Pi = Q\begin{pmatrix}R_{11} & R_{12}\\ 0 & R_{22}\end{pmatrix}$ 에서 $R_{22} \approx 0$ 으로 버리고 푼다.

• 기본 해(basic solution): $R_{11} y = (Q^{\mathsf T}b)_{1:r}$ 를 풀고 나머지 $n - r$ 성분은 0. 피벗이 고른 $r$ 개 열만으로 $b$ 를 설명하는 해다. MATLAB 의 직사각 A\b 가 랭크 결핍 경고와 함께 내는 해가 이것이다.
• 최소 노름 해: 기본 해는 일반적으로 노름 최소가 아니다. $\begin{pmatrix}R_{11} & R_{12}\end{pmatrix}$ 를 오른쪽에서 직교 변환으로 한 번 더 삼각화하면(완전 직교 분해, $A = Q\,T\,Z^{\mathsf T}$) 무어-펜로즈 유사역행렬 해가 나온다. LAPACK xGELSY 가 이 경로다. 반면 NumPy lstsq 는 SVD 기반 xGELSD 를 쓴다.

어느 쪽을 원하는지 모르고 섞어 쓰면, 같은 데이터에 두 라이브러리가 둘 다 맞는 다른 답을 내서 디버깅 지옥이 열린다.²

5. 의외의 주력 — 열 고르기

열 피벗 QR 의 순열 $\Pi$ 자체가 “가장 서로 독립적인 열 $k$ 개”에 대한 탐욕적 답이다. 분해보다 이 순열이 목적인 응용이 늘었다.

• 보간 분해·CUR. $A \approx A_{:,J}\,X$ 처럼 원래 열 몇 개로 행렬을 근사한다. 열이 실제 데이터(실제 센서, 실제 문서, 실제 유전자)라 해석이 쉽다. 랜덤 스케치 $\Omega A$ 에 열 피벗 QR 을 걸어 큰 행렬에서도 싸게 고른다.
• Q-DEIM. 축소차수모델의 비선형항을 몇 개의 격자점에서만 평가하려면 그 점을 골라야 한다. 적합직교분해 기저 $U$ 의 전치 $U^{\mathsf T}$ 에 열 피벗 QR 을 걸어 처음 $k$ 개 피벗을 보간점으로 쓰는 것이 드르마치·구게르친(2016)의 Q-DEIM 이다. 원래 DEIM 의 탐욕 절차보다 오차 상수가 좋다.
• 센서 배치. 같은 원리로 POD 모드에서 “측정하면 모드 계수를 가장 잘 복원하는 위치”를 고른다. 유동장·기후장 재구성에서 쓰인다.
• 랜덤 QRCP. 매 단계 노름 비교가 통신 병목이라, 작은 랜덤 스케치에서 피벗 블록을 미리 고르고 원 행렬엔 BLAS-3 로 한 번에 적용하는 방식이 대형 병렬 환경에서 쓰인다.³

6. 관련 문서

• QR 분해 · 특이값 분해 · 절단 특이값 분해 · 무어-펜로즈 유사역행렬
• 최소자승법 · LAPACK · 조건수 · 가우스 소거법
• 축소차수모델 · 적합직교분해
• RRQR · 보간 분해 · CUR 분해

7. Footnotes

1. 이 사건 이후 “LAPACK 이 틀릴 리 없다”는 믿음은 “LAPACK 도 사람이 짰다”로 살짝 수정됐다. 그래도 고쳐진 속도와 투명성을 보면 여전히 이 바닥에서 가장 믿을 만한 코드인 건 맞다. ↩

2. 대표 증상: MATLAB 과 NumPy 로 같은 회귀를 돌렸는데 계수가 다르고 잔차는 똑같다. 둘 다 정답이다. 이걸 이해하기 전까지는 둘 중 하나가 버그라고 확신하며 반나절을 날린다. ↩

3. 피벗 하나 고르려고 수천 노드가 전원 통신하는 건 낭비라는, 지극히 현생적인 동기다. 수학적으로 더 우아해서가 아니라 네트워크가 느려서 생긴 알고리즘이 HPC 에는 의외로 많다. ↩