本内容为【科研私家菜】R语言机器学习与临床预测模型系列课程

你想要的R语言学习资料都在这里， 快来收藏关注【科研私家菜】

01 Logistic回归

逻辑斯蒂回归就是以对数发生比为响应变量进行线性拟合，即log(P(Y)/1 - P(Y)) =B0+B1x。 这里的系数是通过极大似然估计得到的，而不是通过OLS。极大似然的直观意义就是，我们要找到一对B0和B1的估计值，使它们产生
的对观测的预测概率尽可能接近Y的实际观测结果，这就是所谓的似然性。
R中的glm()函数可以拟合广义线性模型，这是一系列模型，其中包括逻辑斯蒂回归模型。其中一个重点是必须在函数中使用family = binomial这个参数。该参数告诉R运行逻辑斯蒂回归。

library(MASS)
data(biopsy)
str(biopsy)
biopsy$ID = NULL
names(biopsy) = c("thick", "u.size", "u.shape", "adhsn", 
                  "s.size", "nucl", "chrom", "n.nuc", "mit", "class")
names(biopsy)
biopsy.v2 <- na.omit(biopsy)
y <- ifelse(biopsy.v2$class == "malignant", 1, 0)
library(reshape2)
library(ggplot2)
biop.m <- melt(biopsy.v2, id.var = "class")
ggplot(data = biop.m, aes(x = class, y = value)) + 
  geom_boxplot() +
  facet_wrap(~ variable, ncol = 3)
library(corrplot)
bc <- cor(biopsy.v2[ ,1:9]) #create an object of the features
corrplot.mixed(bc)

set.seed(123) #random number generator
ind <- sample(2, nrow(biopsy.v2), replace = TRUE, prob = c(0.7, 0.3))
train <- biopsy.v2[ind==1, ] #the training data set
test <- biopsy.v2[ind==2, ] #the test data set
str(test) #confirm it worked
table(train$class)
table(test$class)

full.fit <- glm(class ~ ., family = binomial, data = train)
summary(full.fit)
confint(full.fit)
exp(coef(full.fit))
library(car)
vif(full.fit)

train.probs <- predict(full.fit, type = "response")
train.probs[1:5] #inspect the first 5 predicted probabilities
contrasts(train$class)

library(InformationValue)
trainY <- y[ind==1]
testY <- y[ind==2]
confusionMatrix(trainY, train.probs)
# optimalCutoff(trainY, train.probs)
misClassError(trainY, train.probs)
confusionMatrix(trainY, train.probs)

test.probs <- predict(full.fit, newdata = test, type = "response")
misClassError(testY, test.probs)
confusionMatrix(testY, test.probs)

通过summary()函数，我们可以查看各个预测变量的系数及其p值。可以看到，只有两个特征的p值小于0.05（ thickness和nuclei）。使用confint()函数可以对模型进行95%置信区间的检验。
对于逻辑斯蒂回归模型的系数，你不能解释为“当X改变1个单位时Y会改变多少”。这时，优势比的作用就体现出来了。对数函数log(P(Y)/1 - P(Y)) =
B0 + B1x的B系数可以通过exponent(beta)指数函数转化为优势比。
可以使用下面的exp(coef())函数形式：
exp(coef(full.fit))
优势比可以解释为特征中1个单位的变化导致的结果发生比的变化。如果系数大于1，就说明当特征的值增加时，结果的发生比会增加。反之，系数小于1就说明，当特征的值增加时，结果的发生比会减小。在本例中，除u.size之外的所有特征都会增加对数发生比。
同线性回归一样，潜在的多重共线性在逻辑斯蒂回归中也可以算出VIF统计量
library(car)
vif(full.fit)
下一步需要评价模型在训练集上执行的效果，然后再评价它在测试集上的拟合程度。快速实现评价的方法是生成一个混淆矩阵。在后面的章节中，我们使用的混淆矩阵是由caret包实现的，InformationValue包也可以实现混淆矩阵。

效果如下：

02 使用交叉验证的逻辑回归

交叉验证的目的是提高测试集上的预测正确率，以及尽可能避免过拟合。 K折交叉验证的做法是将数据集分成K个相等的等份，每个等份称为一个K子集（ K-set）。算法每次留出一个子集，使用其余K -1个子集拟合模型，然后用模型在留出的那个子集上做预测。将上面K次验证的结果进行平均，可以使误差最小化，并且获得合适的特征选择。你也可以使用留一交叉验证方法，这里的K等于N。模拟表明， LOOCV可以获得近乎无偏的估计，但是会有很高的方差。所以，大多数机器学习专家都建议将K的值定为5或10。
bestglm包可以自动进行交叉验证，这个包依赖于我们在线性回归中使用过的leaps包。

library(bestglm)
X <- train[, 1:9]
Xy <- data.frame(cbind(X, trainY))
bestglm(Xy = Xy, IC = "CV", CVArgs = list(Method = "HTF", K = 10, REP = 1), 
        family=binomial)
reduce.fit <- glm(class ~ thick + u.size + nucl, family = binomial, data = train)

test.cv.probs = predict(reduce.fit, newdata = test, type = "response")
misClassError(testY, test.cv.probs)
confusionMatrix(testY, test.cv.probs)

bestglm(Xy = Xy, IC = "BIC", family = binomial)
bic.fit <- glm(class ~ thick + adhsn + nucl + n.nuc, 
               family = binomial, data = train)
test.bic.probs = predict(bic.fit, newdata = test, type = "response")
misClassError(testY, test.bic.probs)
confusionMatrix(testY, test.bic.probs)

关注R小盐，关注科研私家菜(VX_GZH: SciPrivate)，有问题请联系R小盐。让我们一起来学习 R语言机器学习与临床预测模型