세상에는 단순한 비례 관계를 지닌 현상은 그렇게 많지 않습니다. 특히 자연 법칙이 아니라 상당히 여러 가지 요소가 작용하는 사회 현상, 질병 발생의 경우 그렇습니다. 따라서 이런 경우 회귀 모형에서 뭔가 어긋나는 경우가 생길 수밖에 없습니다. 예를 들어 BMI 1이 증가할 때 마다 당뇨가 정확히 몇 %씩 증가하거나 혹은 소득이 10% 증가할 때마다 자녀의 대학 진학률이 몇 % 씩 선형적으로 증가하기는 어렵다는 것입니다. 당뇨나 대학 진학률에 영향을 미치는 요인은 다양하며 각 요소들은 당연히 정확한 선형 비례 관계를 지니고 있지 않습니다. 더구나 사실 자료의 분포도 정규 분포가 아니라 한 쪽으로 치우친 분포인 경우가 많습니다.
변수 변환은 이렇게 한쪽으로 치우친 변수를 변환해 가능한 회귀 모형을 적합하게 만드는 것입니다. 대표적인 방법은 로그 변환 (log(x)), 역수 변환 (1/x), 루트 변환 (√x), 제곱 변환 (x^2) 등입니다. 로지스틱 회귀 분석에서 쓰이는 로짓 변환도 가능합니다. 앞서 보스턴 자료에서 대기질과 집값과의 관계를 다시 보겠습니다.
아무래도 대기 질 분포가 0.4 - 0.5 사이에 몰려 있는 것으로 보입니다. 이를 해결할 수 있는 가장 좋은 방법은 사실 사분위수로 나눠서 등급으로 바꾸는 것일 수도 있습니다. 아예 범주형으로 바꾸면 분포에 대한 가정은 할 필요가 없으니까요. 아무튼 변수 변환을 통해서 좀 더 회귀 모형을 적합하게 바꿀 수 있을지 살펴보겠습니다. 물론 결과 변수가 아니라 독립 변수를 변환하며 궁극적인 목표는 변수가 아니라 잔차를 적합하게 바꾸는 것입니다. 앞서 설명한 것처럼 중심극한 정리를 통해 충분히 큰 샘플은 정규 분포에 근사한다고 생각할 수 있습니다.
model1=lm(medv~nox, data=data)
summary(model1)
par(mfrow=c(2,2))
plot(model1)
nox에 로그를 한번 씌워 보겠습니다.
nox2=log(data$nox)
> model1=lm(medv~nox2, data=data)
> summary(model1)
Call:
lm(formula = medv ~ nox2, data = data)
Residuals:
Min 1Q Median 3Q Max
-13.289 -5.106 -2.008 2.787 31.529
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 10.537 1.179 8.939 <2e-16 span="">
nox2 -19.665 1.835 -10.717 <2e-16 span="">
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 8.308 on 504 degrees of freedom
Multiple R-squared: 0.1856, Adjusted R-squared: 0.184
F-statistic: 114.9 on 1 and 504 DF, p-value: < 2.2e-16
미미한 차이는 있는지 모르겠지만, 사실 대동소이해 보입니다. 변수 변환에서 어려운 부분은 어떻게 변환을 해야 좋을지 알기 어려운 경우도 많다는 점입니다. 이 경우 car 패키지가 도움이 될 수 있습니다. car 패키지에는 회귀 모형을 적합하기 위한 여러 가지 도구가 있습니다. powerTransform은 변수를 정규화할 수 있는 가장 가능성 있는 지수를 제시합니다.
summary(car::powerTransform(data$nox))
bcPower Transformation to Normality
Est Power Rounded Pwr Wald Lwr Bnd Wald Upr Bnd
data$nox -0.9156 -1 -1.366 -0.4653
Likelihood ratio test that transformation parameter is equal to 0
(log transformation)
LRT df pval
LR test, lambda = (0) 16.38843 1 5.1599e-05
Likelihood ratio test that no transformation is needed
LRT df pval
LR test, lambda = (1) 73.92699 1 < 2.22e-16
Est Power 가 가능한 지수 변환으로 -0.9156이 나왔습니다. 물론 이렇게 지수를 적용하기는 불편하기 때문에 반올림해서 -1 정도가 적당합니다. 마지막 부분은 변수 전환이 필요한지를 보여주는데 필요하지 않을 가능성이 1 < 2.22e-16로 극도로 낮다고 하네요. -1이란 물론 역수를 의미합니다.
nox3=1/data$nox
model1=lm(medv~nox3, data=data)
summary(model1)
par(mfrow=c(2,2))
plot(model1)
이전보다 약간 펴진 것 같다는 생각은 들어도 여전히 잔차도가 100% 마음에 들지는 않습니다. 참고로 이 경우 표준 잔차도 이외에도 앞서 소개한 residplot이 ( https://blog.naver.com/jjy0501/221252592503 참조) 잔차의 분포 변화를 확인하는데 더 좋을 수 있습니다.
한편 car 패키지는 선형성 가정을 위반한 경우의 변형 방식에 대한 정보도 제공합니다. boxTidwell이 그것입니다. 참고로 car 뒤의 :: 표시는 여러 함수가 있을 경우 car 패키지의 것을 사용하라는 지시입니다.
car::boxTidwell(medv~nox, data=data)
> car::boxTidwell(medv~nox, data=data)
MLE of lambda Score Statistic (z) Pr(>|z|)
-0.031947 1.5234 0.1277
iterations = 8
결과 해석에 의하면 선형성을 개선시키기 위해서는 -0.031947라는 매우 독특한 변수 전환이 제시되었습니다. 하지만 중요한 사실은 P 값이 0.1277로 변수 변환이 필요하지 않다라는 점입니다. P 값이 작을 수록 변수 변환이 의미있다는 이야기입니다.
등분산 가정에 위반되는 경우 ncvTest와 spreadLevelPlot 이 도움이 될 수 있습니다.
> ncvTest(model1)
Non-constant Variance Score Test
Variance formula: ~ fitted.values
Chisquare = 3.053339, Df = 1, p = 0.080571
> spreadLevelPlot(model1)
Suggested power transformation: 0.6346467
지수 변환은 0.6이지만 P 값은 역시 필요가 없다는 점을 시사하고 있습니다. 이 함수는 여러 개의 변수를 포함한 모델에서도 동일하게 작동하므로 각각 테스트할 수 있습니다.
다만 애시당초 상관 관계 자체가 비선형인 자료를 자꾸 변환하고 이상치를 제거해서 무리하게 선형으로 끼워맞추면 이상한 결론에 도달할 수 있습니다. 더구나 변수를 변환할 경우 그 결과에 대한 해석이 상당히 복잡해집니다. 따라서 변수 변환 역시 신중하게 결정해야 합니다. 사실 변수 변환이 필요하지 않은 경우가 더 많으며 자동화된 툴 역시 가능성을 설명해 주는 것 뿐이지 꼭 해야 한다는 의미는 아닙니다. 확실하게 모형을 개선할 수 있다는 판단이 서지 않는다면 무리하게 변수 변화를 시도할 이유는 없을 것입니다.
댓글
댓글 쓰기