R语言之Logistic回归

发表于   |   分类于   |  

背景

Logistic回归通常情况下会用在分类问题中，一般是二分类问题。Logistic回归直接将预测限定在（0， 1）区间内，因此也可以将它理解为线性回归的归一化。

数据

数据集为一个医学数据，主要用来预测新生儿是否需要额外的医学护理。数据在这里

R语言实现过程

加载数据

首先加载数据集，并划分训练集和测试集。

load("D:/zmPDSwR-master/CDC/NatalRiskData.rData")
train <- sdata[sdata$ORIGRANDGROUP <= 5,]
test <- sdata[sdata$ORIGRANDGROUP > 5,]

构建模型

由于数据集中变量很多，先构建一个公式来方便建模：

complications <- c("ULD_MECO", "ULD_PRECIP", "ULD_BREECH")
riskfactors <- c("URF_DIAB", "URF_CHYPER", "URF_PHYPER", "URF_ECLAM")
y <- "atRisk"
x <- c("PWGT",
    "UPREVIS",
    "CIG_REC",
    "GESTREC3",
    "DPLURAL",
    complications,
    riskfactors)
fmla <- paste(y, paste(x, collapse = "+"), sep = "~")

> print(fmla)
[1] "atRisk~PWGT+UPREVIS+CIG_REC+GESTREC3+DPLURAL+ULD_MECO+ULD_PRECIP+ULD_BREECH+URF_DIAB+URF_CHYPER+URF_PHYPER+URF_ECLAM"

建模并预测：

model <- glm(fmla, data = train, family = binomial(link = "logit"))

train$pred <- predict(model, newdata = train, type = "response")
test$pred <- predict(model, newdata = test, type = "response")

评价模型

先通过图形来看一下打分情况：

library(ggplot2)
ggplot(train, aes(x = pred, color = atRisk, linetype = atRisk)) + geom_density() + scale_colour_hue()

id

打分的分布并不是分离的，发现模型的效果不是很好。这个时候需要选择一个阈值，高于该阈值的打分为正例，反之为负例。在选择阈值时，主要是在准确率和召回率之间寻找平衡。下面是用ROC包来做这件事：

library(ROCR)
library(grid)
predObj <- prediction(train$pred, train$atRisk)
precObj <- performance(predObj, measure = "prec")
recObj <- performance(predObj, measure = "rec")

precision <- (precObj@y.values)[[1]]
prec.x <- (precObj@x.values)[[1]]
recall <- (recObj@y.values)[[1]]

rocFrame <- data.frame(threshold = prec.x, precision = precision, recall = recall)

nplot <- function(plist) {
    n <- length(plist)
    grid.newpage()
    pushViewport(viewport(layout = grid.layout(n, 1)))
    vplayout = function(x, y) {
        viewport(layout.pos.row = x, layout.pos.col = y)
    }

    for(i in 1:n) {
        print(plist[[i]], vp = vplayout(i, 1))
    }
}

pnull <- mean(as.numeric(train$atRisk))

p1 <- ggplot(rocFrame, aes(x = threshold)) + geom_line(aes(y = precision / pnull)) + coord_cartesian(xlim = c(0, 0.05), ylim = c(0, 10))
p2 <- ggplot(rocFrame, aes(x = threshold)) + geom_line(aes(y = recall)) + coord_cartesian(xlim = c(0, 0.05))

nplot(list(p1, p2))

id

选择一个阈值0.02来看看结果：

> ctab.test <- table(pred = test$pred > 0.02, atRisk = test$atRisk)
> ctab.test
    atRisk
pred    FALSE TRUE
FALSE  9487   93
TRUE   2405  116
> precision <- ctab.test[2, 2] / sum(ctab.test[2, ])
> precision
[1] 0.04601349
> recall <- ctab.test[2, 2] / sum(ctab.test[, 2])
> recall
[1] 0.5550239
> enrich <- precision / mean(as.numeric(test$atRisk))
> enrich
[1] 2.664159

再来看看模型的概要：

> summary(model)

Call:
glm(formula = fmla, family = binomial(link = "logit"), data = train)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-0.9732  -0.1818  -0.1511  -0.1358   3.2641  

Coefficients:
                        Estimate Std. Error z value Pr(>|z|)    
(Intercept)              -4.412189   0.289352 -15.249  < 2e-16 ***
PWGT                      0.003762   0.001487   2.530 0.011417 *  
UPREVIS                  -0.063289   0.015252  -4.150 3.33e-05 ***
CIG_RECTRUE               0.313169   0.187230   1.673 0.094398 .  
GESTREC3< 37 weeks        1.545183   0.140795  10.975  < 2e-16 ***
DPLURALtriplet or higher  1.394193   0.498866   2.795 0.005194 ** 
DPLURALtwin               0.312319   0.241088   1.295 0.195163    
ULD_MECOTRUE              0.818426   0.235798   3.471 0.000519 ***
ULD_PRECIPTRUE            0.191720   0.357680   0.536 0.591951    
ULD_BREECHTRUE            0.749237   0.178129   4.206 2.60e-05 ***
URF_DIABTRUE             -0.346467   0.287514  -1.205 0.228187    
URF_CHYPERTRUE            0.560025   0.389678   1.437 0.150676    
URF_PHYPERTRUE            0.161599   0.250003   0.646 0.518029    
URF_ECLAMTRUE             0.498064   0.776948   0.641 0.521489    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 2698.7  on 14211  degrees of freedom
Residual deviance: 2463.0  on 14198  degrees of freedom
AIC: 2491

Number of Fisher Scoring iterations: 7

要点总结

  • Logistic回归是二元分类首选的统计建模方法。首先尝试Logistic回归，如果适用，再尝试其他方法；
  • 处理具有大量变量或者大量类别值得类别型变量的问题时，Logistic回归将遇到麻烦；
  • Logistic可以被很好地校准：它再现数据的边缘概率；
  • 即便在变量相关的情况下，Logistic回归也能很好地进行预测，但相关的变量会降低模型给出的建议质量；
  • 特别大的系数量级、系数估计上得到特别大的标准差、系数上的错误符号，都可能是输入变量具有相关性的征兆；
  • 太多的菲舍尔迭代次数，或者具有特别大标准差的特大系数，可能标志着一个输入变量或输入变量组合与所要应答的一个子集完全相关。这可能需要通过分割数据来处理这个问题；
  • glm()能提供好的诊断，但利用测试数据对模型进行检查仍然是最有效的诊断方法；
  • 伪R-平方是关于拟合优度的一个有用的启发式度量指标。

参考文献

[1][数据科学 理论、方法与R语言实践]