앞서 피마 인디언 예제에서 BMI와 당뇨의 발생을 알아봤는데, 사실 BMI 를 연속 변수로 취급할수도 있지만, 범주형으로 바꿔서 볼수도 있습니다. 사실 이쪽이 더 유용한 결과를 제시하기도 합니다. BMI가 1단위 늘어날때마다 당뇨 위험도가 계속 같은 비율로 상승하지는 않을 것이기 때문입니다. OR값은 대략적인 평균 추정치에 지나지 않는데, 비만, 과체중, 정상 체중으로 그룹을 나누면 더 유용한 결과를 얻을지도 모릅니다.
즉 정상 체중인 사람에 비해 비만이나 과체중인 경우 당뇨 위험도가 얼마나 상승하는지 알 수 있는 것이죠. 실제 BMI가 이런식으로 사용된다는 점을 생각하면 더 유용성이 높고 정규 분포에서 벗어난 결과 변수의 경우 범주형으로 바꾸는 편이 더 정확할 수 있다는 점도 생각할 수 있습니다. 여기서는 25와 30을 기준으로 정상, 과체중, 비만으로 나누겠습니다. obesity라는 새로운 범주형 변수를 만듭니다.
library(mlbench)
data(PimaIndiansDiabetes)
pima <- pimaindiansdiabetes="" span="">
summary(pima)
pima$obesity[pima$mass>=30]=2
pima$obesity[pima$mass<30 span="">
pima$obesity[pima$mass<25 span="">
table(pima$obesity)
> table(pima$obesity)
0 1 2
117 179 472
연구 대상 중 정상 체중은 별로 없고 비만 비율이 높다는 것을 알 수 있습니다. 로지스틱 회귀 분석에 넣기 위해서는 factor라는 사실을 인식하기 위해 factor(obesity)라고 해야 합니다.
library(moonBook)
out=glm(diabetes~factor(obesity),family=binomial,data=pima)
summary(out)
extractOR(out)
> summary(out)
Call:
glm(formula = diabetes ~ factor(obesity), family = binomial,
data = pima)
Deviance Residuals:
Min 1Q Median 3Q Max
-1.1168 -1.1168 -0.7112 1.2393 2.2649
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -2.4849 0.3469 -7.162 7.94e-13 ***
factor(obesity)1 1.2393 0.3906 3.173 0.00151 **
factor(obesity)2 2.3406 0.3590 6.520 7.05e-11 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 993.48 on 767 degrees of freedom
Residual deviance: 905.53 on 765 degrees of freedom
AIC: 911.53
Number of Fisher Scoring iterations: 5
> extractOR(out)
OR lcl ucl p
(Intercept) 0.08 0.04 0.16 0.0000
factor(obesity)1 3.45 1.61 7.43 0.0015
factor(obesity)2 10.39 5.14 20.99 0.0000
정상 체중에 비해 과체중인 경우 오즈비가 3.45, 비만이면 10.39라는 결과가 나오는데, 그만큼 비만이 당뇨에 미치는 영향이 크다는 결과입니다. BMI를 연속변수로 취급했을 때 보다 더 이해하기 쉬운 결과입니다.
하지만 우리가 OR값을 그대로 받아들이기에는 몇 가지 문제가 있습니다. 가장 큰 문제는 당뇨에 영향을 미치는 인자가 BMI 하나만이 아니라는 점입니다. 가능성 있는 변수를 넣은 다중 로지스틱 회귀 분석 (multiple logistic regression analysis)가 사실 더 널리 사용되는 방법이라고 하겠습니다. 피마 인디언에서 구한 주요 변수를 모두 넣은 모델을 만들어 보겠습니다.
out=glm(diabetes~factor(obesity)+pressure+triceps+insulin+glucose+age+pedigree,family=binomial,data=pima)
summary(out)
extractOR(out)
> summary(out)
Call:
glm(formula = diabetes ~ factor(obesity) + pressure + triceps +
insulin + glucose + age + pedigree, family = binomial, data = pima)
Deviance Residuals:
Min 1Q Median 3Q Max
-2.4705 -0.7428 -0.3837 0.7104 3.0096
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -7.4833720 0.6912973 -10.825 < 2e-16 ***
factor(obesity)1 1.1652222 0.4314899 2.700 0.00692 **
factor(obesity)2 2.2445087 0.4077293 5.505 3.69e-08 ***
pressure -0.0101593 0.0051619 -1.968 0.04905 *
triceps 0.0010748 0.0067944 0.158 0.87430
insulin -0.0015264 0.0009212 -1.657 0.09753 .
glucose 0.0349252 0.0037025 9.433 < 2e-16 ***
age 0.0317551 0.0082582 3.845 0.00012 ***
pedigree 0.9370598 0.2979485 3.145 0.00166 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 993.48 on 767 degrees of freedom
Residual deviance: 725.79 on 759 degrees of freedom
AIC: 743.79
Number of Fisher Scoring iterations: 5
> extractOR(out)
OR lcl ucl p
(Intercept) 0.00 0.00 0.00 0.0000
factor(obesity)1 3.21 1.38 7.47 0.0069
factor(obesity)2 9.44 4.24 20.98 0.0000
pressure 0.99 0.98 1.00 0.0491
triceps 1.00 0.99 1.01 0.8743
insulin 1.00 1.00 1.00 0.0975
glucose 1.04 1.03 1.04 0.0000
age 1.03 1.02 1.05 0.0001
pedigree 2.55 1.42 4.58 0.0017
비슷한 결과지만 OR 값은 약간 줄어들었다는 점을 알 수 있습니다. 하지만 여기서 주의할 점이 데이터를 그냥 넣으면 안된다는 점입니다. 사실은 결측값이 0으로 들어간 것 때문에 잘못된 결과이며 먼저 결측치를 제거해야 합니다. 만약 결측치를 제거하지 않을 경우 사실은 단순 로지스틱 회귀 분석과 다중 로지스틱 회귀 분석에 들어간 인구 집단이 서로 달라질 수 있습니다. R에서는 결측치는 무조건 제외하고 분석을 진행합니다.
pima<-subset pima="" pressure="">0)
pima<-subset pima="" triceps="">0)
pima<-subset mass="" pima="">0)
pima<-subset insulin="" pima="">0)
pima<-subset glucose="" pima="">0)
out=glm(diabetes~factor(obesity),family=binomial,data=pima)
summary(out)
extractOR(out)
out=glm(diabetes~factor(obesity)+pressure+triceps+insulin+glucose+age+pedigree,family=binomial,data=pima)
summary(out)
extractOR(out)
> out=glm(diabetes~factor(obesity),family=binomial,data=pima)
> summary(out)
Call:
glm(formula = diabetes ~ factor(obesity), family = binomial,
data = pima)
Deviance Residuals:
Min 1Q Median 3Q Max
-1.0435 -1.0435 -0.6899 1.3175 2.4954
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -3.0681 0.7233 -4.242 2.22e-05 ***
factor(obesity)1 1.7537 0.7705 2.276 0.022833 *
factor(obesity)2 2.7447 0.7340 3.739 0.000185 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 498.10 on 391 degrees of freedom
Residual deviance: 460.58 on 389 degrees of freedom
AIC: 466.58
Number of Fisher Scoring iterations: 5
> extractOR(out)
OR lcl ucl p
(Intercept) 0.05 0.01 0.19 0.0000
factor(obesity)1 5.78 1.28 26.15 0.0228
factor(obesity)2 15.56 3.69 65.58 0.0002
>
> out=glm(diabetes~factor(obesity)+pressure+triceps+insulin+glucose+age+pedigree,family=binomial,data=pima)
> summary(out)
Call:
glm(formula = diabetes ~ factor(obesity) + pressure + triceps +
insulin + glucose + age + pedigree, family = binomial, data = pima)
Deviance Residuals:
Min 1Q Median 3Q Max
-2.7631 -0.6768 -0.3561 0.6357 2.6944
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -9.8526920 1.3011873 -7.572 3.67e-14 ***
factor(obesity)1 1.2443047 0.8266796 1.505 0.132277
factor(obesity)2 1.9356374 0.8147969 2.376 0.017520 *
pressure 0.0041776 0.0117286 0.356 0.721701
triceps 0.0156416 0.0154807 1.010 0.312310
insulin -0.0007836 0.0013274 -0.590 0.554948
glucose 0.0371823 0.0057950 6.416 1.40e-10 ***
age 0.0473885 0.0139856 3.388 0.000703 ***
pedigree 1.1158455 0.4205538 2.653 0.007971 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 498.10 on 391 degrees of freedom
Residual deviance: 343.38 on 383 degrees of freedom
AIC: 361.38
Number of Fisher Scoring iterations: 6
> extractOR(out)
OR lcl ucl p
(Intercept) 0.00 0.00 0.00 0.0000
factor(obesity)1 3.47 0.69 17.54 0.1323
factor(obesity)2 6.93 1.40 34.21 0.0175
pressure 1.00 0.98 1.03 0.7217
triceps 1.02 0.99 1.05 0.3123
insulin 1.00 1.00 1.00 0.5549
glucose 1.04 1.03 1.05 0.0000
age 1.05 1.02 1.08 0.0007
pedigree 3.05 1.34 6.96 0.0080
>
결측치를 제거하고 보니 OR값이 3.47과 6.93으로 줄어들었습니다. 그런데 문제는 샘플 사이즈가 작아지면서 95% 신뢰구간이 넓게 벌어졌다는 것입니다. 그 결과 과체중 (obesity1)의 통계적 유의성이 사라졌습니다. 이 문제를 해결할 방법을 고민해야 합니다. 가장 간단한 해결책은 결측치가 많고 결과에 영향을 주는 변수가 아닌 triceps와 insulin을 제거하는 것입니다.
pima <- pimaindiansdiabetes="" span="">
pima$obesity[pima$mass>=30]=2
pima$obesity[pima$mass<30 span="">
pima$obesity[pima$mass<25 span="">
table(pima$obesity)
pima<-subset pima="" pressure="">0)
pima<-subset mass="" pima="">0)
pima<-subset glucose="" pima="">0)
out=glm(diabetes~factor(obesity),family=binomial,data=pima)
summary(out)
extractOR(out)
out=glm(diabetes~factor(obesity)+pressure+glucose+age+pedigree,family=binomial,data=pima)
summary(out)
extractOR(out)
> out=glm(diabetes~factor(obesity),family=binomial,data=pima)
> summary(out)
Call:
glm(formula = diabetes ~ factor(obesity), family = binomial,
data = pima)
Deviance Residuals:
Min 1Q Median 3Q Max
-1.094 -1.094 -0.709 1.263 2.306
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -2.5867 0.3919 -6.600 4.12e-11 ***
factor(obesity)1 1.3339 0.4329 3.081 0.00206 **
factor(obesity)2 2.3874 0.4031 5.922 3.18e-09 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 931.94 on 723 degrees of freedom
Residual deviance: 855.40 on 721 degrees of freedom
AIC: 861.4
Number of Fisher Scoring iterations: 5
> extractOR(out)
OR lcl ucl p
(Intercept) 0.08 0.03 0.16 0.0000
factor(obesity)1 3.80 1.62 8.87 0.0021
factor(obesity)2 10.88 4.94 23.99 0.0000
>
> out=glm(diabetes~factor(obesity)+pressure+glucose+age+pedigree,family=binomial,data=pima)
> summary(out)
Call:
glm(formula = diabetes ~ factor(obesity) + pressure + glucose +
age + pedigree, family = binomial, data = pima)
Deviance Residuals:
Min 1Q Median 3Q Max
-2.8133 -0.7259 -0.3803 0.7219 2.5514
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -8.223152 0.831018 -9.895 < 2e-16 ***
factor(obesity)1 1.225678 0.471464 2.600 0.009330 **
factor(obesity)2 2.174498 0.441818 4.922 8.58e-07 ***
pressure -0.001569 0.008389 -0.187 0.851656
glucose 0.034531 0.003597 9.600 < 2e-16 ***
age 0.032288 0.008562 3.771 0.000163 ***
pedigree 0.984527 0.305211 3.226 0.001257 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 931.94 on 723 degrees of freedom
Residual deviance: 676.98 on 717 degrees of freedom
AIC: 690.98
Number of Fisher Scoring iterations: 5
> extractOR(out)
OR lcl ucl p
(Intercept) 0.00 0.00 0.00 0.0000
factor(obesity)1 3.41 1.35 8.58 0.0093
factor(obesity)2 8.80 3.70 20.91 0.0000
pressure 1.00 0.98 1.01 0.8517
glucose 1.04 1.03 1.04 0.0000
age 1.03 1.02 1.05 0.0002
pedigree 2.68 1.47 4.87 0.0013
그 결과 샘플 사이즈가 유지되면서 더 안정적인 결과를 얻었습니다. 이렇게 다중 로지스틱 회귀 분석에서 모델을 결정하는 일이 분석에서 매우 중요합니다. 다음에 이 문제를 가지고 더 이야기해 보겠습니다.
댓글
댓글 쓰기