앤디 필드의 유쾌한 R 통계학 챕터 8 - 로지스틱 회귀분석

독립변수와 종속변수가 선형 관계가 아닐 때 로지스틱 회귀분석 사용을 고려한다.

범주형 데이터인지 확인해서 as.factor로 적절하게 변경해야 한다.

 

library(car);library(mlogit)

eelData <- read.delim('datafiles/eel.dat', header = TRUE)
summary(eelData)

# 로지스틱 회귀분석에 사용하기 위해 factor로 변경한다.
eelData$Cured <- as.factor(eelData$Cured)
eelData$Intervention <- as.factor(eelData$Intervention)

# 보기 쉽게 하기 위해 기저 범주를 변경한다.
eelData$Cured <- relevel(eelData$Cured, 'Not Cured')
eelData$Intervention <- relevel(eelData$Intervention, 'No Treatment')

# 로지스틱 회귀 분석은 glm 함수를 활용한다.
eelModel.1 <- glm(Cured ~ Intervention, eelData, family = binomial())
eelModel.2 <- glm(Cured ~ Intervention + Duration, eelData, family = binomial())
summary(eelModel.1) # Null deviance는 기저 모형의 이탈도, Residual deviance는 모형의 이탈도이다.
summary(eelModel.2) # 이탈도가 낮을 수록 모델은 적합하다.

# 모형의 유의성을 검증하기 위해 카이제곱을 계산한다
modelChi <- eelModel.1$null.deviance - eelModel.1$deviance # 이탈도 계산
chidf <- eelModel.1$df.null - eelModel.1$df.residual # 자유도 계산
chisq.prob <- 1 - pchisq(modelChi, chidf) # 카이제곱 통계량 계산

# 오즈비와 신뢰구간 계산하기. 1보다 크면 예측변수가 증가함에 따라 결과 변수도 1로 처리될 승산이 높은 것이다.
exp(eelModel.1$coefficients)
exp(confint(eelModel.1)) # 신뢰구간에 1이 포함되면 신뢰할 수 없다

# 두 모델 비교하기(anova 함수로)
anova(eelModel.1, eelModel.2)

# 잔차 진단을 위한 통계량
eelData$predicted.probabilities <- fitted(eelModel.1)
eelData$standardized.residuals <- rstandard(eelModel.1)
eelData$studentized.residuals <- rstudent(eelModel.1)
eelData$dfbeta <- dfbeta(eelModel.1)
eelData$dffit <- dffits(eelModel.1)
eelData$leverage <- hatvalues(eelModel.1)
head(eelData[, c('Cured', 'Intervention', 'Duration', 'predicted.probabilities')])
eelData[, c('leverage', 'studentized.residuals', 'dfbeta')]

# 다중공선성 진단
penaltyData <- read.delim('datafiles/penalty.dat', header = TRUE)
penaltyData$Scored <- as.factor(penaltyData$Scored)
penaltyModel.2 <- glm(Scored ~ Previous + PSWQ + Anxious, penaltyData, family = binomial())
vif(penaltyModel.2)

# 선형성 검증하기
penaltyData$logPSWQInt <- log(penaltyData$PSWQ)*penaltyData$PSWQ
penaltyData$logAnxiousInt <- log(penaltyData$Anxious)*penaltyData$Anxious
penaltyData$logPreviousInt <- log(penaltyData$Previous)*penaltyData$Previous
head(penaltyData)

penaltyTest.1 <- glm(Scored ~ Previous + PSWQ + Anxious + logPSWQInt + logAnxiousInt
                     + logPreviousInt, penaltyData, family = binomial())
summary(penaltyTest.1)

# 다중 로지스틱 회귀분석
chatData <- read.delim('datafiles/Chat-up Lines.dat', header = TRUE)
chatData$Success <- as.factor(chatData$Success)
chatData$Gender <- as.factor(chatData$Gender)
chatData$Gender <- relevel(chatData$Gender, ref = 2)

# 다중 로지스틱 회귀 분석에 맞게 형식 변경하기
mlChat <- mlogit.data(chatData, choice = 'Success', shape = 'wide')
head(mlChat)

# 다중 로지스틱 회귀분석 실행
chatModel <- mlogit(Success ~ 1 | Good_Mate + Funny + Gender + Sex + Gender:Sex + Funny:Gender, mlChat)
summary(chatModel)
data.frame(exp(chatModel$coefficients))
exp(confint(chatModel))