Press "Enter" to skip to content

机器学习 – 线性回归与逻辑回归(实践部分)

本站内容均来自兴趣收集,如不慎侵害的您的相关权益,请留言告知,我们将尽快删除.谢谢.

之前对线性回归和逻辑回归的理论部分做了较为详细的论述,下面通过一些例子再来巩固一下之前所学的内容。

 

需要说明的是,虽然我们在线性回归中都是直接通过公式推导求出w和b的精确值,但在实际运用中基本上都会采用梯度下降法作为首选,因为用代码表示公式会比较繁琐,而梯度下降法只需要不断对参数更新公式进行迭代即可,用代码表示非常简单,且迭代次数基本上不会对性能有太大的拖累,所以下面我们也将采用梯度下降的方法。

 

使用线性回归模拟数据分布

 

了解数据初始分布情况

 

现有一个 数据集 (提取码:yim9),描述了一系列点的分布情况。我们先将初始的分布情况绘制出来,便于对数据有个整体的了解:

 

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from sklearn.linear_model import LinearRegression
data = np.loadtxt('./data/linear_regression_data.txt', delimiter=',')
print(data[:10])  # 查看部分样本数据
data.shape        # 查看样本规模
# 将数据集转换成 y = w*x 的向量形式
x = np.c_[np.ones(data.shape[0]), data[:, 0]]  # 使用np.c_转换为列向量
y = np.c_[data[:, 1]]  # 1是给b留的位置
# 展示样本分布情况
plt.scatter(x[:, 1], y, color='green', marker='x', linewidth=1)
plt.xlabel('X')
plt.ylabel('Y')

 

可以看到,点在X轴左侧分布得比较密集,但是越往右,越呈现出按照直线分布的趋势,所以可以尝试使用线性回归进行模拟。

 

使用最小二乘法构造目标函数

 

我们需要求出目标函数,以便求出梯度,并记录目标函数在每轮迭代中的值

 

# 目标函数
def target_func(x, y, w=[[0],[0]]):
    # 这里我们设置w和b的初始值为0
    h = x.dot(w)                                   # h表示预测值
    m = y.size                                     # 通常会使用样本大小m平均整体误差的值
    target = 1.0/(2*m) * (np.sum(np.square(h-y)))  # 使用最小二乘法得到目标函数
    return target

 

计算梯度及参数更新公式

 

# 计算梯度
def grad(x, y, w):
    h = x.dot(w)
    m = y.size
    grad = (1.0/m)*(x.T.dot(h-y))  # 对目标函数求导即是梯度
    return grad
# 使用梯度循环更新参数
def grad_descent(x, y, nums_iter=1000, n=0.01, w=[[0],[0]]):
    # 默认迭代1000次,学习率为0.01
    target_his = np.zeros(nums_iter)    # 存放目标函数每轮迭代后的值
    for i in np.arange(nums_iter):
        # 循环更新w,并记录目标函数的值
        w = w - n*grad(x, y, w)
        target_his[i] = target_func(x, y, w)
    return w, target_his
w, target_his = grad_descent(x, y)
# 展示迭代效果
plt.plot(target_his)
plt.xlabel('Iter Times')
plt.ylabel('Loss Function')

可以看出,迭代1000次,性能逐渐趋于平稳。

 

最后与sklearn的结果相对比:

 

# 样本
plt.scatter(x[:, 1], y, color='red', marker='x', linewidth=1)
# 自主
xx = np.arange(5, 23)
yy = w[0] + xx*w[1]
plt.plot(xx, yy, label='LD(GR)')
# sklearn
regr = LinearRegression()
# reshape转换成列向量
regr.fit(x[:, 1].reshape(-1, 1), y.ravel())
yy2 = regr.intercept_ + regr.coef_ * xx
plt.plot(xx, yy2, label='LD(sklearn)')
plt.xlim(4, 24)
plt.xlabel('X')
plt.ylabel('Y')
plt.legend(loc=4)

 

可以看出效果基本一致。

 

使用逻辑回归分类

 

现有一个 考试成绩数据集 (提取码:7n04),该数据集为某班级英语和数学两个科目的考试成绩,前两列分别表示英语和数学成绩,最后一列表示最终是否通过。

 

加载数据,查看数据样本及可视化

 

import pandas as pd
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
from scipy.optimize import minimize
# Step 1. 加载数据
def loaddata(data):
     return np.loadtxt(data, delimiter=',')
data = loaddata('./data/sample_1.txt')
# 查看数据
print(data.shape)
print(data[:5, :])
# 运行结果:
(100, 3)
[[34.62365962 78.02469282  0.        ]
 [30.28671077 43.89499752  0.        ]
 [35.84740877 72.90219803  0.        ]
 [60.18259939 86.3085521   1.        ]
 [79.03273605 75.34437644  1.        ]]

 

# Step 2. 整理训练样本
# np.c_是按行连接两个矩阵,就是把两矩阵左右相加,要求行数相等
X = np.c_[np.ones((data.shape[0],1)), data[:,0:2]]
y = np.c_[data[:,2]]
print(X[:5])
print(X.shape)
print(y[:5])
# 运行结果:
[[ 1.         34.62365962 78.02469282]
 [ 1.         30.28671077 43.89499752]
 [ 1.         35.84740877 72.90219803]
 [ 1.         60.18259939 86.3085521 ]
 [ 1.         79.03273605 75.34437644]]
(100, 3)
[[0.]
 [0.]
 [0.]
 [1.]
 [1.]]

 

X的第一列是给偏置b预留的,其余两列为特征值。

 

# Step 3. 可视化
def plotData(data, label_x, label_y, label_pos, label_neg, axes=None):
    # 获得正负样本的下标(正样本表示考试通过的,负样本表示考试未通过)
    neg = data[:,2] == 0
    pos = data[:,2] == 1
    if axes == None:
        axes = plt.gca() # 获取图的坐标信息
    # 依次绘制正负样本,使用不同颜色区分
    axes.scatter(data[pos][:,0], data[pos][:,1], marker='*', c='r', s=30, linewidth=2, label=label_pos)
    axes.scatter(data[neg][:,0], data[neg][:,1], c='c', s=30, label=label_neg)
    axes.set_xlabel(label_x)
    axes.set_ylabel(label_y)
    axes.legend(frameon= True, fancybox = True);
    
plotData(data, 'English', 'Math', 'Pass', 'Fail')

 

从样本分布可以看出,正负样本分布的还算整齐,可以使用决策边界进行划分。

 

定义Sigmoid函数、损失函数(即目标函数)、梯度函数

 

# Step 4. 训练前准备
# Sigmoid函数
def sigmoid(z):
    return(1 / (1 + np.exp(-z)))
# 损失函数
def lossFunction(w, X, y):
    m = y.size
    h = sigmoid(X.dot(w.reshape(-1, 1)))
    # 套公式得到损失函数
    J = -1.0*(1.0/m)*(np.log(h).T.dot(y)+np.log(1-h).T.dot(1-y))
    if np.isnan(J[0]):
        return(np.inf)
    return J[0]
# 梯度计算函数
def gradient(X, y, w):
    grad = np.zeros(w.shape)
    error = (sigmoid(X.dot(w.reshape(-1, 1))).reshape(-1, 1) - y).ravel()
    for j in range(len(w.ravel())):
        # 循环更新w,得到梯度
        term = np.multiply(error, X[:,j])
        grad[j] = np.sum(term) / len(X)
    return grad
# 初始化w
initial_w = np.zeros(X.shape[1])
loss = lossFunction(initial_w, X, y)
grad = gradient(X, y, initial_w)
print('Loss: 
', loss)
print('Grad: 
', grad)
# 输出结果:
Loss: 
 [0.69314718]
Grad: 
 [ -0.1        -12.00921659 -11.26284221]

 

关于损失函数,代码里使用的是下面这个公式,而不是之前推导到最后的公式,因为用代码更好实现:

使用梯度下降法训练

 

# step 5. 使用梯度下降法训练
# 梯度下降函数
def gradientDescent(X, y, w, alpha=0.01, num_iters=1000):
    m = y.size
    J_history = np.zeros(num_iters)
    for iter in np.arange(num_iters):
        h = X.dot(w.reshape(-1, 1))
        grad = gradient(X, y, w)
        w = w - alpha*(1.0/m)*grad
        J_history[iter] = lossFunction(w, X, y)
    return(w, J_history)
# 画出每一次迭代和损失函数变化
w , Cost_J = gradientDescent(X, y, initial_w, alpha=0.15, num_iters=100000)
plt.plot(Cost_J)
plt.ylabel('Loss J')
plt.xlabel('Iter')

通过将学习率设置为0.15,并经过10W次的迭代后,损失函数计算出来的误差控制在了0.32左右。

 

下面对迭代出来的w进行可视化。

 

可视化决策边界

 

# Step 6. 可视化
plotData(data, 'English', 'Math', 'Pass', 'Fail')
# np.linspace(x1_min, x1_max) 生成的是x轴英语成绩的列向量
# np.linspace(x2_min, x2_max) 生成的是y轴数学成绩的列向量
# xx1, xx2为meshgrid生成的坐标矩阵
x1_min, x1_max = X[:,1].min(), X[:,1].max(),
x2_min, x2_max = X[:,2].min(), X[:,2].max(),
xx1, xx2 = np.meshgrid(np.linspace(x1_min, x1_max), np.linspace(x2_min, x2_max))
# 将预测值作为等高线绘制出来,即决策边界
# 将xx1和xx2作为x的特征值,得到预测值h
h = sigmoid(np.c_[np.ones((xx1.ravel().shape[0],1)), xx1.ravel(), xx2.ravel()].dot(w.T))
h = h.reshape(xx1.shape)
plt.contour(xx1, xx2, h, 1, linewidths=1, colors='b');

 

 

使用scipy.optimize.minimize实现

 

scipy.optimize.minimize可以根据现有的损失函数和梯度,优化整个梯度下降的过程,可以自动选择学习率,从而免去了不断手动调试的麻烦。

 

from scipy.optimize import minimize
res = minimize(lossFunction, initial_w, args=(X,y), jac=gradient, options={'maxiter':400})
print(res)

 

传入相应的参数,并设置迭代轮数(这里仅设置了400次),就可得到优化后的w。从输出可以看到,损失函数的误差已经减小到了0.2,比之前手动实现的要小很多,而且更加省时。

 

# Step 5. 可视化
plotData(data, 'English', 'Math', 'Pass', 'Fail')
x1_min, x1_max = X[:,1].min(), X[:,1].max(),
x2_min, x2_max = X[:,2].min(), X[:,2].max(),
xx1, xx2 = np.meshgrid(np.linspace(x1_min, x1_max), np.linspace(x2_min, x2_max))
h = sigmoid(np.c_[np.ones((xx1.ravel().shape[0],1)), xx1.ravel(), xx2.ravel()].dot(res.x.T))
h = h.reshape(xx1.shape)
plt.contour(xx1, xx2, h, 1, linewidths=1, colors='b');

 

从训练结果来看,手动和调包实现的效果好像差不多,但是明显后者更加省时省力。

Be First to Comment

发表评论

您的电子邮箱地址不会被公开。 必填项已用*标注