机器学习（三）：逻辑回归之从理论到实践

刺骨的言语ヽ痛彻心扉 2022-05-23 13:17 248阅读 0赞

文中的代码和数据集下载地址：  
[https://github.com/TimePickerWang/MachineLearningInAction][https_github.com_TimePickerWang_MachineLearningInAction]

## 一、基本理论 ##

逻辑回归是一个二值型的分类器，它是利用Sigmoid函数来进行分类的，Sigmoid函数的表达式如下：

####  g(z)=11\+e−z  g ( z ) = 1 1 + e − z  ####

其函数图像如下图：  
![这里写图片描述][70]  
结合表达式可以知道：  
当z=0时，g(z)=0.5  当 z = 0 时 ， g ( z ) = 0.5   
当z>0时，g(z)>0.5  当 z > 0 时 ， g ( z ) > 0.5 ，且随着z的增大，g(z)会越来越接近1  
当z<0时，g(z)<0.5  当 z < 0 时 ， g ( z ) < 0.5 ，且随着z的减小，g(z)会越来越接近0

利用这一特性，我们可以进行这样一种分类方式：  
当输入z>0.5  z > 0.5 时，把样本分为1类  
当输入z<0.5  z < 0.5 时，把样本分为0类

那么z  z 具体是什么呢？假设每一个样本有n个特征：\{x1,x2,x3......xn x 1,x 2,x 3. . . . . .x n\}，则：

z=w0x0\+w1x1\+w2x2\+...\+wnxn  z = w 0 x 0 + w 1 x 1 + w 2 x 2 + . . . + w n x n   
　　  
　　其中，  x0=1  x 0 = 1  是一个常数，表示一个偏差量。  w0,w1,w2,...wn  w 0 , w 1 , w 2 , . . . w n  可以理解为每一个特征的权重，我们的目的就是寻求一组  w0,w1,w2,...wn  w 0 , w 1 , w 2 , . . . w n  ，使得分类正确的概率尽可能的大。

那么问题来了，怎么表示分类正确的概率呢？假设有一个样本i  i ，其类标签是y(y=0或1) y ( y = 0或1 )，则可以利用上面提到的“g(z)  g ( z ) ”表示样本i  i 为类1的概率，由于Sigmoid的值域为（０,１），所以可以令”1−g(z) 1 − g ( z )”表示样本为类0的概率。则样本i  i 分类正确的概率可以表示如下：

p=g(z)y(1−g(z))1−y p = g ( z) y( 1 − g ( z ))1 − y

上式中，当y=1  y = 1 时，p=g(z)  p = g ( z ) ；当y=0  y = 0 时，p=1−g(z)  p = 1 − g ( z ) 。所以无论样本i  i 的类别是0还是1，我们总是希望p p的值越大越好。但是，求上式的最大值是很不方便的，我们可以利用对数似然函数来对上式取对数，即：

p=ylog(g(z))\+(1−y)log(1−g(z))  p = y l o g ( g ( z ) ) + ( 1 − y ) l o g ( 1 − g ( z ) )   
所以当有m个样本时，样本分类正确的概率可以这样表示：

L(W)=∑mi=1y(i)log(g(z(i)))\+(1−y(i))log(1−g(z(i)))  L ( W ) = ∑ i = 1 m y ( i ) l o g ( g ( z ( i ) ) ) + ( 1 − y ( i ) ) l o g ( 1 − g ( z ( i ) ) )   
其中  y(i)  y ( i )  表示第  i  i  个样本的类标签，  z(i)=w0x(i)0\+w1x(i)1\+w2x(i)2\+...\+wnx(i)n z( i )=w 0x 0( i )\+w 1x 1( i )\+w 2x 2( i )\+ . . . \+w nx n( i ) ,这里\{  x(i)0,x(i)1,x(i)2,......x(i)n  x 0 ( i ) , x 1 ( i ) , x 2 ( i ) , . . . . . . x n ( i )  \}表示第  i  i  个样本的特征。由于  L L 是关于权重  w0,w1,w2,...wn  w 0 , w 1 , w 2 , . . . w n  的函数，所以这里用  L(W)  L ( W )  表示样本分类正确的概率。

以上内容可以用向量可以表示如下：

W=\[w0,w1,w2,...wn\]  W = \[ w 0 , w 1 , w 2 , . . . w n \]

X=⎡⎣⎢⎢⎢⎢⎢⎢⎢x(1)0x(2)0⋮x(m)0x(1)1x(2)1⋮x(m)1⋯⋯⋯⋯⋯⋯⋱⋯⋯⋯x(1)nx(2)n⋮x(m)n⎤⎦⎥⎥⎥⎥⎥⎥⎥  X = \[ x 0 ( 1 ) x 1 ( 1 ) ⋯ ⋯ ⋯ x n ( 1 ) x 0 ( 2 ) x 1 ( 2 ) ⋯ ⋯ ⋯ x n ( 2 ) ⋮ ⋮ ⋱ ⋮ x 0 ( m ) x 1 ( m ) ⋯ ⋯ ⋯ x n ( m ) \]

Z=XWT  Z = X W T

为了求L(W)  L ( W ) 的最大值，下面介绍一个优化算法：梯度上升算法。  
函数y=−2x2\+5  y = − 2 x 2 + 5 的图像如下入所示：  
![这里写图片描述][70 1]  
我们知道其最大值点为（0， 5）。为了求其最大值，用如下方式更新x的值：

x=x\+α∂y∂x(\*)  (\*) x = x + α ∂ y ∂ x

其中α>0  α > 0 为学习速率，是一个常量。∂y∂x  ∂ y ∂ x 表示y对x的偏导数，这里可以理解为y在x处的斜率。从上图可以看出，在x从-3增长到-1时，每个点所在的切线斜率越来越小，即∂y∂x  ∂ y ∂ x 的值越来越小。结合图像和(∗)  ( ∗ ) 式可知：  
当x<0时，∂y∂x  ∂ y ∂ x >0，且其值随着x的增长越来越小，越来越趋近于0。根据(∗)  ( ∗ ) 式，x增长的越缓慢，越来越趋近于0  
当x>0时，∂y∂x  ∂ y ∂ x <0，且其值随着x的减小越来越大，越来越趋近于0。根据(∗)  ( ∗ ) 式，x减小的越缓慢，越来越趋近于0  
总之最后∂y∂x  ∂ y ∂ x 是越来越趋近于0的，而我们知道，∂y∂x=0  ∂ y ∂ x = 0 的点是函数的极值点，这里就是其最大值点。在这里学习速率α  α 的选取是很重要的，如果太大，则x一次可能增加或减少很大的值，有可能越过极值点；如果太小，则需要更新更多的次数才能到达极致点。

那么之前的L(W)  L ( W ) 怎么运用梯度上升了?由于L(W)  L ( W ) 是关于权重w0,w1,w2,...wn  w 0 , w 1 , w 2 , . . . w n 的函数，所以需要对wi  w i 求导。先从只有一个样本的情况来看，因为有：

#####  p=ylog(g(z))\+(1−y)log(1−g(z))  p = y l o g ( g ( z ) ) + ( 1 − y ) l o g ( 1 − g ( z ) )  #####

#####  g(z)=11\+e−z  g ( z ) = 1 1 + e − z  #####

#####  z=w0x0\+w1x1\+w2x2\+...\+wnxn  z = w 0 x 0 + w 1 x 1 + w 2 x 2 + . . . + w n x n  #####

则有：

###  ∂p∂wi=dpdg(z)dg(z)dz∂z∂wi  ∂ p ∂ w i = d p d g ( z ) d g ( z ) d z ∂ z ∂ w i  ###

而

#####  dpdg(z)=yg(z)−1−y1−g(z)  d p d g ( z ) = y g ( z ) − 1 − y 1 − g ( z )  #####

#####  dg(z)dz=e−z(1\+e−z)2=11\+e−z(1−11\+e−z)=g(z)(1−g(z))  d g ( z ) d z = e − z ( 1 + e − z ) 2 = 1 1 + e − z ( 1 − 1 1 + e − z ) = g ( z ) ( 1 − g ( z ) )  #####

#####  ∂z∂wi=xi  ∂ z ∂ w i = x i  #####

所以：

###  ∂p∂wi=(y−g(z))xi  ∂ p ∂ w i = ( y − g ( z ) ) x i  ###

对每个权重wi  w i ，运用梯度上升为：

###  wi=wi\+α(y−a)xi  w i = w i + α ( y − a ) x i  ###

其中，a=g(z)  a = g ( z ) ，即是分类的预测值。

## 二、向量化 ##

接下来，我们把这一过程向量化。  
之前有

W=\[w0,w1,w2,...wn\]  W = \[ w 0 , w 1 , w 2 , . . . w n \]

Z=XWT  Z = X W T   
m个样本的标签向量：  
 Y=\[y(1),y(2),...,y(m)\]T  Y = \[ y ( 1 ) , y ( 2 ) , . . . , y ( m ) \] T   
m个样本的分类预测：  
 A=g(Z)  A = g ( Z )

以上各个矩阵的维度分别是：

W−(1,n\+1)  W − ( 1 , n + 1 )   
 X−(m,n\+1)  X − ( m , n + 1 )   
 Z,A,Y−(m,1)  Z , A , Y − ( m , 1 )

因此，以上的梯度上升用矩阵表示如下：

W=W\+α(Y−A)T∗X  W = W + α ( Y − A ) T ∗ X

## 三、用python3实现逻辑回归算法 ##

首先获取数据，代码如下：  
主要代码如下：

# sigmoid函数
    def sigmoid(x):
        return 1.0/(1 + np.exp(-x))
    
    
    # 梯度上升算法,data_mat样本集，class_label样本类别标签，iter_num迭代次数，alpha为学习率
    def grad_ascent(data_mat, class_label, iter_num=500, alpha=0.001):
        weights = np.ones((data_mat.shape[1], 1))
        for i in range(iter_num):
            hyp = sigmoid(np.dot(data_mat, weights))
            error = class_label - hyp
            weights += alpha * np.dot(data_mat.transpose(), error)
        return weights
    
    # 分类函数，inx为待分类的样本，weights为权重向量
    def classify(inx, weights):
        prob = sigmoid(np.sum(np.dot(inx, weights)))
        return 1 if prob > 0.5 else 0

以上代码中sigmoid()对应之前的g(z)，grad\_ascent()是梯度下降的主要代码，classify()是跟新样本分类的函数，当预测值大于0.5时，分为类1，否则分为类0。参数的含义在注释中已经介绍了。上述代码中的权重向量weights是一个列向量，而之前理论部分W  W 是一个行向量，在代码中只需稍微变形就可以了。

接下来导入数据，代码如下：

# 获取测试数据
    def load_data():
        data = pd.read_csv("./Data/testSet.txt", sep='\t', header=None)
        data.insert(0, 'Ones', 1)  # 第0行插入1
        data_mat = np.array(data.iloc[:, 0: -1])
        class_label = np.array(data.iloc[:, -1]).reshape(-1, 1)  # 转成列向量
        return data_mat, class_label

在平面上显示为如下图：  
![这里写图片描述][70 2]  
由于数据只有2个特征，所以我们需要往数据中加入偏置量，即之前的x0 x 0,这个在上述代码块的第4行  
接下来我们用梯度上升算法跟新权重W  W ，并画出跟新后的决策边界，代码如下：

# 画决策边界
    def plot_best_fit(data_mat, class_label, weights):
        data_zero = data_mat[class_label.ravel() == 0, :]  # 类标签为0的数据
        data_one = data_mat[class_label.ravel() == 1, :]  # 类标签为1的数据
        plt.scatter(data_zero[:, 1], data_zero[:, 2], c='R', alpha=0.5)
        plt.scatter(data_one[:, 1], data_one[:, 2], c='G', alpha=0.5)
        x = np.linspace(-3, 3, 2)  # -3与3之间通过2个点作出决策边界
        y = (-weights[0] - weights[1] * x)/weights[2]
        plt.plot(x, y)
        plt.xlabel("X1")
        plt.ylabel("X2")
        plt.show()
    
    
    data_mat, class_label = load_data()
    weights = lr.grad_ascent(data_mat, class_label)  # 梯度上升算法
    print(weights)
    plot_best_fit(data_mat, class_label, weights)

结果如下：  
![这里写图片描述][70 3]

以上过程中，学习速率默认为0.001，进行了500次迭代。最后，我们可以看下学习曲线，也就是代价函数的情况。之前没有提到过代价函数，这里理解为−1mL(W) −1 mL ( W )即可，它表示预测类别和真实类别的差异，其值越小越好，其图像如下：  
![这里写图片描述][70 4]

在书上还介绍了一种随机梯度上升算法，代码如下：

# 随机梯度上升算法
    def stoc_grad_ascent(data_mat, class_label, iter_num=500):
        m, n = data_mat.shape
        weights = np.ones(n)
        for j in range(iter_num):
            data_index = [i for i in range(m)]  # 样本的索引值
            for i in range(m):
                alpha = 4/(1.0 + j + i) + 0.01
                rand_index = random.randint(0, len(data_index) - 1)  # 随机索引值
                h = sigmoid(np.dot(data_mat[rand_index], weights))
                error = class_label[rand_index] - h
                weights += alpha * error * data_mat[rand_index]
                del data_index[rand_index]
        return weights

其原理就是每次只用一个样本来进行梯度上升，跟新权重，这里不再讨论。

以上理论部分一、二都是书上没有的，属本人总结归纳，公式的推导部分是我用markdown格式的文本打出来的，应该算写的比较详细了。希望对大家能有帮助。源码和数据集在开头已经标出，欢迎下载和点星星。。。

[https_github.com_TimePickerWang_MachineLearningInAction]: https://github.com/TimePickerWang/MachineLearningInAction
[70]: /images/20220523/a306fbdcaa2348719f47bdea853a7898.png
[70 1]: /images/20220523/883de3206d6843b98e9153249b6c5cd7.png
[70 2]: /images/20220523/d0aa6e60ce8c4d5690e09223fbd5b451.png
[70 3]: /images/20220523/46670d295841411eafc0d29efcc6d969.png
[70 4]: /images/20220523/f530dcd3684d4b9ea61dd6c29a65aff8.png