ADP 실기 — 시계열, 정상성부터 확인하고 들어간다

시계열은 거의 매 회차에 나온다. 17회 코로나, 19회 Traffic EPS, 20회 온도, 23회 코로나, 25회 유입관광객, 26회 은 가격. 문제 형태는 매번 같다 — 탐색 → 정상성 확인 → 처리 → 모형 적합 → 평가.

이론은 ARIMA 를 세 글자로 읽기에 정리해 두었다. 이 글은 코드다.

먼저 그린다 — 롤링 평균과 표준편차

ts = pd.Series(ts_data, index=pd.date_range(start='2021-01-01', periods=1000, freq='D'))

rolling_mean = ts.rolling(30).mean()
rolling_std  = ts.rolling(30).std()

fig, ax = plt.subplots(figsize=(10, 6))
ax.plot(ts, color='gray', label='Original')
ax.plot(rolling_mean, color='blue', label='Rolling Mean')
ax.plot(rolling_std, color='orange', label='Rolling Std')
ax.legend(loc='best')
ax.set_title('Rolling Mean and Standard Deviation')
plt.show()

이 그림 한 장이 정상성 판정의 절반이다. 롤링 평균이 우상향하면 추세가 있는 것이고, 롤링 표준편차가 커지면 분산이 시간에 따라 변하는 것이다. 19회는 “시계열 데이터의 정규성과 이분산성을 설명하기 위한 시각화”를 대놓고 요구했는데, 답이 바로 이 주황색 선이다.

26회의 은 가격 문제는 여기서 끝났다 — 원계열과 3-window 롤링 평균을 한 그래프에 그리고, 1월 대비 9월 증가율을 계산하는 것.

ts.rolling(3).mean()
(ts['2023-09'].mean() / ts['2023-01'].mean() - 1) * 100

정상성 검정 — ADF

그림만으로 쓰지 말고 검정을 같이 붙인다.

from statsmodels.tsa.stattools import adfuller

result = adfuller(ts)
print('ADF Statistic:', result[0])
print('p-value:', result[1])

ADF 는 귀무가설이 “단위근이 있다 = 비정상”이다. 그래서 방향이 다른 검정들과 반대다 — $p < 0.05$ 여야 정상이다. 정규성 검정과 헷갈려서 반대로 쓰면 결론이 통째로 뒤집힌다.

차분

diff_ts = ts.diff(periods=1)
diff_ts.dropna(inplace=True)

차분하면 첫 값이 NaN 이 되므로 dropna() 를 반드시 붙인다. 안 붙이면 뒤에서 ACF 가 에러를 낸다.

로그를 먼저 씌우면 분산도 같이 잡힌다. 값이 커질수록 변동폭도 커지는 데이터(매출·확진자)에 쓴다.

ts_log = np.log(ts).dropna()

분해 — 추세 · 주기 · 계절성

import statsmodels.api as sm

decomposition = sm.tsa.seasonal_decompose(ts_log, period=30)

trend       = ts_log.rolling(window=12).mean().dropna()
cycle       = decomposition.resid.dropna()
seasonality = decomposition.seasonal.dropna()

fig, ax = plt.subplots(nrows=3, figsize=(12, 10))
ax[0].plot(ts_log, label='Original'); ax[0].plot(trend, label='Trend'); ax[0].legend()
ax[1].plot(cycle, label='Cycle'); ax[1].legend()
ax[2].plot(seasonality, label='Seasonality'); ax[2].legend()
plt.show()

period 를 반드시 준다. 일별 데이터의 주 단위 계절성이면 7, 월별 데이터의 연 단위면 12다. 주기를 모르는 채로 분해하면 아무 의미가 없다 — 그래서 분해 전에 데이터 간격이 무엇인지부터 확인한다.

model='additive' 가 기본이고, 계절 변동폭이 수준에 비례해 커지면 model='multiplicative' 다.

날짜에서 성분을 뽑을 때는 .dt 접근자를 쓴다.

tdf['date'].dt.year
tdf['date'].dt.month
tdf['date'].dt.day
tdf['date'].dt.day_name()
tdf['date'].dt.quarter

19회가 “매년 분기별로 작성된 20년치” 데이터였다. .dt.quarter 로 접는 것이 첫 단계였다.

ACF 와 PACF — p, q 를 읽는다

from statsmodels.graphics.tsaplots import plot_acf, plot_pacf

fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(8, 8))
plot_acf(diff_ts,  ax=ax1, lags=50)
plot_pacf(diff_ts, ax=ax2, lags=50)
plt.show()

차분한 계열에 그린다. 원계열에 그리면 추세 때문에 ACF 가 천천히 감소하는 모양만 나오고, 거기서는 아무것도 못 읽는다.

	ACF	PACF
AR(p)	천천히 감소	lag p 이후 뚝 끊긴다
MA(q)	lag q 이후 뚝 끊긴다	천천히 감소
ARMA	둘 다 천천히 감소	둘 다 천천히 감소

끊기는 쪽(절단)을 보고 차수를 읽는다. 파란 신뢰띠 밖으로 나온 마지막 lag 가 차수다.

ARIMA

from statsmodels.tsa.arima.model import ARIMA

model = ARIMA(ts, order=(1, 1, 1))
result = model.fit()
print(result.summary())

forecast = result.predict(start='2022-01-01', end='2022-12-31')

fig, ax = plt.subplots(figsize=(10, 6))
ax.plot(ts, label='Original')
ax.plot(forecast, color='green', label='Forecast')
ax.legend(loc='best')
plt.show()

예전 메모에는 from statsmodels.tsa.arima_model import ARIMA 와 .fit(disp=0) 으로 적어 두었다. 그 경로는 statsmodels 0.12 에서 폐기 예고되고 0.13 에서 빠졌다. 시험 환경이 0.13.2 이므로 statsmodels.tsa.arima.model 을 쓰고, disp 인자도 넣지 않는다. 새 경로는 옛 버전에서도 도니 그냥 이쪽만 외우는 편이 낫다.

차수를 손으로 못 정하겠으면 pmdarima 가 있다.

import pmdarima as pm

model = pm.auto_arima(ts, seasonal=True, m=12, trace=True,
                      error_action='ignore', suppress_warnings=True)
print(model.summary())

AIC 를 기준으로 자동 탐색한다. 19회의 “여러 파라미터를 적용해 보고 가장 성능이 좋은 것을 제시”에 그대로 쓸 수 있다. 다만 답안에는 자동 탐색 결과만 붙이지 말고, ACF/PACF 로 읽은 차수와 비교해서 한 문단을 쓴다.

SARIMA

계절성이 있으면 이쪽이다. 25회의 “계절성 있는 시계열 모델을 적합하기”가 이것이다.

from statsmodels.tsa.statespace.sarimax import SARIMAX

model = SARIMAX(data, order=(1, 1, 1), seasonal_order=(1, 1, 1, 12))
result = model.fit()

predictions = result.get_prediction(start=pd.to_datetime('2022-01-01'),
                                    end=pd.to_datetime('2023-01-01'),
                                    dynamic=False)
predicted_values = predictions.predicted_mean

plt.plot(data, label='Actual')
plt.plot(predicted_values, color='red', linestyle='--', label='Predicted')
plt.legend(); plt.show()

$\text{SARIMA}(p,d,q)(P,D,Q)_s$ 에서 뒤 괄호가 주기 $s$ 에 대한 같은 구조다. 월별 데이터면 $s = 12$, 분기별이면 4, 일별의 주 단위면 7이다.

dynamic=False 는 매 시점 예측에 실제 과거값을 쓴다는 뜻이다(한 스텝 예측). dynamic=True 면 예측값을 다시 입력으로 넣는다. 검증할 때는 False, 진짜 미래를 볼 때는 True 쪽에 가깝다.

평가

from sklearn.metrics import mean_squared_error, r2_score

mse  = mean_squared_error(ts[start:end], forecast)
rmse = np.sqrt(mse)
print('RMSE:', rmse)
print('R2:', r2_score(ts[start:end], forecast))

분포를 나란히 그려 놓으면 “예측이 실제보다 폭이 좁다”같은 말을 할 수 있다.

fig, ax = plt.subplots(1, 2, figsize=(12, 4))
sns.histplot(ts[start:end], ax=ax[0], color='gray'); ax[0].set_title('Original')
sns.histplot(forecast, ax=ax[1], color='green');     ax[1].set_title('Forecast')
plt.show()

잔차 진단 — 19회의 4번 문제

“위 모델의 잔차와 잡음 시각화 및 분석”이 따로 배점이었다. 좋은 시계열 모형의 잔차는 백색잡음이어야 한다 — 평균 0, 일정한 분산, 자기상관 없음.

result.plot_diagnostics(figsize=(12, 8))    # 잔차·히스토그램·QQ·ACF 를 한 번에
plt.show()

from statsmodels.stats.diagnostic import acorr_ljungbox

print(acorr_ljungbox(result.resid, lags=[10], return_df=True))

Ljung-Box 는 귀무가설이 “자기상관이 없다”이므로 $p > 0.05$ 여야 한다. $p$ 가 작으면 잔차에 아직 구조가 남았다는 뜻이고, 차수를 다시 봐야 한다.

plot_diagnostics 한 줄이면 그림 네 장이 나온다. 시간이 없을 때 이것만 붙여도 배점을 꽤 가져간다.

시계열을 회귀로 푸는 경우

20회는 온도 데이터에 지연값(lag) 이 이미 열로 들어 있었고, 랜덤포레스트와 SVM 으로 예측하라고 했다. ARIMA 가 아니라 지도학습 문제였다.

df['lag1'] = df['temp'].shift(1)
df['lag2'] = df['temp'].shift(2)
df['rolling3'] = df['temp'].rolling(3).mean().shift(1)
df = df.dropna()

shift(1) 이 반드시 붙는다. 롤링 평균에 현재 값이 들어가면 예측하려는 값을 입력으로 쓰는 셈이라 성능이 비현실적으로 좋아진다. 시계열에서 가장 흔한 누수다.

분할도 다르다. 시계열은 무작위로 섞어 나누지 않는다.

split = int(len(df) * 0.7)
train, test = df[:split], df[split:]

또는 TimeSeriesSplit 을 쓴다.

from sklearn.model_selection import TimeSeriesSplit
tscv = TimeSeriesSplit(n_splits=5)

20회가 “train/test set 나눌 방법을 설명하라”고 한 이유가 여기 있다. train_test_split 에 shuffle=True 를 그대로 쓰면 미래 데이터로 과거를 예측하게 되고, 그 한 줄로 문제 전체가 틀린 답이 된다.