在日常學習中會經常使用線性回歸模型對一些數據進行預測：房租價格，二手房價格，保溫杯評分的影響因素（R線性回歸（一））等，這些響應變數都屬於連續型變數。在許多情況下，假設因變數為正態分布並不合理，例如結果變數是類別型的，是/否，通過/未通過，差/良好/優秀這些都不是正態分布。這時候就需要用logistic回歸。

首先捋一下邏輯回歸、線性回歸、廣義線性模型這其中的關係。

回歸分析模型中包括線性回歸和邏輯回歸，不同的模型是基於不同類型的因變數，廣義線性模型是為了克服線性回歸模型的缺點出現的，是線性回歸模型的推廣。

回歸分析中的各種變體

  • 簡單線性：用一個量化的解釋變數預測一個量化的響應變數
  • 多項式：用一個量化的解釋變數預測一個量化的響應變數，模型的關係式N階多項式
  • 多元線性：用兩個或多個量化的解釋變數預測一個量化的響應變數
  • 多變數：用一個或多個解釋變數預測多個響應變數
  • Logistic回歸：用一個或多個解釋變數預測一個類別型響應變數
  • 泊松回歸：用一個或多個解釋變數預測一個代表頻數的響應變數
  • 時間序列回歸：對誤差項想關的時間序列數據建模

廣義線性模型中，自變數可以是離散的，也可以是連續的。

與線性回歸模型相比較，有以下推廣：

（1）隨機誤差項不一定服從正態分布，可以服從二項、泊松、負二項、正態、伽馬、逆高斯等分布，這些分布被統稱為指數分布族。

（2）引入聯接函數。因變數和自變數通過聯結函數產生影響，聯接函數滿足單調，可導。

根據不同的數據，可以自由選擇不同的模型。大家比較熟悉的Logit模型就是使用Logit聯接、隨機誤差項服從二項分布得到的模型。

廣義線性模型中模型形式基本相似，不同的就是因變數：

  • 如果是連續的，就是多重線性回歸
  • 如果是二項分布，就是Logistic回歸
  • 如果是Poisson分布，就是poisson回歸
  • 如果是負二項分布，就是負二項回歸

Logistic回歸的因變數可以是二分類的，也可以是多分類的，但是二分類的更為常用，也更容易理解。所以實際中最常用的就是二分類的logistic回歸。

一、logistic回歸模型概述

廣義線性回歸是探索「響應變數的期望」與「自變數」的關係，以實現對非線性關係的某種擬合。這裡有一個「連接函數」和一個「誤差函數」，「響應變數的期望」經過連接函數作用後，與「自變數」存在線性關係。選取不同的「連接函數」與「誤差函數」可以構造不同的廣義回歸模型。當誤差函數取「二項分布」而連接函數取「logit函數」時，就是常見的「logistic回歸模型」,在0-1響應的問題中的得到了大量的應用。（摘抄自：如何在R語言中使用Logistic回歸模型 - Little_Rookie - 博客園）

線性回歸的公式如下：

z= heta_{0}+ heta_{1}x_{1}+ heta_{2}x_{2}+...... heta_{n}x_{n}= heta^{T}x

對於Logistic Regression來說，其思想也是基於線性回歸（Logistic Regression屬於廣義線性回歸模型），廣義線性回歸模型公式如下：

h{ heta}(x)=frac{1}{1+e^{-z}}=frac{1}{1+e^{- heta^{T}x}}

其中

y=frac{1}{1+e^{-x}}

被稱作sigmoid函數，Logistic Regression演算法是將線性函數的結果映射到了sigmoid函數中。

sigmoid的函數圖形如下：

可以看到，sigmoid的函數輸出是介於（0,1）之間的，中間值是0.5，因為 h_{ heta}(x) 輸出是介於（0,1）之間，也就表明了數據屬於某一類別的概率。例如：

h_{ heta}(x) <0.5，則說明當前數據屬於A類；

h_{ heta}(x) >0.5，則說明當前數據屬於B類。

所以我們可以將sigmoid函數看成樣本數據的概率密度函數。

二、Logistic回歸實例

數據來源《預測分析R語言實現》中邏輯回歸部分。數據傳送門：Statlog (Heart) Data Se。

該數據包含了270個可能心臟有問題的病人的觀測數據。其中120個病人表現出心臟有問題。我們根據病人的基本情況和一系列的醫學檢查來預測病人是否患有心臟病。

數據中的輸入和輸出特徵：

data<-read.table("E://Data For R/RData/heart.txt")
head(data)
names(data)<-c("AGE","SEX","CHESTPAIN","RESTBP","CHOL","SUGAR","ECG","MAXHR","ANGINA","DEP","EXERCISE",
"FLUOR","THAL","OUTPUT")
head(data)
####我們在訓練數據之前都需要對數據進行預處理
#觀察數據是否有缺失值
colSums(is.na(data))
AGE SEX CHESTPAIN RESTBP CHOL SUGAR ECG
0 0 0 0 0 0 0
MAXHR ANGINA DEP EXERCISE FLUOR THAL OUTPUT
0 0 0 0 0 0 0
##觀察是否具有缺失值時，我們也可以通過可視化的手段去判斷。
install.packsges(Amelia)
library(Amelia)
missmap(data,main="Missing values vs observed")##見圖一
str(data)
#在處理數值型數據時常遇見的一個問題，某個特徵實際上是分類變數而非數值變數。本列中CHESTPAIN
ECG、THAL、EXERCISE特徵都是分類特徵。對於SEX、SUGAR、ANGINA都是二元類別，這裡不需要處理。
data$CHESTPAIN<-as.factor(data$CHESTPAIN)
data$ECG<-as.factor(data$ECG)
data$THAL<-as.factor(data$THAL)
data$EXERCISE<-as.factor(data$EXERCISE)
#對於輸出變數OUTPUT變數，類別1對應沒有心臟病，類別2對應患有心臟病。我們將其改成熟悉的0-1
data$OUTPUT<-data$OUTPUT-1
###需要將數據分成訓練集和測試集，採取隨機抽樣的方法，訓練集：測試集=4:1
set.seed(125424)
aa<-sample(1:nrow(data),size=round(4/5*nrow(data),0))
train_data<-data[aa,]
test_data<-data[-aa,]
###數據準備完成
heart_model<-glm(OUTPUT~.,data=data,family=binomial("logit"))
summary(heart_model)
Call:
glm(formula = OUTPUT ~ ., family = binomial("logit"), data = train_data)

Deviance Residuals:
Min 1Q Median 3Q Max
-2.76261 -0.42245 -0.09795 0.29330 2.50841

Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -8.207e+00 3.801e+00 -2.159 0.03083 *
AGE -2.175e-02 3.102e-02 -0.701 0.48325
SEX 1.929e+00 7.310e-01 2.638 0.00833 **
CHESTPAIN2 2.251e+00 1.171e+00 1.922 0.05463 .
CHESTPAIN3 1.577e+00 1.003e+00 1.572 0.11600
CHESTPAIN4 2.866e+00 9.815e-01 2.920 0.00350 **
RESTBP 3.308e-02 1.478e-02 2.238 0.02522 *
CHOL 4.663e-03 5.010e-03 0.931 0.35198
SUGAR -1.373e+00 8.955e-01 -1.533 0.12528
ECG1 -1.228e+01 1.455e+03 -0.008 0.99327
ECG2 3.519e-01 4.891e-01 0.719 0.47188
MAXHR -1.966e-02 1.262e-02 -1.558 0.11926
ANGINA 9.421e-01 5.395e-01 1.746 0.08077 .
DEP 5.425e-01 2.861e-01 1.897 0.05788 .
EXERCISE2 1.132e+00 6.123e-01 1.849 0.06442 .
EXERCISE3 3.973e-01 1.083e+00 0.367 0.71375
FLUOR 1.485e+00 3.415e-01 4.348 1.37e-05 ***
THAL6 3.079e-02 9.074e-01 0.034 0.97293
THAL7 1.535e+00 5.282e-01 2.906 0.00366 **
---
Signif. codes: 0 『***』 0.001 『**』 0.01 『*』 0.05 『.』 0.1 『 』 1

(Dispersion parameter for binomial family taken to be 1)

Null deviance: 297.19 on 215 degrees of freedom
Residual deviance: 124.84 on 197 degrees of freedom
AIC: 162.84

Number of Fisher Scoring iterations: 14
#我們根據摘要，看到FLOUR、CHESTPAIN4、THAL7對於心臟病是最強的特徵預測因子。很多輸入特徵有較高的
p值。這表明在有其他特徵存在的情況下，他們可能不是很好的心臟病指標。一般情況下，我們認為心臟病病發的可能性
應該隨著年齡的增長而增長，但是年齡的回歸係數是負數，所以我們肯定在特徵里存在某種程度的共線性。
heart_model2<-glm(OUTPUT~AGE,data=train_data,family=binomial("logit"))
summary(heart_model2)
Call:
glm(formula = OUTPUT ~ AGE, family = binomial("logit"), data = train_data)

Deviance Residuals:
Min 1Q Median 3Q Max
-1.5757 -1.0612 -0.8208 1.1695 1.7225

Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -3.04165 0.88451 -3.439 0.000584 ***
AGE 0.05187 0.01589 3.264 0.001097 **
---
Signif. codes: 0 『***』 0.001 『**』 0.01 『*』 0.05 『.』 0.1 『 』 1

(Dispersion parameter for binomial family taken to be 1)

Null deviance: 297.19 on 215 degrees of freedom
Residual deviance: 285.80 on 214 degrees of freedom
AIC: 289.8

Number of Fisher Scoring iterations: 4
#重新訓練一個只有AGE變數的邏輯回歸模型，我們會得到一個正的回歸係數和一個較低的p值，這兩者
都支持我們對於存在共線性的判斷。

##運用predict()函數來計算模型的輸出，這個輸出就是輸入屬於類別1的概率。通過運用一個閾值來進行
二元分類，通過針對訓練和測試數據進行預測，並運用我們的預期輸出與之比較，來衡量其分類的精確度。
train_predictions<-predict(heart_model,newdata=train_data,type="response")
train_class_predictions<-as.numeric(train_predictions>0.5)
mean(train_class_predictions==train_data$OUTPUT)
[1] 0.8842593
test_predictions<-predict(heart_model,newdata=test_data,type="response")
test_class_predictions<-as.numeric(test_predictions>0.5)
mean(test_class_predictions==test_data$OUTPUT)
[1] 0.8518519
##訓練集合測試集上的分類精確度很相似，接近90%。模型尚可。

圖一：缺失數據可視化

參考：

邏輯回歸 - CSDN博客

淺談線性、非線性和廣義線性回歸模型 - 火星十一郎 - 博客園

如何在R語言中使用Logistic回歸模型 - Little_Rookie - 博客園

Logistic Regression（邏輯回歸）原理及公式推導 - CSDN博客

R機器學習之二：邏輯回歸 - CSDN博客（ROC曲線解釋比較詳細）


推薦閱讀：
相关文章