회귀모형의 평가
회귀분석은 확률에 기반한 것으로 그 모형에 의한 추정값은 관측값과 차이를 발생시킵니다. 모형의 평가는 그 차이의 수준을 용인할 수 있는지에 대한 평가로서 앞서 소개한 분산분석을 적용합니다. 분산분석은 여러개의 그룹 (변수)들 사이에 일어나는 각각의 변동(분산)을 비교하여 일반적으로 일어날 수 있는 수준인지를 판단하는 것입니다.
그림 1에 나타낸 것과 같이 관측치 y의 불편추정치(unbiased estimator)로 그 값들의 평균 $\bar{y}$이 사용됩니다. 평균값이 회귀모형에 의한 예측치 $\hat{y}$와 일치한다면 회귀분석의 의미는 없어집니다. 즉, 회귀모델이 적합하다면 평균과 추정치 사이에 차이가 발생하며 추정치와 관측치 사이에 오차가 발생됩니다. 적합한 회귀모형에 의한 반응변수의 평균과 예측값 그리고 관측값(y) 사이의 관계는 식 1과 같이 정의할 수 있습니다.
| $$(\bar{y}-y)^2=(\bar{y}-\hat{y})^2+(\bar{y}-y)^2+\alpha$$ | (식 1) |
x=np.linspace(-1, 2, 100)
y=x+0.5
plt.figure(figsize=(4,3))
plt.plot(x, y, color="g", label="regression")
plt.hlines(1.7, -1, 2, color="k", ls="--", label="mean line")
plt.scatter(0.25, 1.7, s=20, color="k", label=r"$\bar{y}$")
plt.scatter(0.25, 0.75, s=20, color="r", label=r"$\hat{y}$")
plt.scatter(0.25, 0, s=20, color="b", label=r"$y_{obs}$")
plt.vlines(0.25, 0.75, 1.7, color="r", ls="--")
plt.vlines(0.25, 0, 0.75, color="b", ls="--")
plt.vlines(0.18, 0, 1.7, color="k", ls="--", alpha=0.7)
plt.text(0.3, 1.2, r"SSReg=$(\bar{y}-\hat{y})^2$", color="r", weight="bold", size=8)
plt.text(0.3, 0.25, r"SSE=$(y_{obs}-\hat{y})^2$", color="b", weight="bold", size=8)
plt.text(-0.9, 0.7, r"SST=$(y_{obs}-\bar{y})^2$", color="k", weight="bold", size=8)
plt.xlabel("x")
plt.ylabel("y")
plt.xlim(-1.1, 2.5)
plt.legend(loc="lower right", prop={'size':8}, frameon=False)
plt.show()
식 1은 관측치 y에 대해 반응변수의 평균과 y에 대응하는 예측치 사이의 편차를 계산한 것으로 식 2와 같이 그 변수의 모든 샘플에 대해 편차의 합을 계산할 수 있습니다.
| \begin{align}\sum^n_{i=1}(\bar{y}-y)^2&\approx \sum^n_{i=1}\sum^n_{i=1}(\bar{y}-\hat{y})^2+\sum^n_{i=1}(\bar{y}-y)^2\\\text{SST}&\approx \text{SSReg}+\text{SSE} \\ \text{총제곱합}&\approx \text{설명될 수 있는 부분}+\text{설명될 수 없는 부분}\\\bar{y},\; \hat{y}:&\; \text{평균, 추정치}\\\end{align} | (식 2) |
| SST ≈ SSReg + SSE | |
| 총제곱합 ≈ 설명될수 있는 부분 + 설명될 수 없는 부분 |
식 2의 좌항은 총제곱합(Sum square total, SST)이고 우항은 회귀제곱합(Sum square Regression, SSReg)과 잔차제곱합(Sum square error, SSE)의 합입니다. 회귀제곱합(SSReg)은 모델에 의한 예측치와 반응변수의 평균값(고정된 값)과의 차이로 회귀모형에 의존됩니다. 즉, 이 오차의 원인은 회귀계수에 기인하는 것으로 설명이 가능한 오차이지만 잔차제곱합(SSE)는 확률변수인 관측값과 회귀모형 모두에 의존되므로 설명할 수 없는 오차입니다.
예 1)
kospi 지수의 일일 주가 자료중 시가(Open)을 설명변수로 하여 종가(Close)를 추정하는 회귀모델을 작성하고 MSE를 계산합니다.
| Open | Close | |
|---|---|---|
| 0 | 2874.50 | 2944.45 |
| 1 | 2943.67 | 2990.57 |
| 2 | 2993.34 | 2968.21 |
| ⋮ | ⋮ | ⋮ |
다음 코드는 분석을 위한 자료를 호출하기 위한 것입니다.
st=pd.Timestamp(2021,1, 1)
et=pd.Timestamp(2024, 5, 10)
kos=fdr.DataReader('KS11',st, et)[["Open","Close"]]
kos.index=range(len(kos))
kos.head(3).round(2)
| Open | Close | ||
|---|---|---|---|
| 0 | 2201.21 | 2175.17 | |
| 1 | 2192.58 | 2176.46 | |
| 2 | 2154.97 | 2155.07 |
위 자료를 표준화합니다.
X=kos.values[:,0].reshape(-1,1) y=kos.values[:,1].reshape(-1,1) from sklearn.preprocessing import StandardScaler #독립변수 정규화(표준화) xScaler=StandardScaler().fit(X) X_n=xScaler.transform(X) #반응변수 정규화(표준화) yScaler=StandardScaler().fit(y) y_n=yScaler.transform(y)
자료의 회귀모델을 구현하기 위해 sklearn.liniear_model.linearRegression(fit_intercept=True) 클래스 적용합니다.
from sklearn.linear_model import LinearRegression mod = LinearRegression() mod.fit(X_n, y_n) pre=mod.predict(X_n) pre[:3].round(3)
array([[0.559],
[0.788],
[0.953]])
n=len(y_n) mu=y_n.mean() pre=mod.predict(X_n) ssreg=np.sum((pre-mu)**2) round(ssreg, 2)
822.19
sse=np.sum((pre-y_n)**2) round(sse, 2)
4.81
sst=np.sum((y_n-mu)**2) round(sst, 2)
827.0
round((ssreg+sse)-sst,5)
0.0
위 결과는 statmodels.api.OLS() 클래스 적용한 모델 reg의 속성 ess과 ssr로 확인할 수 있습니다. 이 속성은 각각 sse, ssreg를 나타냅니다.
import statsmodels.api as sm X_n0=sm.add_constant(X_n) X_n0.shape, y_n.shape reg=sm.OLS(y_n, X_n0).fit()
print("sse: %.2f, rss: %.2f" %(reg.ess, reg.ssr))
sse: 822.19, rss: 4.81
모형의 적합성은 SSE가 감소할수록 높아집니다. 즉, 식 3에서 나타낸 것과 같이 전체 오차 중에 설명할 수 있는 오차의 비가 증가할수록 모형의 적합도는 증가합니다. 이 결과를 결정계수(coefficient of determinant, R2)라고 하며 회귀모형에 의한 설명력을 나타내는 지표로 사용됩니다.
| \begin{align}R^2&= \frac{\text{SSReg}}{\text{SST}}\\ & = \frac{\text{SST}-\text{SSE}}{\text{SST}}\\& = 1- \frac{\text{SSE}}{\text{SST}}\\& =1-\frac{\sum^n_{i=1}(y-\hat{y})^2}{\sum^n_{i=1}(y-\bar{y})^2}\end{align} | (식 3) |
결정계수는 r의 제곱과 같습니다. 이 결정계수는 [0, 1] 사이의 값으로 1에 근접할 수록 모형의 추정치가 관찰치에 근접함을 나타냅니다.
r2=ssreg/sst round(r2, 4)
0.9942
statmodels.api.OLS()에 의해 생성된 모델의 rsquared 속성과 sklearn.linear_model.LinearRegression으로 생성된 모델의 score() 메소드에 의해 계산할 수 있습니다.
print(f'결정계수 by OLS() : {round(reg.rsquared, 3)}')
결정계수 by OLS() : 0.994
R2=mod.score(X_n, y_n)
print(f'결정계수 by LinearRegression() : {round(R2, 3)}')
결정계수 by LinearRegression() : 0.994
설명할 수 있는 오차제곱(SSReg)과 설명할 수 없는 오차제곱(SSE)을 각각의 자유도(표 1)로 나눈 값들을 MSReg와 MSE라고 하며 각 오차분포의 분산이 됩니다.
| SS | df | MS |
|---|---|---|
| SST | 자료수 - 1 | $\text{MST} = \frac{\text{SST}}{\text{자료수-1}}$ |
| SSReg | 설명변수의 수 - 1 | $\text{MSReg} = \frac{\text{SSReg}}{\text{설명변수의 수-1}}$ |
| SSE | 자료수 - 설명변수의 수 | $\text{MSE} = \frac{\text{SSE}}{\text{자료수 - 설명변수의 수}}$ |
분산들의 비는 F분포를 따르며 그 통계량에 따라 두 분포의 동일성에 대한 검정을 실시할 수 있습니다. 즉, 식 4와 같이 계산된 F 통계량에 대응하는 p-value와 신뢰구간에 따라 두 분포의 유효한 차이의 유무를 판단할 수 있습니다.
| F 통계량 = | MSReg | (식 4) |
| MSE | ||
| H0: b1 = 0, | H1: b1 ≠ 0 |
위 모델 mod에 대한 f 검정을 실시합니다.
mst=sst/(len(y_n)-1)
print(f'mst: {round(mst, 3)}')
mst: 1.001
msreg=ssreg/1
print(f'msreg: {round(msreg, 3)}')
msreg: 822.193
mse=sse/(len(y_n)-2)
print(f'mse: {round(mse, 3)}')
mse: 0.006
fratio=msreg/mse
print(f'F 통계량: {round(fratio, 3)}')
F 통계량: 141094.943
pval=stats.f.sf(msreg/mse, 1, len(y_n)-2)
print(f'p-val: {round(pval, 3)}')
p-val: 0.0
ci=stats.f.interval(0.95, 1, len(y_n)-2)
print(f'하한: {round(ci[0], 3):}, 상한: {round(ci[1], 3):}')
하한: 0.001, 상한: 5.042
위 결과는 유의수준 0.05를 기준으로 신뢰구간과 유의확률 기준으로 귀무가설을 채택할 수 없음을 나타냅니다.위 계산은 sklearn.feature_selection.f_regression 함수로 대신할 수 있습니다. 이 함수는 검정통계량과 유의확률을 반환합니다.
from sklearn.feature_selection import f_regression
re=f_regression(pre, y_n.ravel())
print(f'통계량: {np.around(re[0], 3)}, p-value: {np.around(re[1], 3)}')
통계량: [141094.943], p-value: [0.]
댓글
댓글 쓰기