행렬 제곱근

1. 개요[편집]

행렬 제곱근 Matrix square root
문제	$X^{2} = A$ 를 만족하는 $X$
해의 개수	정칙 · 서로 다른 고유값 $s$ 개면 1차 제곱근 $2^{s}$ 개
주 제곱근	고유값이 $\mathbb{R}^{-}$ 에 없으면 유일 (스펙트럼이 우반평면)
표준 알고리즘	슈어 방법 (Björck–Hammarling, 1983)
반복법	DB (Denman–Beavers) · 뉴턴–슐츠
소프트웨어	MATLAB sqrtm, SciPy scipy.linalg.sqrtm

$\sqrt{4} = \pm2$ 라는 초등학교 수준의 다가성이, 행렬로 올라가면 연속체만큼 많은 해와 “아예 해가 없는 행렬”을 낳는다.

행렬 제곱근(matrix square root)은 주어진 정사각행렬 $A$ 에 대해 $X^{2} = A$ 를 만족하는 행렬 $X$ 다. 스칼라와 달리 해의 개수부터가 사연이 있고, 그중 어느 것을 답이라고 부를지를 정하는 것이 문제의 절반이다. 실무에서 "$A^{1/2}$" 라고 쓰면 거의 언제나 주 제곱근(principal square root), 즉 스펙트럼이 열린 우반평면에 놓인 유일한 해를 뜻한다.

행렬함수의 일반 이론이 그대로 적용되지만, 제곱근은 유독 알고리즘이 갈라지는 함수다. 우아한 뉴턴 반복이 수치적으로 못 쓸 물건이고, 대신 슈어 분해 위에서 삼각행렬을 손으로 채우는 재귀가 표준이 되었다. 그 사연이 이 문서의 중심이다.

2. 해가 몇 개나 있나[편집]

$A$ 가 정칙이고 서로 다른 고유값이 $s$ 개라면, 각 고유값에서 $\pm\sqrt{\lambda_i}$ 중 하나를 고르는 만큼 — 즉 정확히 $2^{s}$ 개의 1차(primary) 제곱근이 있다. 1차라는 것은 “같은 고유값에는 반드시 같은 가지를 배정한다”는 뜻이고, 이 조건이 붙어야 $X$ 가 $A$ 의 다항식으로 표현되고 $A$ 와 교환한다.

가지를 고유값마다 제각각 배정하면 비1차 제곱근이 나오는데, 이건 고유값이 겹치면서 비유도(derogatory), 즉 같은 고유값이 조르당 블록 두 개 이상에 나타날 때만 가능하다. $A = I_2$ 가 교과서 예다. 1차 제곱근은 $\pm I$ 둘뿐인데, 임의의 정칙 $Z$ 에 대해

$$Z\begin{bmatrix}1 & \\ & -1\end{bmatrix}Z^{-1}$$

도 전부 제곱하면 $I$ 다. 연속체만큼 많다. 이들은 $I$ 와 교환하기는 하지만(단위행렬이니까) 일반의 비1차 제곱근은 $A$ 와 교환하지 않고, 어떤 알고리즘도 이들을 “구하려” 하지 않는다.

반대로 제곱근이 아예 없는 행렬도 있다. 고전적인 예가

$$N = \begin{bmatrix}0 & 1 \\ 0 & 0\end{bmatrix}$$

이다. $X^2 = N$ 이면 $X^4 = 0$ 이라 $X$ 는 멱영이고 $2\times2$ 멱영행렬은 $X^2 = 0 \ne N$ 이므로 모순. 일반적으로 특이행렬의 제곱근 존재 여부는 영고유값에 붙은 조르당 블록 크기들이 짝지어지는지로 결정된다. 정칙이면 언제나 존재한다.¹

3. 주 제곱근[편집]

$A$ 의 고유값이 닫힌 음의 실축 $(-\infty, 0]$ 에 하나도 없으면, 모든 고유값이 우반평면에 오도록 가지를 고른 제곱근이 유일하게 존재한다. 이것이 주 제곱근 $A^{1/2}$ 다.

• $A$ 가 실행렬이면 $A^{1/2}$ 도 실행렬이다.
• $A$ 가 에르미트 양정부호면 $A^{1/2}$ 도 에르미트 양정부호이고, 촐레스키 분해와 달리 유일하며 대칭이다. 이 성질 때문에 통계·양자화학이 촐레스키 대신 제곱근을 고집하는 경우가 생긴다.
• 조건수는 얌전한 편이다. 에르미트 양정부호에서 상대 조건수가 $\tfrac12\,\kappa_2(A)^{1/2}$ 로, 원 문제 조건수의 제곱근 규모다. 제곱근을 취하는 연산이 조건수도 제곱근으로 눌러 준다.

실 제곱근은 따로 따져야 한다. 실행렬 $A$ 가 실 제곱근을 가지려면 모든 음의 실고유값에 대해 각 크기의 조르당 블록이 짝수 개씩 있어야 한다. 그래서 $-1$($1\times1$)은 실 제곱근이 없고, $-I_2$ 는 회전행렬

$$\begin{bmatrix}0 & 1\\ -1 & 0\end{bmatrix}$$

이라는 실 제곱근을 갖는다. 짝이 맞으면 복소 켤레쌍처럼 묶어 실수 산술로 처리할 수 있다는 뜻이고, 짝이 안 맞으면 답이 복소수로 나온다. 실 데이터를 넣었는데 sqrtm 이 복소 배열을 뱉는 사고의 대부분이 여기서 온다.

4. 슈어 방법 — 실무 표준[편집]

뵈르크–함마를링(1983)의 방법은 단순하고 강력하다. $A = QTQ^{H}$ 로 슈어 분해하면 $U^2 = T$ 인 상삼각 $U$ 를 찾는 문제로 바뀌고, $A^{1/2} = QUQ^{H}$ 다. 상삼각행렬의 제곱은 다시 상삼각이므로 구조가 보존된다.

$(U^{2})_{ij} = t_{ij}$ 를 성분별로 쓰면 대각은 즉시 $u_{ii} = \sqrt{t_{ii}}$ 이고, 비대각은

$$u_{ij} \;=\; \frac{t_{ij} - \displaystyle\sum_{k=i+1}^{j-1} u_{ik}u_{kj}}{u_{ii} + u_{jj}}$$

로 대각선에서 바깥으로 한 겹씩 채워 나간다. 삼각 부분의 비용은 $n^3/3$ flops 뿐이고 대부분의 시간은 슈어 분해가 먹는다(전체 $\approx 28n^3$ 규모). 하이엄(1987)이 켤레쌍을 $2\times2$ 블록으로 다루는 실 슈어 판을 만들어, 실행렬은 복소 산술 없이 처리한다.

분모 $u_{ii} + u_{jj}$ 를 보자. 주 가지를 고르면 $\operatorname{Re}u_{ii} \ge 0$ 이므로 이 합이 0이 되려면 두 값이 모두 순허수이고 부호가 반대여야 하는데, 그건 $t_{ii}, t_{jj}$ 가 음의 실축 위에 있다는 뜻이다. 주 제곱근의 존재 조건이 곧 이 점화식이 나눗셈 사고를 내지 않을 조건이라는 게 깔끔하다. 다만 분모가 0은 아니어도 작으면 조건이 나빠지므로, 가까운 고유값을 한 블록에 모아 블록 단위로 실베스터 방정식을 푸는 슈어-파레 방법 스타일의 블록화 판을 쓴다.

안정성에는 각주 하나가 붙는다. 계산된 $\hat X$ 는 $\hat X^{2} = A + \Delta A$, $\lVert\Delta A\rVert \lesssim c\,u\,\lVert\hat X\rVert^{2}$ 를 만족한다 — 잔차가 $\lVert A\rVert$ 가 아니라 $\lVert X\rVert^{2}$ 에 비례한다. 두 값의 비 $\alpha = \lVert X\rVert^{2}/\lVert A\rVert$ 는 항상 1 이상이고, 비정규 행렬에서는 꽤 커질 수 있다. “안정적”이라는 말이 언제나 $\lVert A\rVert$ 기준이라는 보장은 없다는 것을 상기시키는 사례.

5. 뉴턴 반복은 왜 못 쓰나[편집]

$X^2 = A$ 에 뉴턴-랩슨법을 걸면 그림처럼 예쁜 식이 나온다.

$$X_{k+1} = \tfrac12\left(X_k + X_k^{-1}A\right), \qquad X_0 = A$$

수학적으로는 2차 수렴하고 $A^{1/2}$ 로 간다. 그런데 수치적으로 불안정하다. 고정점 $X = A^{1/2}$ 근처에서 오차 전파를 선형화하면

$$E_{k+1} \approx \tfrac12\left(E_k - X^{-1}E_kX\right)$$

이고, $X$ 의 고유값 $\mu_i = \lambda_i^{1/2}$ 로 대각화하면 $(i,j)$ 성분이 매 반복 $\tfrac12\left(1 - \mu_j/\mu_i\right)$ 배가 된다. 반복이 안정하려면 모든 $i,j$ 에 대해

$$\left\lvert \frac{\lambda_j^{1/2}}{\lambda_i^{1/2}} - 1 \right\rvert \;\le\; 2$$

여야 한다. 고유값 제곱근의 비가 대략 3배를 넘으면 그 순간 반올림 오차가 매 반복 증폭된다는 뜻이고, 실전 행렬에서 이 조건이 지켜지는 경우는 흔치 않다. 겉보기엔 잘 수렴하다가 어느 순간 발산하거나 엉뚱한 값에 정착하는 전형적인 증상이 나온다. 행렬 부호 함수의 뉴턴 반복이 항상 안정한 것과 정확히 반대라, 두 반복이 형제처럼 생겼다는 사실이 오히려 함정이 된다.²

6. DB 반복과 뉴턴–슐츠[편집]

원인은 반복이 $X_k^{-1}A$ 라는 비대칭 결합을 만든다는 데 있고, 처방은 $A^{1/2}$ 와 $A^{-1/2}$ 를 동시에 굴려 결합을 끊는 것이다. 덴먼–비버스(1976) 반복:

$$Y_0 = A,\quad Z_0 = I,\qquad Y_{k+1} = \tfrac12\left(Y_k + Z_k^{-1}\right),\quad Z_{k+1} = \tfrac12\left(Z_k + Y_k^{-1}\right)$$

$Y_k \to A^{1/2}$, $Z_k \to A^{-1/2}$ 로 2차 수렴하며 안정하다. 반복당 역행렬이 두 번이라 슈어 방법보다 비싸지만, 연산이 전부 역행렬과 곱이라 병렬·GPU 이식이 쉽다. 수렴 가속에는 부호 함수와 같은 행렬식 스케일링 $\mu_k = \lvert\det Y_k \det Z_k\rvert^{-1/(2n)}$ 을 쓴다.

DB 가 어디서 왔는지는 알고 보면 허탈하다. 항등식

$$\operatorname{sign}\begin{bmatrix}0 & A \\ I & 0\end{bmatrix} = \begin{bmatrix}0 & A^{1/2} \\ A^{-1/2} & 0\end{bmatrix}$$

의 좌변에 행렬 부호 함수의 뉴턴 반복을 적용하고 블록으로 풀어 쓰면 그대로 DB 다. 제곱근을 부호 함수 문제로 승격시켰더니 안정해졌다는 것이고, 반복의 안정성이 문제의 표현 방식에 달려 있음을 보여 주는 교과서적 사례다. 곱의 형태로 정리해 반올림에 더 강하게 만든 product-form DB, 반복마다 증분만 갱신하는 incremental Newton 등의 변종이 있다.

역행렬조차 피하고 싶으면 뉴턴–슐츠의 결합 판을 쓴다.

$$Y_{k+1} = \tfrac12 Y_k\left(3I - Z_kY_k\right), \qquad Z_{k+1} = \tfrac12\left(3I - Z_kY_k\right)Z_k$$

행렬곱만으로 이루어져 있어 텐서 코어에 그대로 얹히지만, $\lVert I - A\rVert < 1$ 일 때만 수렴하므로 사전 스케일링이 필수다. 딥러닝 쪽에서 공분산 화이트닝 층이나 $A^{-1/2}$ 전처리를 GPU 에서 미분 가능하게 구현할 때 이 반복이 재발견되어 널리 쓰인다.

7. 어디에 쓰나[편집]

• 화이트닝과 직교화. 공분산 $\Sigma$ 에 대해 $\Sigma^{-1/2}x$ 로 변수를 표준화하는 것이 화이트닝이다. 촐레스키 인수를 써도 되지만 대칭 제곱근을 쓰면 원 좌표축과의 정렬이 최대로 유지된다(ZCA). 양자화학의 뢰딘 대칭 직교화 $S^{-1/2}$ 도 같은 계산이며, 겹침행렬을 “가장 덜 비틀어” 정규직교기저로 만든다는 변분적 의미가 붙는다.
• 표본 생성과 제곱근 필터. $\Sigma^{1/2}\xi$ 로 다변량 정규 표본을 만들고, 앙상블 칼만 필터의 ETKF 류는 공분산을 인수 형태로만 들고 다니며 제곱근을 갱신한다. 가우시안 프로세스 표본 생성도 같은 자리.
• 로그의 역스케일링-제곱. $\log A = 2^{s}\log\!\left(A^{1/2^{s}}\right)$ 를 이용해 제곱근을 반복해서 $A$ 를 $I$ 쪽으로 밀어 놓고 파데 근사를 쓴다. 행렬 로그 알고리즘의 내부 루프가 통째로 제곱근이라, logm 의 정확도는 sqrtm 의 정확도에 종속된다.
• 극분해와 기하평균. $A = UH$ 의 에르미트 인수는 $H = (A^{H}A)^{1/2}$ 다. 두 양정부호 행렬의 기하평균 $A \,\#\, B = A^{1/2}(A^{-1/2}BA^{-1/2})^{1/2}A^{1/2}$ 도 제곱근 세 번이며, 확산텐서 영상이나 공분산 보간에서 쓰인다.
• 분수 거듭제곱. $A^{p}$ 를 $p$ 의 이진 전개로 쪼개면 반복 제곱근이 필요하고, 마르코프 연쇄의 “반년 전이행렬”($P^{1/2}$) 같은 통계 문제가 여기 걸린다. 참고로 확률행렬의 제곱근이 다시 확률행렬이라는 보장은 없어서, 음수 성분이 나오는 유명한 골칫거리가 있다.
• 전처리. 대칭 양정부호 전처리기 $M$ 을 $M^{1/2}$ 로 쪼개 좌우 대칭으로 거는 것은 전처리기 붙은 켤레기울기법의 이론 전개에 늘 등장하지만, 실제로 $M^{1/2}$ 를 만드는 구현은 없다. 대칭 분할이 존재한다는 사실만 쓰고 계산은 하지 않는, 순수하게 증명용으로 소비되는 행렬 제곱근이다.³

8. 관련 문서[편집]

• 행렬함수 · 행렬 부호 함수 · 행렬 지수함수 · 슈어-파레 방법
• 슈어 분해 · 조르당 표준형 · 비정규 행렬 · 극분해
• 촐레스키 분해 · 특이값 분해 · 실베스터 방정식 · 파데 근사
• 화이트닝 · 주성분 분석 · 앙상블 칼만 필터 · 가우시안 프로세스
• 조건수 · 후진 오차 해석 · 뉴턴-랩슨법 · 전처리기

9. Footnotes[편집]

1. 스칼라에서는 $0$ 의 제곱근이 $0$ 하나로 아무 문제가 없는데, 행렬에서는 “영고유값이 어떻게 생겼는지”가 존재성을 좌우한다. 크기 2짜리 조르당 블록 $J_2(0)$ 하나가 외톨이로 있으면 제곱근이 없고, 두 개면 있다. 대수학의 사소해 보이는 구조가 존재성을 통째로 뒤집는 지점이라, 처음 배우는 사람이 “그럴 리가”라고 말하는 순간 예제를 손으로 확인시키는 것이 이 바닥의 국룰. ↩

2. $\tfrac12(X + X^{-1})$ 은 안정하고 $\tfrac12(X + X^{-1}A)$ 는 불안정하다. 두 식의 차이는 $A$ 하나뿐인데 운명이 갈린다. 이유는 앞쪽 반복의 오차 사상 $\tfrac12(E - SES)$ 가 멱등 사영자인 반면 뒤쪽은 고유값 비에 따라 증폭 인자가 1을 넘길 수 있다는 것. “수식이 예쁘면 안정할 것”이라는 직관이 배신당하는 대표 사례다. ↩

3. 켤레기울기법 수렴 증명에는 $M^{-1/2}AM^{-1/2}$ 의 조건수가 등장하지만, 코드에는 $M^{-1}r$ 을 푸는 줄만 있다. 두 세계를 잇는 것은 “대칭 분할이 존재하므로 스펙트럼이 같다”는 한 문장이고, 그 한 문장 덕분에 아무도 만들지 않는 행렬이 모든 교과서에 등장한다. 계산하지 않는 대상을 정의해 두는 것이 왜 유용한가에 대한 좋은 답변. ↩