회귀분석에 사용되는 데이터는 그 자체로 사용하기 보다는 스케일링이나 함수 변환 등의 전처리 과정을 거치는 경우가 많다. 전처리 과정은 공분산 행렬의 조건을 향상시키거나 데이터 간의 관계를 선형 모형에 맞게 바꾸기 위해 사용된다.
조건수(conditional number)는 공분산 행렬 $X^TX$의 가장 큰 고유치와 가장 작은 고유치의 비율을 뜻한다.
$$ \text{condition number} = \dfrac{\lambda_{\text{max}}}{\lambda_{\text{min}}} $$조건수가 크면 역행렬을 계산할 때 오차가 미치는 영향이 커진다. 여기에서는 다음 연립방정식을 예로 들어 설명을 하겠다.
$$ Ax = b $$행렬 $A$가 단위 행렬이면 조건수는 가장 작은 경우로 조건수의 값이 1이다.
$$ \text{cond}(I) = 1 $$A = np.eye(4)
이 행렬 $A$와 곱해져서 상수 벡터 $b$가 되는 벡터 $x$를 역행렬 $A^{-1}$을 사용하여 계산할 수 있다. 이 예에서는 상수 벡터 $b$가 1-벡터이다.
b = np.ones(4)
sp.linalg.solve(A, b)
만약 상수 벡터에 약간의 오차가 있었다면 연립방정식의 해에도 동일한 수준의 오차가 발행한다.
sp.linalg.solve(A + 0.0001 * np.eye(4), b)
다음과 같은 행렬을 생각하자. 이 행렬은 4차 힐버트 행렬로 조건수가 15000이 넘는다. 이렇게 연립방정식을 이루는 행렬의 조건수가 커지면 상수항 오차가 작은 경우에도 해의 오차가 커지게 된다.
A = sp.linalg.hilbert(4)
A
np.linalg.cond(A)
이 행렬과 곱해져서 상수 벡터가 되는 벡터를 역행렬을 사용하여 찾으면 다음과 같다.
sp.linalg.solve(A, b)
조건수가 크면 약간의 오차만 있어도 해가 전혀 다른 값을 가진다. 따라서 조건수가 크면 회귀분석을 사용한 예측값도 오차가 커지게 된다.
sp.linalg.solve(A + 0.0001 * np.eye(4), b)
회귀분석에서 조건수가 커지는 경우는 크게 두 가지가 있다.
보스턴 집값 데이터의 경우 회귀분석을 하면 조건수가 15,000 정도로 크게 나온다.
from sklearn.datasets import load_boston
boston = load_boston()
dfX = pd.DataFrame(boston.data, columns=boston.feature_names)
dfy = pd.DataFrame(boston.target, columns=["MEDV"])
df = pd.concat([dfX, dfy], axis=1)
model1 = sm.OLS.from_formula("MEDV ~ " + "+".join(boston.feature_names), data=df)
result1 = model1.fit()
print(result1.summary())
그 이유는 각 독립 변수들이 0.1 단위부터 수백 단위까지 제각각의 크기를 가지고 있기 때문이다.
dfX.describe().iloc[[3, 7], :]
이렇게 조건수가 크면 모수 추정 오차가 증폭될 가능성이 커진다. 이 효과를 확실하게 보기 위하여 일부러 다음처럼 조건수를 더 크게 만들어 보았다. 이 상태로 선형 회귀분석을 하면 성능이 급격하게 떨어진다.
summary 리포트에서 성능을 나타내는 지표는 R-squared라는 이름의 결정계수다. 결정계수에 대해서는 분산분석에서 자세히 설명한다. 결정계수 값이 원래의 모형에서는 0.741이었는데 조건이 나빠지면서 0.332로 감소한 것을 볼 수 있다.
dfX2 = dfX.copy()
dfX2["TAX"] *= 1e13
df2 = pd.concat([dfX2, dfy], axis=1)
model2 = sm.OLS.from_formula("MEDV ~ " + "+".join(boston.feature_names), data=df2)
result2 = model2.fit()
print(result2.summary())
StatsModels에서는 모형지정 문자열에서 scale 명령을 사용하여 스케일링을 할 수 있다. 이 방식으로 스케일을 하면 스케일링에 사용된 평균과 표준편차를 저장하였다가 나중에 predict 명령을 사용할 때도 같은 스케일을 사용하기 때문에 편리하다. 더미 변수인 CHAS는 스케일을 하지 않는다는 점에 주의한다.
feature_names = list(boston.feature_names)
feature_names.remove("CHAS")
feature_names = ["scale({})".format(name) for name in feature_names] + ["CHAS"]
model3 = sm.OLS.from_formula("MEDV ~ " + "+".join(feature_names), data=df2)
result3 = model3.fit()
print(result3.summary())
항상 감사합니다. 잘 보고 배우고 있습니다.
질문이 있습니다.
model3 = sm.OLS.from_formula("np.log(MEDV) ~ "
"scale(CRIM) + scale(I(CRIM ** 2)) + "
"scale(ZN) + scale(I(ZN ** 2)) + "
"scale(INDUS) + "
"scale(NOX) + scale(RM) + scale(AGE) + "
"scale(np.log(DIS)) + scale(RAD) + scale(TAX) + "
"scale(np.log(PTRATIO)) + scale(B) + scale(np.log(LSTAT)) + CHAS",
data=df)
여기에서
1. I(CRIM ** 2), I(ZN**2) 의 I()는 어떤 함수인가요?
2. CRIM, ZN 을 제곱하는 이유가 있나요?? 제곱을 해야하는 경우와 log를 씌워줘야 하는 경우의 차이를 알고 싶습니다.
3. CRIM이나 ZN의 경우는 "scale(CRIM) + scale(I(CRIM ** 2)) " 처럼 두 개의 feature를 사용하는데, 이렇게 하는 이유가 궁금합니다.
3. 지수화 또는 로그화를 안해주는 feature들은 어떤 기준으로 scale만 해주는 건가요?
사용자에 의해 삭제되었습니다.