本內(nèi)容為【科研私家菜】R語言機(jī)器學(xué)習(xí)與臨床預(yù)測(cè)模型系列課程

你想要的R語言學(xué)習(xí)資料都在這里种远， 快來收藏關(guān)注【科研私家菜】

01 Logistic回歸

邏輯斯蒂回歸就是以對(duì)數(shù)發(fā)生比為響應(yīng)變量進(jìn)行線性擬合卵迂，即log(P(Y)/1 - P(Y)) =B0+B1x认臊。 這里的系數(shù)是通過極大似然估計(jì)得到的愧膀，而不是通過OLS馆类。極大似然的直觀意義就是歼培，我們要找到一對(duì)B0和B1的估計(jì)值亩码，使它們產(chǎn)生
的對(duì)觀測(cè)的預(yù)測(cè)概率盡可能接近Y的實(shí)際觀測(cè)結(jié)果季率，這就是所謂的似然性。
R中的glm()函數(shù)可以擬合廣義線性模型描沟，這是一系列模型飒泻，其中包括邏輯斯蒂回歸模型。其中一個(gè)重點(diǎn)是必須在函數(shù)中使用family = binomial這個(gè)參數(shù)吏廉。該參數(shù)告訴R運(yùn)行邏輯斯蒂回歸泞遗。

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()函數(shù)，我們可以查看各個(gè)預(yù)測(cè)變量的系數(shù)及其p值席覆∈氛蓿可以看到，只有兩個(gè)特征的p值小于0.05（ thickness和nuclei）佩伤。使用confint()函數(shù)可以對(duì)模型進(jìn)行95%置信區(qū)間的檢驗(yàn)聊倔。
對(duì)于邏輯斯蒂回歸模型的系數(shù)，你不能解釋為“當(dāng)X改變1個(gè)單位時(shí)Y會(huì)改變多少”生巡。這時(shí)耙蔑，優(yōu)勢(shì)比的作用就體現(xiàn)出來了。對(duì)數(shù)函數(shù)log(P(Y)/1 - P(Y)) =
B0 + B1x的B系數(shù)可以通過exponent(beta)指數(shù)函數(shù)轉(zhuǎn)化為優(yōu)勢(shì)比孤荣。
可以使用下面的exp(coef())函數(shù)形式：
exp(coef(full.fit))
優(yōu)勢(shì)比可以解釋為特征中1個(gè)單位的變化導(dǎo)致的結(jié)果發(fā)生比的變化甸陌。如果系數(shù)大于1须揣，就說明當(dāng)特征的值增加時(shí)，結(jié)果的發(fā)生比會(huì)增加邀层。反之返敬，系數(shù)小于1就說明，當(dāng)特征的值增加時(shí)寥院，結(jié)果的發(fā)生比會(huì)減小劲赠。在本例中，除u.size之外的所有特征都會(huì)增加對(duì)數(shù)發(fā)生比秸谢。
同線性回歸一樣凛澎，潛在的多重共線性在邏輯斯蒂回歸中也可以算出VIF統(tǒng)計(jì)量
library(car)
vif(full.fit)
下一步需要評(píng)價(jià)模型在訓(xùn)練集上執(zhí)行的效果，然后再評(píng)價(jià)它在測(cè)試集上的擬合程度估蹄∷芗澹快速實(shí)現(xiàn)評(píng)價(jià)的方法是生成一個(gè)混淆矩陣。在后面的章節(jié)中臭蚁，我們使用的混淆矩陣是由caret包實(shí)現(xiàn)的最铁，InformationValue包也可以實(shí)現(xiàn)混淆矩陣。

效果如下：

02 使用交叉驗(yàn)證的邏輯回歸

交叉驗(yàn)證的目的是提高測(cè)試集上的預(yù)測(cè)正確率垮兑，以及盡可能避免過擬合冷尉。 K折交叉驗(yàn)證的做法是將數(shù)據(jù)集分成K個(gè)相等的等份，每個(gè)等份稱為一個(gè)K子集（ K-set）系枪。算法每次留出一個(gè)子集雀哨，使用其余K -1個(gè)子集擬合模型，然后用模型在留出的那個(gè)子集上做預(yù)測(cè)私爷。將上面K次驗(yàn)證的結(jié)果進(jìn)行平均绽昼，可以使誤差最小化贷币，并且獲得合適的特征選擇截驮。你也可以使用留一交叉驗(yàn)證方法啄糙，這里的K等于N。模擬表明工秩， LOOCV可以獲得近乎無偏的估計(jì)嘉栓，但是會(huì)有很高的方差。所以拓诸，大多數(shù)機(jī)器學(xué)習(xí)專家都建議將K的值定為5或10侵佃。
bestglm包可以自動(dòng)進(jìn)行交叉驗(yàn)證，這個(gè)包依賴于我們?cè)诰€性回歸中使用過的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)

關(guān)注R小鹽馋辈，關(guān)注科研私家菜(VX_GZH: SciPrivate)，有問題請(qǐng)聯(lián)系R小鹽倍谜。讓我們一起來學(xué)習(xí) R語言機(jī)器學(xué)習(xí)與臨床預(yù)測(cè)模型