앞서 대기오염과 집값의 관계를 살펴봤습니다. 하지만 집값에 영향을 미치는 인자가 워낙 여러개이다보니 단순히 두 가지 변수만으로는 제대로 관계를 설명하기 어렵습니다. 여기서는 여러 가지 변수를 포함한 다중 선형회귀 모델을 만들고 모델의 비교와 평가에 대해서 알아보겠습니다. 물론 실제적으로 그 과정은 상당히 복잡할 수 있습니다. 역시 보스턴 데이터를 불러와 연습 예제로 삼습니다.
library(MASS)
data("Boston", package = "MASS")
data<-boston span="">-boston>
model1=lm(medv~nox, data=data)
summary(model1)
par(mfrow=c(2,2))
plot(model1)
기본 모델과의 비교를 위해 일단 nox 변수 하나만 넣은 모델을 model1이라고 하겠습니다. 여기에 집값에 영향을 미칠 수 있는 몇 가지 변수를 더 넣어 보겠습니다. 범죄율, 가구당 방수, 노후화 주택 비율, 학생당 교사수를 더 넣어 보겠습니다.
model2=lm(medv~nox+crim+rm+age+tax+ptratio, data=data)
summary(model2)
plot(model2)
> model2=lm(medv~nox+crim+rm+age+tax+ptratio, data=data)
> summary(model2)
Call:
lm(formula = medv ~ nox + crim + rm + age + tax + ptratio, data = data)
Residuals:
Min 1Q Median 3Q Max
-13.957 -3.190 -0.524 1.804 40.138
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 4.841171 4.620672 1.048 0.295275
nox -9.402249 3.935328 -2.389 0.017256 *
crim -0.126495 0.036642 -3.452 0.000603 ***
rm 6.929609 0.402242 17.227 < 2e-16 ***
age -0.013521 0.013561 -0.997 0.319210
tax -0.001990 0.002494 -0.798 0.425304
ptratio -0.999508 0.143792 -6.951 1.14e-11 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 5.727 on 499 degrees of freedom
Multiple R-squared: 0.6168, Adjusted R-squared: 0.6122
F-statistic: 133.9 on 6 and 499 DF, p-value: < 2.2e-16
기본 모델이 adjusted R squared 값이 0.181에 불과했던 것과 비교해서 0.6122로 높아지고 잔차도도 많이 개선되었습니다. 몇 개 이상치만 제거하면 그럭저럭 맞는 모델이 만들어질 것 같습니다. 여기에 고급 주택과 연관이 있는 주택당 25000평방피트 (약 700평) 면적 비율을 집어넣으면 어떻게 될까요.
model3=lm(medv~nox+crim+rm+age+tax+ptratio+zn, data=data)
summary(model3)
plot(model3)
> model3=lm(medv~nox+crim+rm+age+tax+ptratio+zn, data=data)
> summary(model3)
Call:
lm(formula = medv ~ nox + crim + rm + age + tax + ptratio + zn,
data = data)
Residuals:
Min 1Q Median 3Q Max
-14.018 -3.111 -0.655 1.931 40.182
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 6.123873 4.880091 1.255 0.210115
nox -10.313773 4.090789 -2.521 0.012006 *
crim -0.125244 0.036685 -3.414 0.000692 ***
rm 6.963393 0.404482 17.216 < 2e-16 ***
age -0.016862 0.014165 -1.190 0.234453
tax -0.001552 0.002551 -0.609 0.543122
ptratio -1.043181 0.153397 -6.801 3e-11 ***
zn -0.012090 0.014755 -0.819 0.412958
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 5.729 on 498 degrees of freedom
Multiple R-squared: 0.6173, Adjusted R-squared: 0.612
F-statistic: 114.8 on 7 and 498 DF, p-value: < 2.2e-16
별로 차이 없는 결과가 나왔습니다. 그런데 이렇게만 봐서는 넣어도 좋은지 아닌지 객관적으로 이해하기 어렵습니다. 모델을 비교하는 방법 중 기본이 되는 것은 두 모델이 결과값에 차이가 있는지를 보는 것입니다. ANOVA를 이용해서 모델간 비교가 가능합니다.
anova(model2, model3)
> anova(model2, model3)
Analysis of Variance Table
Model 1: medv ~ nox + crim + rm + age + tax + ptratio
Model 2: medv ~ nox + crim + rm + age + tax + ptratio + zn
Res.Df RSS Df Sum of Sq F Pr(>F)
1 499 16368
2 498 16346 1 22.037 0.6714 0.413
역시 두 모델간 비교는 별 차이가 없는 것으로 나왔습니다. 그렇다면 굳이 더 설명 변수를 포함할 이유가 없습니다. 현상을 설명할 수 있는 가장 단순한 모델이 가장 좋은 모델이라고 생각되기 때문입니다. 불필요한 변수를 자꾸 넣으면 모델이 복잡해지고 상호 작용의 가능성이 커지면서 오히려 왜곡된 결과를 얻게 될 수 있습니다. 이번엔 큰 주택 대신 저소득 가구 비율을 모델에 포함시켜 보겠습니다.
model4=lm(medv~nox+crim+rm+age+tax+ptratio+lstat, data=data)
summary(model4)
plot(model4)
> model4=lm(medv~nox+crim+rm+age+tax+ptratio+lstat, data=data)
> summary(model4)
Call:
lm(formula = medv ~ nox + crim + rm + age + tax + ptratio + lstat,
data = data)
Residuals:
Min 1Q Median 3Q Max
-13.7019 -3.2454 -0.7473 1.9018 29.0106
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 20.5985664 4.4527073 4.626 4.76e-06 ***
nox -6.0227360 3.5825578 -1.681 0.0934 .
crim -0.0587656 0.0338477 -1.736 0.0832 .
rm 4.4010468 0.4377443 10.054 < 2e-16 ***
age 0.0328562 0.0130724 2.513 0.0123 *
tax 0.0001715 0.0022703 0.076 0.9398
ptratio -0.9377174 0.1305012 -7.186 2.46e-12 ***
lstat -0.5681351 0.0544016 -10.443 < 2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 5.192 on 498 degrees of freedom
Multiple R-squared: 0.6857, Adjusted R-squared: 0.6813
F-statistic: 155.2 on 7 and 498 DF, p-value: < 2.2e-16
이 모델은 Adjusted R squared 값이 0.6813으로 증가되어 훨씬 설명력이 높아 보입니다. 직접 비교해보면 역시 차이가 있습니다.
> anova(model2, model4)
Analysis of Variance Table
Model 1: medv ~ nox + crim + rm + age + tax + ptratio
Model 2: medv ~ nox + crim + rm + age + tax + ptratio + lstat
Res.Df RSS Df Sum of Sq F Pr(>F)
1 499 16368
2 498 13427 1 2940.6 109.06 < 2.2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
모델을 비교하는 또 다른 기준은 AIC (Akaike information criterion) 값입니다. 이 값이 적을수록 모델이 좋은 것으로 평가합니다.
R에서 AIC로 쉽게 비교가 가능합니다.
> AIC(model2, model4)
df AIC
model2 8 3211.085
model4 9 3112.880
model4의 AIC 값이 더 작으므로 model4가 더 좋다는 결론을 내릴 수 있습니다. 그렇지만 이런 식으로 하나하나 다 비교하는 것은 시간 관계상 매우 어려울 수 있습니다. 이럴 때 자주 사용하는 방법이 단계적 회귀분석(stepwise regression) 입니다. 가장 많은 변수부터 제거해나가는 후진법과 가장 적은 변수부터 추가해가는 전진, 그리고 양방향 방법이 가능합니다.
R에서는 구현은 간단합니다.
fit1 <- .="" data="" lm="" medv="" span="">->
fit2 <- 1="" data="" lm="" medv="" span="">->
stepAIC(fit1,direction="backward")
stepAIC(fit2,direction="forward",scope=list(upper=fit1,lower=fit2))
stepAIC(fit2,direction="both",scope=list(upper=fit1,lower=fit2))
각각의 경우 마지막에 나오는 모델이 선택된 모델이며 서로 다른 방법을 사용해도 결국 같은 모델을 얻게 되는 경우가 많습니다.
하지만 이런 자동화된 모델 선택은 여러 가지 문제점이 있습니다. 불필요하게 많은 변수를 넣거나 혹은 반대로 꼭 보정해야 할 변수를 생략할 수 있습니다. 따라서 이미 잘 알려진 위험 인자나 관련 변수를 넣어 모델을 만들거나 혹은 각 그룹별로 유의한 차이가 있는 변수를 넣는 방법 등 여러 가지 다른 방법이 있습니다. 자동화된 단계적 회귀분석은 모델 선정에 있어 주관이 덜 개입되었다고 주장할 순 있지만, 이 역시 초기에 넣는 변수의 종류에 따라 달라질 수 있어 완전히 주관적이지 않은 것은 아닙니다.
댓글
댓글 쓰기