로지스틱 회귀분석에 대해 질문 올립니다.

SOL 2017.09.15 13:10 Views : 13920

안녕하세요. 교수님.

로지스틱 회귀분석 중 교수님께 여쭙고 싶은 내용이 있어서 질문을 올립니다.

Result of Multiple Logistic Regression
	Estimate	Std. Error	z value	Pr(>\|z\|)	OR	lcl	ucl
(Intercept)	-3.1239	0.3416	-9.15	0.0000	0.04	0.02	0.08
Age	0.0547	0.0112	4.88	0.0000	1.06	1.03	1.08
AT_PN	-1.3352	0.7450	-1.79	0.0731	0.26	0.05	0.98
glm( Recalcitrant ~Age+AT_PN ,data=df,family=binomial)

여기서 종속 변수는 0 또는 1로 코딩되어있으며,

Age 는 연속형 변수, AT_PN 은 0 또는 1로 코딩되어 있는 명목형 변수입니다.

표에서 95% CI 는 0.05-0.98 로 구간이 1 을 포함하지 않고 있는데 Pr(>|z|) value 는 0.073 으로 0.05 보다 큰 것 같습니다.

이 변수 뿐만 아니라 몇몇 변수에서 동일한 결과를 보이고 있는데,

혹시 이 경우에는 결과를 어떻게 해석을 해야 하는지 궁금합니다.

Pr(>|z|) cut-off 를 0.05 가 아닌 다른 수치로 봐야하는 것인지.. 혹은 lcl-ucl 이 95% CI 가 아닌 것인지 궁금합니다.

추가적으로 한가지만 더 답변해주시면 정말 감사드리겠습니다.

위 분석은 "2차 치료까지 필요할 정도로 질환이 recalcitrant 하였는지 여부"를 종속 변수로 둔 것인데,

AT_PN 이라는 finding 여부의 영향을 age 를 covariate 로 하여 adjusting 하고자 합니다.

그렇다면 위 방법처럼 로지스틱 회귀 분석에 AT_PN 및 age 두 변수만 놓고 돌리면 틀린 것은 없는 것인지요?

현재 Stepwise regression 옵션은 풀고 multiple regression 투입 p value 을 1 로 설정하여 두 변수에 의한 회귀 결과가 강제로 나오도록 하여 분석하고 있습니다.

항상 감사드립니다.

I want to

List

Comment 12

스튜던트

2017.09.15 23:10

1.AT_PN 의 회귀계수(-1.3352)와 표준오차(0.7450) 이 맞게 구해졌다는 가정하에, 오즈비의 95% 신뢰구간을 구해보면

exp(-1.3352-1.96*0.7450), exp(-1.3352+1.96*0.7450) -> (0.06,1.13) 으로 나옵니다. 오즈비가 1을 포함하므로

유의한 설명변수가 아니라며, 이는 p-value가 0.05보다 큰 것과 일치 합니다.

아마도 올리주신 표에 오즈비의 신뢰구간이 잘못 기재되었거나, 잘못 계산된것으로 보입니다

2. 네. 위의 회귀식을 해석하면, '나이가 일정할때, AP-PN이 1일때 recalcitrant가 일어날 오즈는 AP-PN이 0일때 비하여 0.26배다'가 됩니다

Comment
cardiomoon

2017.09.15 23:38

웹R에서 OR 계산하는 부분을 다시 검토해보고 글 올리겠습니다.

Comment
천감기

2017.09.15 23:42

1. 번에 대해서는 원자료로 한번 분석해보고 싶어지네요..일반적인 결과가 아니므로 착오 등이 있지는 않은지 궁금합니다.

2. 에 대해서 일반 위의논지는 맞습니다 age를 보정한 APPN의 여향을 보려고 하는.....

그런데 다른 질문이 생기죠?

왜. 성별 다른 변수는 아닌 나이만 분석했느냐?

나이만이 중요하냐...

그러므로 일반적으로 여러 변수를 투입해서 단변수분석부터 하면서 변수들을 제거해가는 방식으로 하죠.

cutㅡoff를 1로했다는 것은 stepwise 등 하지 않고 그냥 enter했다는 뜻인것 같네요. 그렇게 하기도하고 특히 변수가 작을 때 일 때 조금더 유리하죠.

Comment

로지스틱 회귀분석에서 회귀계수의 신뢰구간

로지스틱 회귀분석에서 신뢰구간을 구하는 것은 다음에 따릅니다. 원 자료가 없으므로 survival 패키지의 colon 데이터를 예로 들겠습니다.

require(survival)

Loading required package: survival

colon<-na.omit(colon)
result=glm(status~rx+sex+age+obstruct+perfor,family=binomial,data=colon)
summary(result)


Call:
glm(formula = status ~ rx + sex + age + obstruct + perfor, family = binomial, 
    data = colon)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-1.4485  -1.1949  -0.9553   1.1368   1.4223  

Coefficients:
             Estimate Std. Error z value Pr(>|z|)    
(Intercept)  0.346152   0.262470   1.319    0.187    
rxLev       -0.091277   0.116251  -0.785    0.432    
rxLev+5FU   -0.619551   0.118099  -5.246 1.55e-07 ***
sex         -0.102240   0.096405  -1.061    0.289    
age         -0.002325   0.004052  -0.574    0.566    
obstruct     0.197600   0.122820   1.609    0.108    
perfor       0.367541   0.284735   1.291    0.197    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 2461.7  on 1775  degrees of freedom
Residual deviance: 2423.3  on 1769  degrees of freedom
AIC: 2437.3

Number of Fisher Scoring iterations: 4

오즈비 및 신뢰구간

오즈비 및 신뢰구간은 다음과 같이 구할수 있습니다.

# Odds ratio
exp(coef(result))

(Intercept)       rxLev   rxLev+5FU         sex         age    obstruct 
  1.4136174   0.9127645   0.5381858   0.9028127   0.9976777   1.2184750 
     perfor 
  1.4441789

# confidence Interval 
exp(confint(result))

Waiting for profiling to be done...

                2.5 %    97.5 %
(Intercept) 0.8455328 2.3671160
rxLev       0.7266781 1.1463109
rxLev+5FU   0.4266480 0.6779222
sex         0.7472524 1.0904893
age         0.9897779 1.0056315
obstruct    0.9580619 1.5510076
perfor      0.8311164 2.5534969

웹 R에서 odds ratio의 계산

웹R에서 odds raio의 계산은 위의 공식을 따릅니다. 다음과 같은 함수를 만들어 계산합니다.

require(moonBook)

Loading required package: moonBook

extractOR

function (x, digits = 2) 
{
    suppressMessages(a <- confint(x))
    result = data.frame(exp(coef(x)), exp(a))
    result = round(result, digits)
    result = cbind(result, round(summary(x)$coefficient[, 4], 
        4))
    colnames(result) = c("OR", "lcl", "ucl", "p")
    result
}
<environment: namespace:moonBook>

extractOR(result)

              OR  lcl  ucl      p
(Intercept) 1.41 0.85 2.37 0.1872
rxLev       0.91 0.73 1.15 0.4323
rxLev+5FU   0.54 0.43 0.68 0.0000
sex         0.90 0.75 1.09 0.2889
age         1.00 0.99 1.01 0.5661
obstruct    1.22 0.96 1.55 0.1076
perfor      1.44 0.83 2.55 0.1968

회귀계수로부터 manual 로 계산하기(1)

회귀모형 중 마지막 설명변수인 perfor에 대해 손으로 계산해보겠습니다.

(res=0.3675+c(-1,1)*1.96*0.2847)

[1] -0.190512  0.925512

exp(res)

[1] 0.8265358 2.5231598

손계산이 정확한가?

여기서 회귀계수인 0.3675에서 표준에러의 1.96배를 빼준 것인데 1.96이라는 숫자는 모형의 선형가정이 맞다고 볼때 자료의 갯수가 많아지면 t분포에 따르고 1.96에 수렴하기 때문입니다. sample size가 작을 경우 1.96보다 커집니다. 정확하게는 장도에 따라 달라지는데 이 경우 자유도가 1769이므로 다음과 같이 계산됩니다.

qt(0.975,1769)

[1] 1.961306

즉 1.96보다는 약간 커집니다.

또한 회귀계수를 소수점 4자리에서 반올림하여 계산한 것과 싵제 값은 차이가 있습니다.

회귀계수로부터 manual 로 계산하기(2)

따라서 회귀계수의 계산을 반올림 없이 하면 다음과 같이 됩니다.

res2=summary(result)$coef[7,1]+c(-1,1)*qt(0.975,1769)*summary(result)$coef[7,2]
res2

[1] -0.1909117  0.9259936

exp(res2)

[1] 0.8262055 2.5243751

역시 confint함수로 계산한 것과는 0.03 정도 차이가 나게 됩니다.

샘플 수가 적은 경우(자유도가 적은 경우) confint()함수로 계산한 것과 손으로 반올림해서 계산한 것은 차이가 날수 있을 것 같습니다. 이경우 자유도가 얼마죠?

자료를 올려봐주시면 검토하는데 도움이 되겠습니다. 스튜던트님 천감기님 같이 검토 부탁드립니다.

Comment

SOL

2017.09.16 14:04

파일 업로드 하였습니다. 일부 민감한 정보를 제외한 raw data 를 업로드 하였습니다.

Recalcitrant ~ (Age + ) Type, Recalcitrant ~ (Age + ) AT_PN 에서 동일한 이슈가 발생하였습니다.

많은 선생님들께 감사드립니다.
- Coding_for_upload.csv [File Size:12.6KB/Download:64]
Comment
스튜던트

2017.09.17 00:38

SOL님이 올려주신 95% 신뢰구간을 가지고 역으로 추정해보니,

AT_PN 모델에서는 표준오차에 1.96을 곱한게 아니라 1.76을 곱하였네요..

원래자료는 자유도가 약 270개 정도 되기 때문에

n수가 부족해서 나타나는 현상은 아닌 것 같고,

아마 문교수님이 말씀하신 것 처럼 비선형 회귀분석이기 때문에 발생하는

현상 같아요.

이건 좀 더 찾아봐야 겠습니다. 알게되면 다시 글 올리겠습니다

Comment

일단 올려주신 자료를 R에서 로지스틱회귀분석을 해 보았습니다.

data=read.csv("Downloads/Coding_for_upload.csv",stringsAsFactors = FALSE)
fit=glm(Recalcitrant~Age+AT_PN,family=binomial,data=data)
summary(fit)


Call:
glm(formula = Recalcitrant ~ Age + AT_PN, family = binomial, 
    data = data)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-1.6841  -0.4841  -0.3443  -0.3014   2.7659  

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept) -3.12392    0.34158  -9.145  < 2e-16 ***
Age          0.05467    0.01120   4.882 1.05e-06 ***
AT_PN       -1.33523    0.74500  -1.792   0.0731 .  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 199.26  on 280  degrees of freedom
Residual deviance: 173.00  on 278  degrees of freedom
AIC: 179

Number of Fisher Scoring iterations: 5

회귀계수의 신뢰구간

confint(fit)

Waiting for profiling to be done...

                  2.5 %      97.5 %
(Intercept) -3.84284465 -2.49633775
Age          0.03342183  0.07768483
AT_PN       -2.99254309 -0.01690978

오즈비

confint(fit)

Waiting for profiling to be done...

                  2.5 %      97.5 %
(Intercept) -3.84284465 -2.49633775
Age          0.03342183  0.07768483
AT_PN       -2.99254309 -0.01690978

exp(confint(fit))

Waiting for profiling to be done...

                 2.5 %     97.5 %
(Intercept) 0.02143255 0.08238617
Age         1.03398662 1.08078198
AT_PN       0.05015971 0.98323239

웹 R의 방법

require(moonBook)

Loading required package: moonBook

extractOR(fit)

              OR  lcl  ucl      p
(Intercept) 0.04 0.02 0.08 0.0000
Age         1.06 1.03 1.08 0.0000
AT_PN       0.26 0.05 0.98 0.0731

이상 웹R의 방법이 문제가 없는 것은 확인하였습니다.

여기서 문제는 회귀계수의 신뢰구간이 손으로 계산한것과 다르다는 것입니다. 왜 그럴까요? AT_PN의 경우 회귀계수 -1.33523 이고 표준오차가 0.74500 이므로 다음과 같습니다.
fit$coef[3]+c(-1,1)*qt(0.975,df=278)*summary(fit)$coef[3,2]
[1] -2.801793  0.131329
이 값은 다음과 차이가 많이 납니다. 
confint(fit)
Waiting for profiling to be done...
                  2.5 %      97.5 %
(Intercept) -3.84284465 -2.49633775
Age          0.03342183  0.07768483
AT_PN       -2.99254309 -0.01690978

이 문제는 조금 더 찾아보고 다시 글 올리겠습니다.

Comment

스튜던트

2017.09.17 19:01

답을 찾은 것 같네요. 저도 교수님 처럼 붙여넣기 하고 싶은데

수식이있다 보니 깨져서 나옵니다..

http://rpubs.com/KOSUNGWOO/308829

여기에 가서 보시거나, 첨부 파일을 읽어보시면 됩니다

parameter.html
- parameter.html [File Size:699.4KB/Download:530]
Comment
cardiomoon

2017.09.17 20:56

스튜던트님 잘 보았습니다. 그러면 가능도비로 p값을 구하려면 어떻게 해야 하나요?

Comment
스튜던트

2017.09.17 22:31

model2<-glm(Recalcitrant~Age,data = sov,family = binomial)

# model1과 비교하기위해 AT_PN변수를 뺀 model2를 적합하였음

library(lmtest)
lrtest(model2,model1)

이렇게 구하거나,

anova(model2,model1,test = "LRT")

이렇게 하면 AT_PN 계수에 대한 가능도비 검정의 p-value를 구할 수 있습니다.

만족할만한 대답이 되었는지는 모르겠습니다만 ㅠㅠ

Comment

                Variables in the Equation
                                95.0% C.I.for EXP(B)
       B    S.E.    Wald    df    Sig.    Exp(B)    Lower    Upper
Step 1a Age    .055    .011    23.829    1    .000    1.056    1.033    1.080
   AT_PN    -1.335    .745    3.212    1    .073    .263    .061    1.133
Constant    -3.124    .342    83.636    1    .000    .044
a. Variable(s) entered on step 1: Age, AT_PN.

요건 SPSS에서 구한 것입니다.

일단 R에서 구한 것과는 다른 값을 제시하고 있고

엑셀과 SPSS는 같은 값을 보여 주네요.

Comment

SOL

2017.09.20 11:26

제가 모르는 별천지의 이야기라서 이해를 하지는 못하였지만 많이 찾아봐주시고 도움 주신 선생님들께 감사드립니다.

Comment

No.	Subject	Author	Date
1174	Web R studio 초기화 후 문의드립니다. [2]	화모기	2016.02.23
1173	[정보] 구조방정식 모형 분석 관련 참고자료 소개 [2]	보석	2016.02.24
1172	그래프가 이상하게 나와서 어떻게 하면 좋을까요? [3]	osr33260216	2016.02.25
1171	데이터 업로드가 되지 않습니다.. [3]	jenni	2016.02.28
1170	어제 강의 너무 잘 들었습니다 [5]	마술샘	2016.02.29
1169	lavaan의 공식 구글 토론실입니다	마술샘	2016.02.29
1168	[오류]web-r 구조방정식 [2]	마술샘	2016.02.29
1167	저 이 홈피 체험 1주일 지난것같은데... [3]	osr33260216	2016.03.03
1166	이럴 땐 어떤 통계 기법을 쓰나요? (기초) [1]	송이버섯	2016.03.03
1165	남북한별 인구를 나타내는 막대 그래프를 만들려는데 왜 바가 않나타나죠? [1]	osr33260216	2016.03.04
1164	제가 금시세 표를 만들었는데 왜 업로드 파일은 이상하게 나왔나요? [1]	osr33260216	2016.03.04
1163	Standard error 질문 드립니다. [7]	R까기	2016.03.05
1162	object 'res' not found [3]	민자	2016.03.06
1161	제가 선생님 말씀대로 금시세표를 csv쉼표분리로 변경했는데 이상하게되었습니다. [1]	osr33260216	2016.03.07
1160	인공지능, 인지 컴퓨팅, 머신 러닝, 딥 러닝 관련 IBM 자료 [1]	Ben	2016.03.09
1159	[RStudio Server] 사용 연장 문의 [2]	KB	2016.03.09
1158	descriptive 표만들기 질문입니다. [4]	난둥	2016.03.10
1157	그래프 작성 문의드립니다. [4]	도치맘	2016.03.10
1156	표만들기 관련 질문입니다. [1]	tympa	2016.03.10
1155	머신러닝 레볼루션 행사안내 [1]	mike	2016.03.10

Write Tag

First Page 1 2 3 4 5 6 7 8 9 10 Last Page

	coeff b	s.e.	Wald	p-value	exp(b)	lower	upper
Intercept	-3.12392	0.341588	83.63607	5.95E-20	0.043984
AT_PN	-1.33523	0.745038	3.211856	0.073107	0.263097	0.061086	1.133167
Age	0.054675	0.011201	23.82864	1.05E-06	1.056197	1.033264	1.07964