자기회귀 모형

1. 개요

자기회귀 모형 Autoregressive Model, AR(p)
식	xt = c + φ1xt−1 + … + φpxt−p + εt
정상성	특성다항식 1 − φ1z − … − φpzp 의 근이 모두 단위원 밖
식별 신호	자기상관은 서서히 감쇠, 부분 자기상관은 p 차 이후 0
추정	율-워커(레빈슨-더빈), 최소제곱, 최대우도, 버그
차수 선택	AIC, BIC
기원	율(1927) 태양흑점 AR(2), 워커(1931)

내일의 나는 어제의 나와 그제의 나의 가중평균 + 약간의 운. 시계열 모형계의 “과거는 반복된다”.

자기회귀 모형(autoregressive model, AR)은 시계열의 현재 값을 자기 자신의 과거 $p$ 개 값의 선형결합과 백색 잡음의 합으로 쓰는 확률 모형이다.

$$x_t = c + \sum_{i=1}^{p} \varphi_i\, x_{t-i} + \varepsilon_t, \qquad \varepsilon_t \sim \text{i.i.d. } (0, \sigma^2)$$

이것을 AR($p$) 라 부른다. 1927년 율(G. U. Yule)이 볼퍼 태양흑점 수를 「주기 함수 + 측정 잡음」이 아니라 「잡음이 계속 흔들어 대는 감쇠 진자」로 설명하려고 AR(2) 를 맞춘 것이 시작이고, 워커(1931)가 일반 차수로 정리했다. 계수와 자기공분산을 잇는 연립방정식이 두 사람 이름을 딴 율-워커 방정식이다.

시뮬레이션 쪽에서 AR 모형은 세 군데서 자주 나온다. 첫째, 합성 입력 생성 — 풍하중·지진동·난류 유입 경계조건처럼 “이런 스펙트럼을 가진 무작위 시계열”이 필요할 때. 둘째, 출력 분석 — 몬테카를로·분자동역학 출력의 자기상관을 모형화해 통계 오차를 바로 잡을 때. 셋째, 신호 처리 — 음성의 선형 예측 부호화, 최대 엔트로피 파워 스펙트럼 추정, 위너 필터·LMS 필터의 예측기 구조. 계수를 빠르게 푸는 수치 알고리즘은 레빈슨-더빈 알고리즘 문서가 따로 다루니, 여기서는 모형으로서의 성질과 맞추는 법에 집중한다.

2. AR(1) — 이산화한 오른슈타인-울렌벡

가장 단순한 $x_t = \varphi x_{t-1} + \varepsilon_t$ 는 $\lvert\varphi\rvert < 1$ 일 때 정상이고

$$\operatorname{Var}(x_t) = \frac{\sigma^2}{1 - \varphi^2}, \qquad \rho_k = \varphi^{k}$$

이다. 자기상관이 지수적으로 감쇠하고, $\varphi \to 1$ 이면 분산이 발산한다. $\varphi = 1$ 이 랜덤 워크(단위근)다.

이 모형은 랑주뱅 동역학의 속도 방정식 같은 오른슈타인-울렌벡 과정 $dX = -\theta X\,dt + s\,dW$ 를 간격 $\Delta$ 로 정확히 샘플링한 것과 같다. 해를 적분하면 $X_{t+\Delta} = e^{-\theta\Delta}X_t + \eta_t$ 이고 $\eta_t$ 는 분산 $s^2(1 - e^{-2\theta\Delta})/2\theta$ 인 정규 잡음이다. 즉 $\varphi = e^{-\theta\Delta}$. 연속시간 확률미분방정식을 오일러-마루야마로 풀면 $\varphi \approx 1 - \theta\Delta$ 가 나오는데, 이건 정확한 지수의 1차 근사일 뿐이다 — 상관 시간을 재야 하는 시뮬레이션에서는 정확한 AR(1) 갱신을 쓰는 것이 국룰이다.

기후 연구에서 AR(1) 은 「적색 잡음」 귀무가설로 쓰인다. 하셀만(1976)의 확률적 기후 모형이 바로 이 구조라, 어떤 기후 시계열의 스펙트럼 봉우리가 의미 있는지는 같은 $\varphi$ 의 AR(1) 스펙트럼 위로 솟았는지로 판정한다.

3. 정상성과 AR(2) 의 진동

일반 AR($p$) 는 지연 연산자 $B x_t = x_{t-1}$ 로 $\varphi(B)x_t = \varepsilon_t$, $\varphi(z) = 1 - \varphi_1 z - \cdots - \varphi_p z^p$ 로 쓴다. $\varphi(z)$ 의 근이 모두 단위원 밖이면 정상이고, 그때 $x_t = \varphi(B)^{-1}\varepsilon_t = \sum_{j\ge0}\psi_j\varepsilon_{t-j}$ 로 무한 이동평균(MA($\infty$)) 표현이 수렴한다. 이것은 선형 점화식 $x_t = \sum\varphi_i x_{t-i}$ 의 안정성 조건과 같은 말이다 — 동반 행렬의 고유값(= $\varphi(z)$ 근의 역수)이 단위원 안에 있어야 한다.

AR(2) 는 정상 영역이 $(\varphi_1, \varphi_2)$ 평면의 삼각형 $\varphi_1 + \varphi_2 < 1$, $\varphi_2 - \varphi_1 < 1$, $\lvert\varphi_2\rvert < 1$ 이고, 그 안에서 $\varphi_1^2 + 4\varphi_2 < 0$ 이면 특성근이 복소수가 되어 잡음에 흔들리는 감쇠 진동자가 된다. 율이 태양흑점에서 본 것이 이것이다¹.

$\varphi_1 = 0.75$, $\varphi_2 = -0.5$, $\sigma = 1$ 로 계산해 보자. 동반 행렬 고유값은 $0.375 \pm 0.599i$, 크기 $\sqrt{0.5} = 0.707$ 이라 정상이다. 고유값 각도 $1.012$ rad 로 「준주기」 $2\pi/1.012 \approx 6.2$ 스텝의 진동이 섞인다. 이론값은

지연 $k$	0	1	2	3	4	5	6
$\rho_k$	1	0.500	−0.125	−0.344	−0.195	0.025	0.117

이고 분산은 $\gamma_0 = \dfrac{1-\varphi_2}{(1+\varphi_2)\bigl[(1-\varphi_2)^2 - \varphi_1^2\bigr]}\sigma^2 = 1.778\sigma^2$ 이다. 자기상관이 부호를 바꾸며 진동 감쇠한다. 한 가지 함정 — 파워 스펙트럼 $\sigma^2/\lvert\varphi(e^{-i\omega})\rvert^2$ 의 봉우리는 $\cos\omega = \varphi_1(\varphi_2 - 1)/(4\varphi_2)$ 에서 $\omega = 0.973$ 으로, 고유값 각도 $1.012$ 와 다르다. 감쇠가 강할수록 둘이 더 벌어지므로 “스펙트럼 봉우리 주기 = 고유 진동 주기”라고 바로 읽으면 안 된다.

4. 식별 — ACF 와 PACF

차수 $p$ 를 데이터에서 짐작하는 고전 도구는 두 개다.

• 자기상관 함수(ACF) $\rho_k$: AR($p$) 에서는 $\varphi(z)$ 근에 따라 지수·진동 감쇠하며 끊기지 않는다.
• 부분 자기상관 함수(PACF) $\alpha_k$: 차수 $k$ 짜리 AR 을 맞췄을 때 마지막 계수 $\varphi_{kk}$. AR($p$) 에서는 $k > p$ 이면 정확히 0 이다. 레빈슨-더빈 점화식이 차수를 올리며 계산하는 반사 계수가 곧 PACF 다.

위 AR(2) 에서 $N = 500$ 개를 생성해(버림 구간 500) 재면 표본 ACF 는 $0.527, -0.098, -0.356, -0.252, -0.014$ 로 이론값을 따라 흔들리고, 표본 PACF 는

$k$	1	2	3	4	5	6
$\hat\alpha_k$	0.527	−0.519	−0.008	−0.027	0.028	−0.002

로 3차부터 신뢰 띠 $\pm 1.96/\sqrt N = \pm 0.088$ 안에 떨어진다. 「PACF 가 $p$ 에서 끊긴다」는 신호가 교과서 그대로 나온다. 이동평균(MA) 모형은 반대로 ACF 가 끊기고 PACF 가 감쇠하므로, 두 그림을 나란히 보는 것이 박스-젠킨스 식별법이다.

5. 추정

5.1. 율-워커

표본 자기공분산 $\hat\gamma_k$ 로 토플리츠 계 $\Gamma_p \boldsymbol\varphi = \boldsymbol\gamma$ 를 푼다. 토플리츠 행렬이라 레빈슨-더빈으로 $O(p^2)$ 에 차수별 해가 한꺼번에 나온다. 편향 표본 자기공분산($1/N$ 정규화)을 쓰면 $\Gamma_p$ 가 항상 양정치라 추정된 모형이 반드시 정상이라는 장점이 있다.

단점은 편향이다. 데이터 양 끝을 0으로 채운 것과 같은 효과가 있어 계수가 0 쪽으로 쏠린다. 단위근 근처에서 특히 심하다.

5.2. 최소제곱(조건부 최대우도)

$x_t$ 를 $(x_{t-1}, \dots, x_{t-p})$ 에 회귀하는 보통 최소자승법이다. 첫 $p$ 개 값을 고정하고 쓴 가우스 우도의 최대화와 같다. 정상성 보장은 없지만 편향이 작다. 정확한 최대우도추정은 첫 $p$ 개 값의 정상 분포까지 우도에 넣는데, 표본이 크면 셋이 거의 같다.

5.3. 실측 — 언제 차이가 나는가

위 AR(2), $N = 500$ 한 표본에서 율-워커는 $(0.800, -0.519)$, 최소제곱은 $(0.805, -0.519)$ 로 차이가 미미하다. 그런데 AR(1), $\varphi = 0.95$, $N = 50$ 으로 2000회 반복하면:

추정법	$\hat\varphi$ 평균	편향
율-워커	0.833	−0.117
최소제곱	0.921	−0.029

율-워커 편향이 네 배다. 짧고 지속성이 강한 시계열(난류 신호의 짧은 구간, 느린 완화를 보이는 MD 관측량)에서는 율-워커를 쓰면 상관 시간을 체계적으로 과소평가한다. 정상성 보장과 작은 편향을 함께 원하면 버그(Burg, 1967) 방법을 쓴다 — 앞·뒤 방향 예측 오차를 함께 최소화하며 반사 계수를 차례로 정해 $\lvert\alpha_k\rvert < 1$ 이 자동으로 지켜진다².

6. 차수 선택 — AIC 와 BIC

차수를 올리면 잔차 분산 $\hat\sigma_p^2$ 은 무조건 줄어든다. 그래서 벌점을 붙인다.

$$\mathrm{AIC}(p) = N\ln\hat\sigma_p^2 + 2p, \qquad \mathrm{BIC}(p) = N\ln\hat\sigma_p^2 + p\ln N$$

아카이케 정보기준은 아카이케(1974)가 바로 이 자기회귀 차수 선택 문제에서 출발해 만든 것이다. 위 표본에서는

$p$	0	1	2	3	4	5
$\hat\sigma_p^2$	2.093	1.513	1.105	1.105	1.104	1.103
AIC	369.4	208.9	53.8	55.7	57.4	59.0
BIC	369.4	213.1	62.2	68.4	74.2	80.0

로 둘 다 $p = 2$ 를 고른다. 하지만 같은 실험을 200번 되풀이해 08차 중에서 고르게 하면 AIC 는 145회(72.5%)만 2를 고르고 나머지는 38 로 과적합하는 반면, 베이즈 정보 기준은 195회(97.5%)를 맞힌다. 이게 이론과 정확히 일치하는 거동이다 — BIC 는 일치 추정량(표본이 커지면 참 차수로 수렴)이고 AIC 는 그렇지 않다. 대신 참 모형이 유한 차수 AR 이 아닐 때 예측 오차를 줄이는 데는 AIC 가 유리하다. 목적이 「모형 구조 식별」이면 BIC, 「예측」이면 AIC 가 국룰이다.

7. 예측과 스펙트럼

7.1. 예측

$h$ 스텝 앞 최적 선형 예측은 점화식을 그대로 돌리며 미래의 잡음을 0으로 두면 된다. 예측 오차 분산은 $\sigma^2\sum_{j<h}\psi_j^2$ 로 커지다가 $h \to \infty$ 에서 과정의 분산 $\gamma_0$ 로 포화한다 — 먼 미래는 평균으로 예측하는 것과 같아진다. 위 AR(2) 의 1스텝 예측 오차 분산은 $\sigma^2 = 1$, 무한 스텝은 $1.778$ 이다. 즉 이 과정에서 과거를 아는 것의 가치는 분산 44% 감소다.

7.2. 최대 엔트로피 스펙트럼

AR 계수로 만든 $\hat S(\omega) = \hat\sigma^2 / \lvert\hat\varphi(e^{-i\omega})\rvert^2$ 는 주어진 처음 $p+1$ 개 자기공분산과 일치하는 스펙트럼 중 엔트로피가 최대인 것이다(버그, 1967). 주기도(피리어도그램)보다 짧은 자료에서 봉우리 분해능이 좋아 지구물리·레이더에서 오래 쓰였다. 다만 차수를 너무 높이면 가짜 봉우리가 갈라지고, 너무 낮으면 봉우리가 뭉개진다. 그래서 차수 선택이 곧 스펙트럼 해상도 선택이다.

8. 시뮬레이션에서 쓰는 법

• 합성 하중 시계열: 목표 파워 스펙트럼에서 자기공분산을 역푸리에로 얻고 율-워커로 AR 계수를 뽑으면, 잡음 한 줄로 원하는 상관 구조의 풍하중·파랑·지진동을 무한히 찍어 낼 수 있다. 같은 일을 FFT 위상 무작위화로도 하지만 AR 은 스트리밍(한 스텝씩 생성)이 된다는 점이 시간 적분 루프에 맞는다.
• 상관된 출력의 오차 막대: MD·MCMC 관측량에 AR(1) 을 맞추면 적분 상관 시간 $\tau = (1+\varphi)/(1-\varphi)$ 와 유효 표본 수 $N/\tau$ 가 바로 나온다. 블로킹·부트스트랩보다 짧은 자료에서 안정적이다.
• 상태공간·필터링: AR($p$) 는 동반 행렬로 $p$ 차원 상태공간 모형이 되므로 칼만 필터에 그대로 들어간다. 결측값 보간·실시간 평활이 이 경로로 처리된다.
• 확장: 이동평균 항을 붙인 ARMA, 차분을 붙인 ARIMA, 여러 변수의 벡터 자기회귀(VAR)와 그 위의 그레인저 인과성 검정, 비선형 판인 임계 AR 등. 순환 신경망은 이 계보의 비선형·고차원 끝판왕이라 볼 수 있다³.

9. 관련 문서

• 레빈슨-더빈 알고리즘 · 토플리츠 행렬 · 최소자승법 · 최대우도추정
• 아카이케 정보기준 · 베이즈 정보 기준 · 부트스트랩
• 파워 스펙트럼 · 위너 필터 · LMS 필터 · 칼만 필터
• 랑주뱅 동역학 · 확률미분방정식 · 순환 신경망
• ARMA · 벡터 자기회귀 · 선형 예측 부호화

10. Footnotes

1. 당시 주류는 흑점 수를 여러 사인파의 합(주기도 분석)으로 설명하려 했고, 율의 요지는 “진자에 콩알을 계속 던지면 주기는 있어 보이지만 위상이 계속 흐트러진다”였다. 주기 성분 피팅이 왜 자꾸 실패하는지를 한 방에 설명해 버린, 꽤 킹받는 반박이었다. ↩

2. 버그의 1967년 결과는 학회 발표(SEG 연례회의)로 나왔고 학위논문(1975)으로 정리됐다. 석유 탐사 회사들이 탄성파 자료 스펙트럼 분석에 갈아넣던 시절 산물이라, 이론보다 현장이 먼저 간 드문 사례다. ↩

3. 실제로 활성화 함수를 항등으로 두고 은닉 상태를 지난 $p$ 개 입력으로 잡은 선형 RNN 은 AR($p$) 와 같은 모형이다. 딥러닝 시계열 논문이 「AR 기준선을 이겼다」고 할 때 그 기준선의 차수를 BIC 로 제대로 골랐는지 확인하는 것이 은근히 중요한 리뷰 포인트다. ↩