邏輯回歸 - 簡書

邏輯回歸
邏輯回歸是一種經(jīng)典的二分類算法，它的決策邊界可以是非線性的（包含高階項(xiàng)）
Sigmoid函數(shù)
公式： $g(z)=\frac{1}{1+{{e}^{-z}}}$ 柿扣，其中自變量取值為任意實(shí)數(shù)
解釋：將任意的輸入映射到了[0,1]區(qū)間
比如在線性回歸中可以得到一個(gè)預(yù)測值（ $\hat{y}$ ）肖方，再將該值映射到Sigmoid函數(shù)，就完成了由值到概率的轉(zhuǎn)換未状，也就是分類任務(wù)
假設(shè)預(yù)測函數(shù) ${{h}_{\mathbf{\theta }}}(x)=g({{\mathbf{\theta }}^{\mathbf{T}}}\mathbf{X})=\frac{1}{1+{{e}^{-{{\mathbf{\theta }}^{\mathbf{T}}}\mathbf{X}}}}$
分類任務(wù)： $\begin{align} & p(y=1|\mathbf{X};\mathbf{\theta })={{h}_{\mathbf{\theta }}}(x) \\ & p(y=0|\mathbf{X};\mathbf{\theta })=1-{{h}_{\mathbf{\theta }}}(x) \\ \end{align}$
相當(dāng)于0-1分布俯画，整合上述式子得 $p(y|\mathbf{X};\mathbf{\theta })={{({{h}_{\mathbf{\theta }}}(x))}^{y}}{{(1-{{h}_{\mathbf{\theta }}}(x))}^{1-y}}$
似然函數(shù)：
$L(\mathbf{\theta })=\prod\limits_{i=1}^{m}{p(y|\mathbf{X};\mathbf{\theta })}=\prod\limits_{i=1}^{m}{{{({{h}_{\mathbf{\theta }}}({{x}_{i}}))}^{{{y}_{i}}}}{{(1-{{h}_{\mathbf{\theta }}}({{x}_{i}}))}^{1-{{y}_{i}}}}}$
取對數(shù)似然函數(shù)：
$\ln L(\mathbf{\theta })=\sum\limits_{i=1}^{m}{{{y}_{i}}\ln ({{h}_{\theta }}({{x}_{i}}))\text{+(1-}}{{y}_{i}})\ln (1-{{h}_{\theta }}({{x}_{i}}))$
現(xiàn)在要求上述對數(shù)似然函數(shù)的最大值，如果要用梯度下降法司草，需要將目標(biāo)函數(shù)轉(zhuǎn)化為求最小值【梯度的方向是值越來越大的方向艰垂，求最小值，應(yīng)該朝著梯度反方向走】埋虹，因此引入
$J(\mathbf{\theta })=-\frac{1}{m}\ln L(\mathbf{\theta })=-\frac{1}{m}\sum\limits_{i=1}^{m}{{{y}_{i}}\ln ({{h}_{\theta }}({{x}_{i}}))\text{+(1-}}{{y}_{i}})\ln (1-{{h}_{\theta }}({{x}_{i}}))$
對上式求偏導(dǎo)
$\begin{align} & \frac{\partial J(\mathbf{\theta })}{\partial {{\theta }_{j}}}=-\frac{1}{m}\sum\limits_{i=1}^{m}{\left[ {{y}_{i}}\frac{1}{{{h}_{\mathbf{\theta }}}({{x}_{i}})}\frac{\partial {{h}_{\mathbf{\theta }}}({{x}_{i}})}{\partial {{\theta }_{j}}}-(1-{{y}_{i}})\frac{1}{1-{{h}_{\mathbf{\theta }}}({{x}_{i}})}\frac{\partial {{h}_{\mathbf{\theta }}}({{x}_{i}})}{\partial {{\theta }_{j}}} \right]} \\ & \text{ =}-\frac{1}{m}\sum\limits_{i=1}^{m}{\left[ {{y}_{i}}\frac{1}{{{h}_{\mathbf{\theta }}}({{x}_{i}})}-(1-{{y}_{i}})\frac{1}{1-{{h}_{\mathbf{\theta }}}({{x}_{i}})} \right]}\frac{\partial {{h}_{\mathbf{\theta }}}({{x}_{i}})}{\partial {{\theta }_{j}}} \\ & \text{ =}-\frac{1}{m}\sum\limits_{i=1}^{m}{\left[ {{y}_{i}}\frac{1}{g({{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{i}})}-(1-{{y}_{i}})\frac{1}{1-g({{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{i}})} \right]}\frac{\partial g({{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{i}})}{\partial {{\theta }_{j}}} \\ & \text{ =}-\frac{1}{m}\sum\limits_{i=1}^{m}{\left[ {{y}_{i}}\frac{1}{g({{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{i}})}-(1-{{y}_{i}})\frac{1}{1-g({{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{i}})} \right]}\frac{{{e}^{-{{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{\mathbf{i}}}}}}{1+{{e}^{-{{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{\mathbf{i}}}}}}x_{i}^{j} \\ & \text{ =}-\frac{1}{m}\sum\limits_{i=1}^{m}{\left[ {{y}_{i}}\frac{1}{g({{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{i}})}-(1-{{y}_{i}})\frac{1}{1-g({{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{i}})} \right]}g({{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{i}})(1-g({{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{i}}))x_{i}^{j} \\ & \text{ =}-\frac{1}{m}\sum\limits_{i=1}^{m}{\left[ {{y}_{i}}(1-g({{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{i}}))-(1-{{y}_{i}})g({{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{i}}) \right]}x_{i}^{j} \\ & \text{ =}-\frac{1}{m}\sum\limits_{i=1}^{m}{\left[ {{y}_{i}}-g({{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{i}}) \right]}x_{i}^{j} \\ & \text{ =}\frac{1}{m}\sum\limits_{i=1}^{m}{\left[ {{h}_{\mathbf{\theta }}}({{x}_{i}})-{{y}_{i}} \right]}x_{i}^{j} \\ \end{align}$
確定了梯度方向之后驹马，利用上式更新參數(shù)秩冈，得到
$\theta _{j}^{'}={{\theta }_{j}}-\alpha \frac{1}{m}\sum\limits_{i=1}^{m}{\left[ {{h}_{\mathbf{\theta }}}({{x}_{i}})-{{y}_{i}} \right]}x_{i}^{j}$
按照上式繼續(xù)迭代
下面是一個(gè)關(guān)于實(shí)現(xiàn)邏輯回歸的小例子
目標(biāo)：建立一個(gè)邏輯回歸模型來預(yù)測一個(gè)大學(xué)生是否被大學(xué)錄取潜腻。
特征：兩次考試成績
導(dǎo)入數(shù)據(jù)

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
pdData=pd.read_csv('C:/Users/lenovo/Desktop/LogiReg_data.txt',header=None,names=['Exam 1','Exam 2','Admitted'])#header=None不讓數(shù)據(jù)自己指定列名
pdData.insert(0,'ones',1) #為了考慮常數(shù)項(xiàng)
pdData.head()

數(shù)據(jù)如以下所示

image.png

通常在做數(shù)據(jù)處理之前遂铡，需要做一些數(shù)據(jù)預(yù)處理，比如數(shù)據(jù)標(biāo)準(zhǔn)化等辣辫，這邊先不做旦事，看看不做標(biāo)準(zhǔn)化會對結(jié)果什么影響。
邏輯回歸的關(guān)鍵是計(jì)算損失值和梯度方向
為了計(jì)算這兩個(gè)函數(shù)急灭，需要引入sigmoid函數(shù)【將值映射到概率】

#定義sigmoid函數(shù)
def sigmoid_2(z):
    return 1/(1+np.exp(-z))
#預(yù)測函數(shù)模型
def model(X,theta):
    return sigmoid_2(np.dot(X,theta.T))

切片姐浮，分出X和y

orig_data=pdData.as_matrix()# 轉(zhuǎn)換成array類型
cols=orig_data.shape[1]
X=orig_data[:,:cols-1]#前三列,切片012
y=orig_data[:,cols-1:]#和y=orig_data[:,cols-1]結(jié)果額維度不一樣
theta=np.zeros([1,3]) #theta初始化是0

注：這里要注意一個(gè)關(guān)于數(shù)組切片的問題

image.png

上面畫圈的兩個(gè)語句一個(gè)返回一維數(shù)組，一個(gè)返回三維數(shù)組

損失函數(shù)

將對數(shù)似然函數(shù)取相反數(shù)
$-\ln L(\mathbf{\theta })=-\sum\limits_{i=1}^{n}{{{y}_{i}}\ln ({{h}_{\theta }}({{x}_{i}}))-\text{(1-}}{{y}_{i}})\ln (1-{{h}_{\theta }}({{x}_{i}}))$
求平均損失
$J(\mathbf{\theta })=\frac{1}{n}\ln L(\mathbf{\theta })=\frac{1}{n}\sum\limits_{i=1}^{n}{\text{-}{{y}_{i}}\ln ({{h}_{\theta }}({{x}_{i}}))\text{-(1-}}{{y}_{i}})\ln (1-{{h}_{\theta }}({{x}_{i}}))$
目的是求平均損失的最小值

def cost(X,y,theta):#對應(yīng)于J(theta)
    left=np.multiply(-y,np.log(model(X,theta)))#np.multiply對應(yīng)位置相乘
    right=np.multiply(1-y,np.log(1-model(X,theta)))
    return np.sum(left-right)/len(X)

計(jì)算梯度

$\frac{\partial J(\mathbf{\theta })}{\partial {{\theta }_{j}}}\text{=}-\frac{1}{n}\sum\limits_{i=1}^{n}{\left[ {{y}_{i}}-g({{\mathbf{\theta }}^{\mathbf{T}}}{{\mathbf{x}}_{i}}) \right]}x_{i}^{j}\text{=}\frac{1}{n}\sum\limits_{i=1}^{n}{\left[ {{h}_{\mathbf{\theta }}}({{x}_{i}})-{{y}_{i}} \right]}x_{i}^{j}$

def gradient(X,y,theta):
    grad=np.zeros(theta.shape)
    error=(model(X,theta)-y).ravel()
   # print(error.shape)
    for j in range(len(theta.ravel())):#求每個(gè)參數(shù)的偏導(dǎo)
        term=np.multiply(error,X[:,j])
        grad[0,j]=np.sum(term)/len(X)
    return grad

參數(shù)更新

首先設(shè)定三種停止策略【迭代次數(shù)葬馋、損失值差距卖鲤、梯度】

def stopCriterion(type,value,threshold):
    if type==STOP_ITER:
        return value>threshold #返回邏輯值#迭代次數(shù)
    elif type==STOP_COST:
        return abs(value[-1]-value[-2])<threshold#兩次迭代結(jié)果的損失值差距較小
    elif type==STOP_GRAD:
        return np.linalg.norm(value)<threshold #梯度較小

為了避免原數(shù)據(jù)存在某種規(guī)律【比如先收集女生成績再收集男生成績】，將原數(shù)據(jù)順序打亂

import numpy.random
#洗牌畴嘶，將數(shù)據(jù)順序打亂
def shuffleData(data):
    np.random.shuffle(data)  #shuffle() 方法將序列的所有元素隨機(jī)排序
    cols=data.shape[1]
    X=data[:,:cols-1]
    y=data[:,cols-1:]
    return X,y

然后進(jìn)行參數(shù)更新【引入time庫是為了比較不同策略的運(yùn)行時(shí)間】

import time
def descent(data,theta,batchSize,stopType,thresh,alpha):
    #batchSize等于總樣本數(shù)蛋逾，即批量梯度下降；batchSize等于1窗悯，即隨機(jī)梯度下降区匣；batchSize取1~n之間的數(shù)據(jù)，即小批量梯度下降
    
    
    #梯度下降求解
    #初始化
    init_time=time.time()  #比較不同梯度下降方法的運(yùn)行速度
    i=0 #迭代次數(shù)蒋院，第0次開始
    k=0 #batch
    X,y=shuffleData(data)
    grad=np.zeros(theta.shape)  #計(jì)算的梯度
    costs=[cost(X,y,theta)]  #計(jì)算損失值
    
    
    
    while True:
        grad=gradient(X[k:k+batchSize],y[k:k+batchSize],theta)
        k += batchSize #取樣本數(shù)據(jù)
        if k >= n:
            k = 0
            X,y=shuffleData(data) #重新洗牌
        theta = theta - alpha*grad  # 參數(shù)更新
        costs.append(cost(X,y,theta))
        i += 1
        
        #停止策略
        if stopType==STOP_ITER:
            value = i
        elif stopType==STOP_COST:
            value=costs
        elif stopType==STOP_GRAD:
            value=grad
        if stopCriterion(stopType,value,thresh):#如果if語句為真亏钩，就跳出整個(gè)循環(huán)
            break
        #print(grad,np.linalg.norm(grad))
    return theta,i-1,costs,grad,time.time()-init_time

注：np.linalg.norm([4,3])表示 $\sqrt{{{3}^{2}}+{{4}^{2}}}$

def runExpe(data,theta,batchSize,stopType,thresh,alpha):
    theta,iter,costs,grad,dur=descent(data,theta,batchSize,stopType,thresh,alpha)
#     name='Original' if (data[:,2]>2).sum()>1 else 'Scaled'
#     name += 'data-learning rate: {} -'.format(alpha)
#     if batchSize==n: strDescType='Gradient'
#     elif batchSize==1: strDescType='Stochastic'
#     else:strDescType='mini-batch {} '.format(batchSize)
#     name += strDescType + 'descent-stop: '
#     if stopType==STOP_ITER: strStop='{} iterations'.format(thresh)
#     elif stopType==STOP_COST: strStop='cost change < {}'.format(thresh)
#     else: strStop ='gradient norm < {}'.format(thresh)
#     name += strStop
    fig,ax=plt.subplots(figsize=(12,4))
    ax.plot(np.arange(len(costs)),costs,c='r')
    ax.set_xlabel('Iteration')
    ax.set_ylabel('Cost')
    #ax.set_title(name.upper())
    print('iter: {}, last cost: {:03.2f}, duration: {:03.2f}s'.format(iter,costs[-1],dur))
    return theta

不同的停止策略

以批量梯度下降為例

停止條件為迭代次數(shù)

n=100
runExpe(orig_data,theta,n,STOP_ITER,thresh=5000,alpha=0.000001)

image.png

看似損失值已經(jīng)穩(wěn)定在最低點(diǎn)0.63

停止條件為損失值
設(shè)定閾值為0.000001，需要迭代110000次左右

runExpe(orig_data,theta,n,STOP_COST,thresh=0.000001,alpha=0.001)

image.png

損失值最低為0.38欺旧，似乎還可以進(jìn)一步收斂

停止條件為梯度大小
設(shè)定閾值0.05姑丑，需要迭代40000次左右

runExpe(orig_data,theta,n,STOP_GRAD,thresh=0.05,alpha=0.001)

image.png

損失值最小為0.49，似乎還可以進(jìn)一步收斂
綜上辞友，基于批量梯度下降方法栅哀，上述三種停止條件得到的損失函數(shù)值為0.63震肮、0.38和0.49，迭代次數(shù)分別為5000次留拾、110000次和40000次戳晌，迭代次數(shù)越多，損失值越小

對比不同的梯度下降方法

停止策略為迭代次數(shù)

隨機(jī)梯度下降

runExpe(orig_data,theta,1,STOP_ITER,thresh=5000,alpha=0.001)

image.png

波動非常大间驮，迭代過程不穩(wěn)定躬厌，這也是隨機(jī)梯度下降的主要缺點(diǎn)
嘗試降低學(xué)習(xí)率為0.000001马昨，增加迭代次數(shù)為15000

runExpe(orig_data,theta,1,STOP_ITER,thresh=15000,alpha=0.000001)

image.png

效果要好一些竞帽，損失值似乎穩(wěn)定在0.63,根據(jù)上面的結(jié)果可知，0.63不算是一個(gè)特別合理的值

小批量梯度下降

#取樣本為16
runExpe(orig_data,theta,16,STOP_ITER,thresh=15000,alpha=0.001)

image.png

上下波動鸿捧，迭代過程不穩(wěn)定
嘗試調(diào)低學(xué)習(xí)率為0.000001

runExpe(orig_data,theta,16,STOP_ITER,thresh=15000,alpha=0.001)

image.png

降低學(xué)習(xí)率之后沒有效果屹篓，迭代過程依舊不穩(wěn)定
因此，可能不是模型本身的問題匙奴，而是數(shù)據(jù)本身的問題堆巧，嘗試著對數(shù)據(jù)做一些變換，此處對數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化泼菌，用標(biāo)準(zhǔn)化后的數(shù)據(jù)求解

#標(biāo)準(zhǔn)化
from sklearn import preprocessing as pp
scaled_data=orig_data.copy()
scaled_data[:,1:3]=pp.scale(scaled_data[:,1:3])
#用標(biāo)準(zhǔn)化后的數(shù)據(jù)求解
runExpe(scaled_data,theta,16,STOP_ITER,thresh=15000,alpha=0.001)

image.png

損失值收斂到0.28谍肤，比0.63好很多
再嘗試一下梯度為停止條件的情況

runExpe(scaled_data,theta,16,STOP_GRAD,thresh=0.004,alpha=0.001)

image.png

迭代次數(shù)由15000增加到60000多，損失值由0.28降低到0.22哗伯，又改善了一步

計(jì)算模型分類結(jié)果的精確率

#設(shè)定閾值為0.5荒揣，大于0.5就可以入學(xué)
def predict(X,theta):
    return [1 if x >= 0.5 else 0 for x in model(X,theta)]

scaled_X=scaled_data[:,:3]
y=scaled_data[:,3]
theta = runExpe(scaled_data,theta,16,STOP_GRAD,thresh=0.004,alpha=0.001)
predictions = predict(scaled_X,theta)
correct = [1 if a==b else 0 for (a,b) in zip(predictions,y)]   #真實(shí)值與預(yù)測值相等，同為1或者同為0
accuracy = sum(correct)/len(correct)

返回accuracy等于0.89
通過這個(gè)例子焊刹，算是對邏輯回歸的基本原理有一個(gè)比較清晰的認(rèn)識了系任。

sigmoid函數(shù)：將值映射到概率的函數(shù)
model：返回預(yù)測結(jié)果值
cost：根據(jù)參數(shù)計(jì)算損失
gradient：計(jì)算每個(gè)參數(shù)的梯度方向
descent：進(jìn)行每個(gè)參數(shù)的更新
accuracy：計(jì)算精度