一、定義

   logistic迴歸是一種廣義線性迴歸（generalized linear model），因此與多重線性迴歸分析有很多相同之處。它們的模型形式基本上相同，都具有 w'x+b，其中w和b是待求引數，其區別在於他們的因變數不同，多重線性迴歸直接將w'x+b作為因變數，即y =w'x+b，而logistic迴歸則通過函式L將w'x+b對應一個隱狀態p，p =L(w'x+b),然後根據p 與1-p的大小決定因變數的值。邏輯迴歸假設資料服從伯努利分佈，通過極大似然估計的方法，運用梯度下降來求解引數，來達到資料二分類的目的。

• 迴歸：假設現在有一些資料點，我們用一條直線對這些點進行擬合（該直線稱為最佳擬合直線），這個擬合的過程就稱為迴歸
  。
• 主要思想：根據現有資料對分類邊界線建立迴歸公式，以此進行分類。
• 邏輯迴歸的一般過程：  

  1. 收集資料：採用任意方法收集資料。

  2. 準備資料：由於需要進行距離計算，因此要求資料型別為數值型。另外，結構化資料格式則最佳。

  3. 分析資料：採用任意方法對資料進行分析。

  4. 訓練演算法：目的是為了找到最佳的分類迴歸係數。

  5. 測試演算法

  6. 使用演算法

二、Sigmoid函式

　　若要處理的是二分類問題，我們期望的函式輸出會是0或1，類似於單位階躍函式，可是該函式是不連續的，不連續不可微。

　　因此我們換用Sigmoid函式，當x=0時，y為0；隨著x的增大，y值趨近於1，隨著x的減小，y趨近於0，當橫座標足夠大時，Sigmoid函式就會看起來像一個階躍函式。

　　python繪製Sigmoid影象：

import matplotlib.pyplot as plt
import numpy as np
 
def sigmoid(z):
    return 1.0/(1.0+np.exp(-z))
 
z=np.arange(-6,6,0.05)
plt.plot(z,sigmoid(z))
plt.axvline(0.0,color='k')
plt.axhline(y=0.0,ls='dotted',color='k')
plt.axhline(y=1.0,ls='dotted',color='k')
plt.axhline(y=0.5,ls='dotted'
 
,color='k')
plt.yticks([0.0,0.5,1.0])
plt.ylim(-0.1,1.1)
plt.xlabel('z')
plt.ylabel('$\phi (z)$')
plt.show()

 三、梯度上升演算法

　　訓練演算法時採用了梯度上升演算法，其思想是：要找到某函式的最大值，最好的方法就是沿著該函式的梯度方向探尋。如果梯度記為∇，則函式f(x,y)的梯度由 下式表示：

　　如圖，梯度上升演算法到達每個點後都會重新估計移動的方向。從P0開始，計算完該點的梯度，函式就根據梯度移動到下一點P1。在P1點，梯度再次被重新計算，並沿新的梯度方向移動到P2。如此迴圈迭代，直到滿足停止條件。迭代的過程中，梯度運算元總是保證我們能選取到最佳的移動方向。

python繪製結果：

#邏輯斯蒂分類演算法  
import numpy as np  

#畫圖的函式  
from matplotlib.colors import ListedColormap  
import matplotlib.pyplot as plt  
def plot_decision_regions(X, y, classifier, test_idx=None, resolution=0.02):  
  
    # setup marker generator and color map  
    markers = ('s', 'x', 'o', '^', 'v')  
    colors = ('red', 'blue', 'lightgreen', 'gray', 'cyan')  
    cmap = ListedColormap(colors[:len(np.unique(y))])#通過ListedColormap來定義一些顏色和標記號，並通過顏色列表生成了顏色示例圖  
    # plot the decision surface  
    #對兩個特徵的最大值最小值做了限定（使用兩個特徵來訓練感知器）  
    x1_min, x1_max = X[:, 0].min() - 1, X[:, 0].max() + 1  
    x2_min, x2_max = X[:, 1].min() - 1, X[:, 1].max() + 1  
    #利用meshgrid函式，將最大值、最小值向量生成二維陣列xx1和xx2  
    xx1, xx2 = np.meshgrid(np.arange(x1_min, x1_max, resolution), np.arange(x2_min, x2_max, resolution))  
    #建立一個與資料訓練集中列數相同的矩陣，以預測多維陣列中所有對應點的類標z  
    Z = classifier.predict(np.array([xx1.ravel(), xx2.ravel()]).T)  
    Z = Z.reshape(xx1.shape)#將z變換為與xx1和xx2相同維度  
    #使用contourf函式，對於網格陣列中每個預測的類以不同的顏色繪製出預測得到的決策區域  
    plt.contourf(xx1, xx2, Z, alpha=0.3, cmap=cmap)  
    plt.xlim(xx1.min(), xx1.max())  
    plt.ylim(xx2.min(), xx2.max())  
  
    for idx, cl in enumerate(np.unique(y)):  
        plt.scatter(x=X[y == cl, 0],   
                    y=X[y == cl, 1],  
                    alpha=0.8,   
                    c=colors[idx],  
                    marker=markers[idx],   
                    label=cl,   
                    edgecolor='black')  
  
    # highlight test samples  
    if test_idx:  
        # plot all samples  
        X_test, y_test = X[test_idx, :], y[test_idx]  
  
        plt.scatter(X_test[:, 0],  
                    X_test[:, 1],  
                    c='',  
                    edgecolor='black',  
                    alpha=1.0,  
                    linewidth=1,  
                    marker='o',  
                    s=100,   
                    label='test set')  

#訓練集與測試集的獲取,採用鳶尾花資料集  
from sklearn import datasets  
iris=datasets.load_iris()  
x=iris.data[:,[2,3]]  
y=iris.target  
  
#對資料集進行劃分  
from sklearn.model_selection import train_test_split
#採用scikit-learn中的cross_validation模組中的train_test_split（）函式，隨機將iris資料特徵矩陣x與類標向量y按照3:7劃分為測試資料集和訓練資料集  
x_train,x_test, y_train, y_test =train_test_split(x,y,test_size=0.3, random_state=0)  
  
#為了優化效能，對特徵進行標準化處理  
from sklearn.preprocessing import StandardScaler  
sc=StandardScaler()  
sc.fit(x_train)#通過fit方法，可以計算訓練資料中每個特徵的樣本均值和方差  
x_train_std=sc.transform(x_train)#通過呼叫transform方法，可以使用前面獲得的樣本均值和方差來對資料做標準化處理  
x_test_std=sc.transform(x_test)  

from sklearn.linear_model import LogisticRegression
lr=LogisticRegression(C=1000.0,random_state=0)
lr.fit(x_train_std,y_train)
print("Training Score:%f"%lr.score(x_train_std,y_train))#返回在(X_train,y_train)上的準確率  
print("Testing Score:%f"%lr.score(x_test_std,y_test))#返回在(X_test,y_test)上的準確率  
  
x_combined_std = np.vstack((x_train_std, x_test_std))#將陣列垂直排列成多個子陣列的列表。  
y_combined = np.hstack((y_train, y_test))# 按水平順序（列）順序堆疊陣列。  
plot_decision_regions(X=x_combined_std, y=y_combined, classifier=lr, test_idx=range(105, 150))  
plt.xlabel('petal length [standardized]')  
plt.ylabel('petal width [standardized]')  
plt.legend(loc='upper left')  
plt.show()  

注意：

原始碼交叉驗證的時候使用：

這樣做會報錯。原因是sklearn中已經廢棄cross_validation，將其中的內容整合到了model_selection。需要按如下修改：

•  梯度下降演算法

　　梯度下降演算法，它與這裡的梯度上升演算法是一樣的，只是公式中的加法需要變成減法。因此，對應的公式可以寫成：

 　　梯度上升演算法用來求函式的最大值，而梯度下降演算法用來求函式的最小值。

 四、Logistic迴歸的優缺點

　　優點：計算代價不高，易於理解和實現

　　缺點：容易欠擬合，分類精度可能不高

　　適用資料型別：數值型和標稱型資料

五、python實現分類器的應用

　　下圖所示為資料集test.txt：

　　資料集有100個樣本點，每個點包含兩個數值型特徵：X1和X2。

　　在此資料集上，我們將通過使用梯度上升法找到最佳迴歸係數，也就是擬合出Logistic迴歸模型的最佳引數。

　　梯度上升法的虛擬碼如下：

每個迴歸係數初始化為1
重複R次：
    計算整個資料集的梯度
    使用alpha*gradient更新迴歸係數的向量
返回迴歸係數

　　①求出迴歸係數[w0,w1,w2]

# 載入資料集
def loadDataSet():
    dataMat = []                                                      # 資料列表
    labelMat = []                                                     # 標籤列表
    fr = open('test.txt')                                          # 開啟檔案
    for line in fr.readlines():                                       # 遍歷每行，讀取資料
        lineArr = line.strip().split()                                # 去除回車符
        dataMat.append([1.0, float(lineArr[0]), float(lineArr[1])])   # 把此行的資料新增到資料列表裡面
        labelMat.append(int(lineArr[2]))                              # 新增標籤進列表
    return dataMat, labelMat
# sigmoid函式
def sigmoid(inX):
    return 1.0 / (1 + exp(-inX))
# 梯度上升演算法
def gradAscent(dataMatIn, classLabels):    
    dataMatrix = mat(dataMatIn)           
    labelMat = mat(classLabels).transpose()            
    m, n = shape(dataMatrix)   # 獲取資料集矩陣的大小，m為行數，n為列數
    alpha = 0.001      # 步長
    maxCycles = 500      # 迭代次數
    weights = ones((n, 1))     # 權重初始化為1
    for k in range(maxCycles):     # 重複矩陣運算
        h = sigmoid(dataMatrix * weights)   
        error = (labelMat - h)      # 計算誤差
        weights = weights + alpha * dataMatrix.transpose() * error 
    return weights
# 獲得迴歸係數
dataMat, labelMat = loadDataSet()
weigths = gradAscent(dataMat, labelMat)
print("w0: %f, w1: %f, W2: %f" % (weigths[0], weigths[1], weigths[2]))

　　②分析資料，畫出決策邊界

# 繪製資料集和Logistic迴歸最佳擬合直線
def plotBestFit(weights):
    dataMat, labelMat = loadDataSet()      
    dataArr = array(dataMat)      
    n = shape(dataArr)[0]     # 獲取資料總數
    xcord1 = []; ycord1 = []       # 存放正樣本
    xcord2 = []; ycord2 = []       # 存放負樣本
    for i in range(n):                                     # 依據資料集的標籤來對資料進行分類
        if int(labelMat[i]) == 1:                          # 資料的標籤為1，表示為正樣本
            xcord1.append(dataArr[i, 1]); ycord1.append(dataArr[i, 2])
        else:      
            xcord2.append(dataArr[i, 1]); ycord2.append(dataArr[i, 2])
    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.scatter(xcord1, ycord1, s=30, c='blueviolet', marker='s')  # 繪製正樣本
    ax.scatter(xcord2, ycord2, s=30, c='pink')      # 繪製負樣本
    x = arange(-3.0, 3.0, 0.1)      
    y = (-weights[0] - weights[1] * x) / weights[2]    
    ax.plot(x, y)
    plt.title('cxr_test')          
    plt.xlabel('X1'); plt.ylabel('X2')

    plt.show()
    dataMat, labelMat = loadDataSet()
    weigths = gradAscent(dataMat, labelMat)
    plotBestFit(weigths.getA())
    print("w0: %f, w1: %f, W2: %f" % (weigths[0], weigths[1], weigths[2]))

　　這個分類器分類效果不錯，從圖上看只分錯了幾個點而已，但是回顧對演算法的原理的解析，我們卻發現這個這個方法需要相當大量的計算，因此下一步我們將對其進行改進，從而使它可以運用在我們現實生活中，處理例子複雜數量大的資料集。

　　③訓練演算法：隨機梯度上升

　　一種改進的方法是一次僅用一個樣本點來更新迴歸係數，即隨機梯度上升法。由於可以在新樣本到來時對分類器進行增量式更新，因此隨機梯度上升法是一個線上學習演算法。

　　隨機梯度上升法可以寫成如下虛擬碼：

所有迴歸係數初始化為1
對資料集中每個樣本
        計算該樣本的梯度
　　    使用alpha*gradient更新回顧係數值
 返回迴歸係數值

• 隨機梯度上升和梯度上升法在程式碼上的區別：

　　（1）後者的變數h和誤差error都是向量。
　　（2）前者沒有矩陣的轉換過程，所有變數的資料型別都是Numpy陣列。

def stocGradAscent0(dataMatrix, classLabels):   
    m, n = shape(dataMatrix)  
    alpha = 0.01     #步長
    weights = ones(n)     #初始化為1
    for i in range(m):    # 重複矩陣運算
        h = sigmoid(sum(dataMatrix[i] * weights))   
        error = classLabels[i] - h     # 計算誤差
        weights = weights + alpha * error * dataMatrix[i]     #更新權重
    return weights
    
    dataMat, labelMat = loadDataSet()
    weigths = stocGradAscent0(array(dataMat), labelMat)
    plotBestFit(weigths)

 　　擬合曲線錯了三分之一的樣本，要進行改進。

　　④改進隨機梯度上升

　　判斷優化演算法的優劣的可靠方法是看它是否收斂，也就是看引數是否穩定。

　　由上圖可知，x2只經過50次迭代就達到了穩定值，但其他兩個需要更多次迭代。

•  有三點改進：

　　（1）alpha每次迭代的時候都會調整。不會減小到0，滿足某個條件時也不是嚴格下降的。
　　（2）隨機選取樣本來更新迴歸係數，這種做法會減少週期性波動。
　　（3）增加了迭代引數的設定。預設為150次，如果給定，將會按照新的引數值進行迭代。

　　對此進行改進的隨機梯度上升演算法的程式碼如下：

def stocGradAscent1(dataMatrix, classLabels, numIter=150):  
    m, n = shape(dataMatrix)           
    weights = ones(n)              
    for j in range(numIter):
        dataIndex = list(range(m))                 
        for i in range(m):
            alpha = 4 / (1.0 + j + i) + 0.0001           
            randIndex = int(random.uniform(0, len(dataIndex)))  
            h = sigmoid(sum(dataMatrix[randIndex] * weights))  
            error = classLabels[randIndex] - h          
            weights = weights + alpha * error * dataMatrix[randIndex]  
            del (dataIndex[randIndex])          
    return weights

    dataMat, labelMat = loadDataSet()
    weigths = stocGradAscent1(array(dataMat), labelMat)
    plotBestFit(weigths)

　　效果很好，而且計算量少。

　　改進後的迴歸係數和迭代次數的關係圖：

 　　可以發現：

　　（1）此關係圖中係數沒有之前的迭代次數圖那樣出現週期性的波動。
　　（2）水平軸比之前的圖短了很多。

 六、 python例項：從疝氣病症狀預測病馬的死亡率

　　①準備資料，處理資料中的缺失值

　　流程：
　　（1）收集資料：給定資料檔案。
　　（2）準備資料：用python解析文字檔案，並填充缺失值。
　　（3）分析資料：視覺化並觀察資料。
　　（4）訓練演算法：使用優化演算法，找到最佳係數。
　　（5）測試演算法：為了量化迴歸的效果，需要觀察錯誤率，根據錯誤率決定是否退到訓練階段，通過改變迭代的次數和步長等引數來得到更好的迴歸係數。

　　資料包含了368個樣本和28個特徵。

　　疝氣病是描述馬胃腸痛的術語。這種病不一定源自馬的腸胃問題，其他問題也可能引發馬疝病。

　　該資料集中包含了醫院檢測馬疝病的一些指標，有的指標比較主觀，有的指標難以測量，例如馬的疼痛級別。

　　預處理資料階段需要做兩件事：
　　（1）所有的缺失值必須用一個實數值來替換，因為使用的NumPy資料型別不允許包含缺失值。這裡選擇實數0來替換所有缺失值，恰好能適用於Logistic迴歸。

　　迴歸係數的更新公式：

　　weights = weights + alpha * error * dataMatrix[randIndex]

　　如果說dataMatrix的某特徵對應值是0，那麼該特徵的係數將不做更新，即：

　　weights = weights

　　由於sigmoid(0)=0.5，即它對結果的預測不懼任何傾向性，因此不會對誤差造成任何影響。

　　（2）如果在測試資料集中發現了一條資料的類別標記已經缺失，那麼簡單的做法就是將該條資料丟棄。

　　②訓練演算法，使用優化演算法，找到最佳係數

from numpy import *

def sigmoid(inX):
    return 1.0/(1+exp(-inX))

def stoGradAscent1(dataMatrix, classLabels, numIter = 150):
    m,n = shape(dataMatrix)
    weights = ones(n)
    for j in range (numIter):
        dataIndex = range(m)
        for i in range(m):
            alpha = 4 / (1.0 + j + i) + 0.01
            randIndex = int(random.uniform(0, len(dataIndex)))
            h = sigmoid(sum(dataMatrix[randIndex] * weights))
            error = classLabels[randIndex] - h
            weights = weights + alpha * error * dataMatrix[randIndex]
            del (list(dataIndex)[randIndex])
    return weights

　　③測試演算法，用Logistic迴歸進行分類

def classifyVector(inX, weights):
    prob = sigmoid(sum(inX * weights))
    #大於0.5 返回 1；否則返回0
    if prob > 0.5:
        return 1.0
    else:
        return 0.0

def colicTest():
    frTrain = open('HorseColicTraining.txt')
    frTest = open('HorseColicTest.txt')
    trainingSet = []
    trainingLabels = []
    
    for line in frTrain.readlines():
        currLine = line.strip().split('\t') #分割
        lineArr = []
        for i in range(21):
            lineArr.append(float(currLine[i]))
        #存入訓練樣本特徵
        trainingSet.append(lineArr)
        #存入訓練樣本標籤
        trainingLabels.append(float(currLine[21]))
    #使用改進後的隨機梯度下降演算法得到迴歸係數
    trainingWeights = stoGradAscent1(array(trainingSet), trainingLabels, 500)

    #統計測試集預測錯誤樣本數量和樣本總數
    errorCount = 0
    numTestVec = 0.0
    for line in frTest.readlines():
        #迴圈一次，樣本數加1
        numTestVec += 1.0
        currLine = line.strip().split('\t') #分割
        lineArr = []
        for i in range(21):
            lineArr.append(float(currLine[i]))
        
        # 利用分類預測函式對該樣本進行預測，並與樣本標籤進行比較
        if int(classifyVector(array(lineArr), trainingWeights)) != int(currLine[21]):
            errorCount += 1
    #計算錯誤率
    errorRate = (float(errorCount) / numTestVec)
    print('the error rate of this test is : %f' % errorRate)
    return errorRate

#呼叫colicTest10次，求平均值
def multiTest():
    numTests = 10
    errorSum = 0.0
    for k in range(numTests):
        errorSum += colicTest()
    print('after %d iterations the average error rete is : %f ' % (numTests,errorSum / float(numTests)))

　　測試結果：

 　　平均值為34%。因為資料集本身有30%的資料缺失。後續可以使用梯度上升法進一步改進。

總結：

• 梯度上升演算法：每次更新迴歸係數所有樣本都參與。

　　優點：分類準確，獲取全域性最優解
　　缺點：當樣本比較多時，訓練速度特別慢
　　適用場合：樣本較少的資料集
　　

• 隨機梯度下降法：每次更新迴歸係數只有一個樣本參與。

　　優點：訓練速度很快
　　缺點：準確率會降低，並不是朝著整體最優方向進行，容易獲取到區域性最優解
　　適用場合：樣本非常多的資料集