机器学习之逻辑回归(Logistic Regression)

最新推荐文章于 2024-08-09 12:39:33 发布

David_Hdw

最新推荐文章于 2024-08-09 12:39:33 发布

阅读量2.8k

点赞数 4

分类专栏：机器学习文章标签：机器学习逻辑回归 logistic regression

本文链接：https://blog.csdn.net/perfect1t/article/details/82966647

版权

机器学习专栏收录该内容

12 篇文章 14 订阅

订阅专栏

一、逻辑回归模型

逻辑回归主要应用于二分类问题，其主要思想是：根据现有数据对分类边界线建立回归公式，以此进行分类。逻辑回归把线性回归 $\textup{z}=w^Tx+b$ 的输出集输入到simoid函数

$y=\frac{1}{1+e^{-\textup{z}}}$

中得到：

$y=\frac{1}{1+e^{-(w^Tx+b)}}$ .

sigmoid函数自变量为实数集R,即自变量的取值范围在负无穷到正无穷之间，值域在0到1之间，函数图像如下：

从图像可以看出来，它把z值转化为接近0或1的值，并且其输出值在z=0附近变化很陡峭。

由公式：

$y=\frac{1}{1+e^{-(w^Tx+b)}}$

我们可以得到：

$ln\frac{y}{1-y}=w^Tx+b$

若将y视为正例的概率，则1-y为反例的概率。则两者的比值

$\frac{y}{1-y}$

称为几率(odds)，反应了作为正例的相对可能性。对几率取对数称为对数几率，也叫logit

$ln\frac{y}{1-y}$ ，

从下面的公式：

$ln\frac{y}{1-y}=w^Tx+b$

我们可以看出，实际上在利用线性回归模型的预测结果去逼近真实标记的对数几率，所以称为逻辑回归模型(logistic regression)，模型为：

$y=\frac{1}{1+e^{-(w^Tx+b)}}$ 。

由 $y=\frac{1}{1+e^{-(w^Tx+b)}}$ ，我们可以得到： $ln\frac{y}{1-y}=w^Tx+b$ 的推导过程：

$y=\frac{1}{1+e^{-(w^Tx+b)}}====>y=\frac{1}{\frac{e^{w^Tx+b}+1}{e^{w^Tx+b}}}$

$y=\frac{e^{w^Tx+b}}{1+e^{w^Tx+b}} = ==>ye^{w^Tx+b}+y=e^{w^Tx+b}$

$y=e^{w^Tx+b}(1-y)====>ln\frac{y}{1-y}=w^Tx+b$

逻辑斯蒂回归的优点：

它直接对分类可能性进行建模，无需事先假设数据分布，避免了假设分布不准确带来的问题。
不仅可以预测出类别，而是可得到近似概率预测，对许多需利用概率辅助决策的任务很有用。
对数几率函数是任意阶可导的凸函数，可以使用许多现有的优化方法进行求解最优解。

二、LR的损失函数 ---- logloss

通过上面，我们知道logistics Regression的模型为：

$y=\frac{1}{1+e^{-(w^Tx+b)}}$

那我们如何求解参数w和b呢？

如上图，显然该二分类问题的分类边界为黑色的粗线，那么打叉的数据为大于零的，即正样本，而圆圈的数据为负样本，分布在分类边界的下方，即带入分类边界函数，结果小于零，那么我们把某个样本带入分类边界函数的结果带进sigmoid函数，有什么样的结果呢？发现当为负样本时，即带入sigmoid的结果为小于0.5，而正样本带入则结果大于0.5，所以我们可以认为集输出的结果为该样本属于正样本的概率，当大于0.5时为正样本，小于0.5时为负样本，当然该阀值是可以调整的。

logLoss (log 损失函数)的函数标准形式：

$L(Y,\ P(Y|X)) = - log\ P(Y|X)$

逻辑回归的P(Y=y|x)的表达式如下：

$P(Y=y|x) = \left\{\begin{matrix} \ \ \ \ \ \ \ h_\Theta (x) = g(f(x)) = \frac{1}{1+e^{-f(x)}}, \ \ \ \ y=1\ \ \ \ \\\\ 1 - h_\Theta (x) = 1 - g(f(x)) = \frac{1}{1+e^{f(x)}}, \ \ y = 0 \end{matrix}\right.$

所以我们可以得到逻辑(斯蒂)回归的损失函数为：

对m个样本有：

$J(\theta)=\frac{1}{m}\sum_{i=1}^{m}Cost(h_{\theta}(x^{(i)}),y^{(i)})$

$=-\frac{1}{m}\sum_{i=1}^{m}[y^{(i)}log\ h_{\theta}(x^{(i)})+(1-y^{(i)})log\ (1-h_{\theta}(x^{(i)}))]$

别忘了正则化项目,，所以损失函数logloss（也叫binary cross entropy，二元交叉熵损失）为：

$J(\theta)=-\frac{1}{m}\sum_{i=1}^{m}[y^{(i)}log\ h_{\theta}(x^{(i)})+(1-y^{(i)})log\ (1-h_{\theta}(x^{(i)}))]\ +\ \frac{\lambda }{2m}\sum_{j=1}^{n}\theta_j^2$

这里先不求解，因为后面我们会发现，该损失函数只是多了一个正则化项，其他与下文通过极大似然方法得到的函数式一模一样。

一般我们采用极大似然法来估计参数w和b，则对于数据集

$D=\left \{ (x_1,y_1),(x_2,y_2),...,(x_m,y_m) \right \}$

LR模型的似然函数为：

$L(w,b)=\prod_{i=1}^{m}p(y_i|x_i;w,b)$

其中 p(y_i|x_i;w,b) 表示给定w,b和 x_i 的情况下，模型预测出来为 y_i 的可能性。所以其对数似然为：

$log(L(w,b))=\sum_{i=1}^{m}log(p(y_i|x_i;w,b))$

也就是：

log(L(w,b))=log(p(Y|X;w,b))

其中，X,Y为数据集，而我们的目的就是求似然的最大，所以相当于求负对数似然函数的最小值，也就是使得模型的损失最小化。所以我们可以得到logistics Regression的损失函数为：

L(Y,P(Y|X))=-log(P(Y|X))

对于一个样本(x,y)，逻辑回归的P(Y=y|x)表达式如下：

$P(Y=y|x)=\left\{\begin{matrix} h_\Theta (x)=g(f(x))=\frac{1}{1+exp(-f(x))},\ \ \ \ \ \ \ \ \ \ \ y=1\\ \\ 1-h_\Theta (x)=1-g(f(x))=\frac{1}{1+exp(f(x))},\ \ \ y=0 \end{matrix}\right.$

2、另外一种理解方法:

三、实战

数据在博主github上：https://github.com/davidHdw/machine-learning

# -*- coding: utf-8 -*-
"""
Created on Mon Oct  8 17:25:33 2018

@author: 87955
"""

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

def loadDataSet(filename):
    '''
    函数说明：读取数据
    参数：
        filename：文件名
    返回：
        pdData : 数据集
    '''
    pdData = pd.read_csv(filename, header=None, names=['Exam 1','Exam 2','Admitted'])
    return pdData

def plotDataSet(pdData):
    '''
    函数说明：对数据集进行可视化
    参数：pdData : DataFrame类型数据集
    '''
    positive = pdData[pdData['Admitted'] == 1]  # 正例
    negative = pdData[pdData['Admitted'] == 0]  # 反例
    
    fig, ax = plt.subplots(figsize=(10, 5))
    ax.scatter(positive['Exam 1'], positive['Exam 2'], s=30, c='b', marker='o', label='Admitted')
    ax.scatter(negative['Exam 1'], negative['Exam 2'], s=30, c='r', marker='x', label='No Admitted')
    ax.legend()
    ax.set_xlabel('Exam 1 Score')
    ax.set_ylabel("Exam 2 Score")
    
    
def sigmoid(z):
    '''
    函数说明：sigmoid函数
    '''
    return (1 / (1+np.exp(-z)))

def logisticModel(X, ws):
    '''
    函数说明：逻辑斯蒂回归模型
    参数：
        X:  数据集（不包括标签）
        ws: 回归参数
    返回：
        当前预测结果
    '''
    return sigmoid(np.dot(X, ws.T))
    
def cost(X, y, ws):
    '''
    函数说明：求平均损失函数
    参数：
        X: 数据集（不含标签）
        y: 标签
    返回：
        平均损失
    '''
    left = np.multiply(-y, np.log(logisticModel(X, ws)))
    right = np.multiply(1-y, np.log(1 - logisticModel(X, ws)))
    return np.sum(left - right) / (len(X))

def gradient(X, y, ws):
    '''
    函数说明：梯度函数
    参数：
        X: 数据集（不含标签）
        y: 标签        
    返回：
        梯度值
    '''
    grad = np.zeros(ws.shape)   
    error = (logisticModel(X, ws) - y).ravel() # 误差
    for j in range(len(ws.ravel())):
        term = np.multiply(error, X[:,j])
        grad[0, j] = np.sum(term) / len(X)
        
    return grad

STOP_ITER = 0
STOP_COST = 1
STOP_GRAD = 2

def stopCriterion(tp, value, threshold):
    '''
    设置三种不同的停止策略
    '''
    if tp == STOP_ITER:
        return value > threshold
    elif tp == STOP_COST:
        return abs(value[-1]-value[-2]) < threshold
    elif tp == STOP_GRAD:
        return np.linalg.norm(value) < threshold
    
def shuffleData(data):
    '''
    函数说明：对数据进行洗牌
    '''
    np.random.shuffle(data)
    cols = data.shape[1]
    X = data[:, 0:cols-1]
    y = data[:, cols-1:]
    return X,y

import time

def graDescent(data, ws, batchSize, stopType, thresh, alpha):
    '''
    函数说明：梯度下降求解
    参数：
        data：数据集
        ws：回归参数
        batchSize：大小
        stopType：停止更新模式标志
        thresh：阀值
        alpha： 学习率
    返回：
    '''
    init_time = time.time()
    i = 0  # 迭代次数
    k = 0  # batch
    X, y = shuffleData(data)
    m = X.shape[0]           # 样本个数
    grad = np.zeros(ws.shape)  # 初始化梯度
    costs = [cost(X, y, ws)]   # 损失值
    
    while True:
        grad = gradient(X[k:k+batchSize], y[k:k+batchSize], ws)
        k += batchSize
        if k >= m:
            k = 0
            X, y = shuffleData(data)  # 重新洗牌
        ws = ws - alpha * grad        # 参数更新
        costs.append(cost(X,y, ws))   # 计算新的损失
        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):
            break
        
    return ws, i-1, costs, grad, time.time() - init_time

def runExpe(data, ws, batchSize, stopType, thresh, alpha):
    ws, iter, costs, grad, dur = graDescent(data, ws, batchSize, stopType, thresh, alpha)
    name = "Original" if (data[:,1]>2).sum() > 1 else "Scaled"
    name += " data - learning rate: {} - ".format(alpha)
    if batchSize == data.shape[0]:
        strDescType = "Gradient"
    elif batchSize == 1:
        strDescType = "Stochastic"
    else: strDescType = "Mini-batch({})".format(thresh)
    
    name += strDescType + " descent - Stop: "
    if stopType == STOP_ITER:
        strStop = "{} iterations".format(thresh)
    elif stopType == STOP_COST:
        strStop = "costs change < {}".format(thresh)
    else: strStop = "gradient norm < {}".format(thresh)
    
    name += strStop
    print("***{}\nTheta: {} - Iter: {} - Last cost: {:03.2f} - Duration: {:03.2f}"
          .format(name, ws, iter, costs[-1], dur))
    fig, ax = plt.subplots(figsize=(12,4))
    ax.plot(np.arange(len(costs)), costs, 'r')
    ax.set_xlabel("Iterations")
    ax.set_ylabel("Cost")
    ax.set_title(name.upper()+" - Error vs. Iteration")
    return ws

#设定阈值
def predict(X, theta):
    return [1 if x >= 0.5 else 0 for x in logisticModel(X, theta)]



if __name__ == "__main__":
    pdData = loadDataSet('LogiReg_data.txt')
    plotDataSet(pdData)
    pdData.insert(0, "ones", 1) # 插入x0 = 1
    orig_data = pdData.as_matrix() # 把DataFrame转成matrix类型
    cols = orig_data.shape[1]
    X = orig_data[:,0:cols-1]
    y = orig_data[:,cols-1:]
    theta = np.zeros([1,3])
    
    runExpe(orig_data, theta, 100, STOP_ITER, thresh=5000, alpha=0.000001)
    runExpe(orig_data, theta, 100, STOP_COST, thresh=0.000001, alpha=0.001)
    runExpe(orig_data, theta, 100, STOP_GRAD, thresh=0.05, alpha=0.001)
    runExpe(orig_data, theta, 1, STOP_ITER, thresh=5000, alpha=0.001)
    runExpe(orig_data, theta, 16, STOP_ITER, thresh=15000, alpha=0.001)
    
    # 对数据进行预处理：标准化（均值为0，方差为1）
    import sklearn.preprocessing as pp
    scaled_data = orig_data.copy()
    scaled_data[:,1:3] = pp.scale(scaled_data[:, 1:3])
    print(scaled_data[:5])
    
    runExpe(scaled_data, theta, 100, STOP_GRAD, thresh=0.02, alpha=0.001)
    theta = runExpe(scaled_data, theta, 1, STOP_GRAD, thresh=0.002/5, alpha=0.001)
    runExpe(scaled_data, theta, 16, STOP_GRAD, thresh=0.002*2, alpha=0.001)
    
    # 正确率验证
    scaled_X = scaled_data[:, :3]
    y = scaled_data[:, 3]
    predictions = predict(scaled_X, theta)
    correct = [1 if ((a == 1 and b == 1) or (a == 0 and b == 0)) else 0 for (a, b) in zip(predictions, y)]
    accuracy = (sum(map(int, correct)) % len(correct))
    print ('accuracy = {0}%'.format(accuracy))