逻辑回归案例练习

2023-11-11

逻辑回归

场景一

在训练的初始阶段，我们将要构建一个逻辑回归模型来预测，某个学生是否被大学录取。设想你是大学相关部分的管理者，想通过查看申请学生的两次测试的评分，来决定他们是否被录取。现在你拥有之前申请学生的可以用于训练逻辑回归的训练样本集。对于每一个训练样本，你有他们两次测试的评分和最后是被录取的结果。为了完成这个预测任务，我们准备构建一个可以基于两次测试评分来评估录取可能性的分类模型。

1.1 查看数据

导包：

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

读取数据：

path = './data/ex2data1.txt'
data = pd.read_csv(path, header=None, names=['Exam 1', 'Exam 2', 'Admitted'])
data.head()

我们接下来创建两个分数的散点图，如果样本是正的（表示被接纳）或负的（表示未被接纳）：

positive = data[data['Admitted'].isin([1])]
negative = data[data['Admitted'].isin([0])]

fig, ax = plt.subplots(figsize=(12,8))
ax.scatter(positive['Exam 1'], positive['Exam 2'], s=50, c='b', marker='o', label='Admitted')
ax.scatter(negative['Exam 1'], negative['Exam 2'], s=50, c='r', marker='x', label='Not Admitted')
ax.legend()
ax.set_xlabel('Exam 1 Score')
ax.set_ylabel('Exam 2 Score')
plt.show()

看起来在两类间，有一个清晰的决策边界。现在我们需要实现逻辑回归，那样就可以训练一个模型来预测结果。

1.2 sigmoid函数

sigmoid函数表达式：
g ( z ) = 1 1 + e x p − z g(z) = \frac{1}{1 + exp^{-z}} g(z)=1+exp−z1
结合 h θ ( x ) = g ( θ T X ) h_{\theta}(x) = g(\theta^TX) hθ(x)=g(θTX) 可得到最终逻辑回归模型的假设函数：
h θ ( x ) = 1 1 + e x p − θ T X h_{\theta}(x) = \frac{1}{1 + exp^{-\theta^TX}} hθ(x)=1+exp−θTX1

使用Python代码实现：

def fn_sigmoid(z):
	return 1 / (1 + np.exp(-z))

让我们做一个快速的测试，来确保它可以工作：

x = np.linspace(-10, 10, 10000)

fig, ax = plt.subplots(figsize=(12,8))
ax.plot(x, sigmoid(x), 'r')
plt.show()

没错，这就是“大S”曲线！现在，我们需要编写代价函数来评估结果。代价函数：

def cost(theta, X, y):
   theta = np.matrix(theta)
    X = np.matrix(X)
    y = np.matrix(y)
    y0 = np.multiply(-y, np.log(fn_sigmoid(X * theta.T)))
    y1 = np.multiply(-(1 - y), np.log(1 - fn_sigmoid(X * theta.T)))
    return np.sum(y0 + y1) / len(X)

现在，我们要对数据做一些预处理，确保数据在交给我们设计的功能函数进行计算时能够正常工作：

# 1.给原数据添加一列Ones，让矩阵乘法更加简单
data.insert(0, 'Ones', 1)

data.head()

# 2.找出特征值X和目标值y
cols = data.shape[1]
X = data.iloc[:,0:cols-1]
y = data.iloc[:,cols-1:cols]  


# 3.将特征值和目标值转换成数组，并且初始化参数θ
X = np.array(X.values)
y = np.array(y.values)
theta = np.zeros(3)

查看特征值、目标值以及初始化的 θ \theta θ 参数：

(X.shape, y.shape, theta.shape)  
# result  
((100, 3), (100, 1), (3,))

接下来让我们计算初始化参数的代价函数( θ \theta θ 的值全为 0)：

cost(theta, X, y)

# result
0.6931471805599453

看起来还不错，接下来我们需要一个函数来计算我们的训练数据、标签和参数 θ \theta θ 的梯度。

1.3 梯度下降（Gradient Descent）

这里我们使用批量梯度下降（Batch Gradient Descent），将批量的数据向量化后再计算：

def gradient(theta, X, y):
    theta = np.matrix(theta)
    X = np.matrix(X)
    y = np.matrix(y)

    parameters = int(theta.ravel().shape[1])
    grad = np.zeros(parameters)

    error = fn_sigmoid(X * theta.T) - y

    for i in range(parameters):
        term = np.multiply(error, X[:,i])
        grad[i] = np.sum(term) / len(X)

    return grad

注意，我们实际上没有在这个函数中执行梯度下降，我们仅仅在计算一个梯度步长。我们可以用 Scipy 的 optimize 命名空间来计算损失和梯度。

我们看看用我们的数据和初始参数 θ \theta θ 为 0 的梯度下降法的结果：

gradient(theta, X, y)

# result
array([ -0.1       , -12.00921659, -11.26284221])

现在可以用 Scipy’s truncated newton（TNC）实现寻找最优化后的参数 θ \theta θ ：

import scipy.optimize as opt
result = opt.fmin_tnc(func=cost, x0=theta, fprime=gradient, args=(X, y))
result  

# result  
(array([-25.16131872,   0.20623159,   0.20147149]), 36, 0)

让我们看看在这个结论下代价函数计算结果是多少：

cost(result[0], X, y)  

# result  
0.20349770158947425

接下来我们需要编写一个函数，用我们所学的参数 θ \theta θ 来为数据集 X X X 输出预测结果。然后，我们可以使用这个函数来给我们的分类器的训练精度打分。逻辑回归模型的假设函数为：
h θ ( x ) = 1 1 + e x p − θ T X h_{\theta}(x) = \frac{1}{1 + exp^{-\theta^TX}} hθ(x)=1+exp−θTX1
当 h θ ( x ) h_{\theta}(x) hθ(x) 大于等于 0.5 时，预测 y=1 ；

当 h θ ( x ) h_{\theta}(x) hθ(x) 小于 0.5 时，预测 y=0 。

def predict(theta, X):
	theta = np.matrix(theta)
    probability = fn_sigmoid(X * theta.T)
    return [1 if x >= 0.5 else 0 for x in probability]

查看准确率：

predictions = predict(result[0], X)
correct = [1 if a == b else 0 for (a, b) in zip(predictions, y)]
accuracy = (sum(map(int, correct)) % len(correct))
print ('accuracy = {0}%'.format(accuracy))

# result
accuracy = 89%

如果使用现在的逻辑回归模型预测一个学生被录取或没有录取，能够达到 89%的精确度。不差！记住，这只是训练集的准确性。我们没有使用交叉验证得到的模型更加准可靠，所以这个准确率有可能高于其真实值。

带正则化的逻辑回归

场景二

1.4 正则化

我们将要通过加入正则项提升逻辑回归算法。正则化是成本函数中的一个术语，它使算法更倾向于选择“更简单”的模型（在这种情况下，模型将会有值更小的系数）。这个理论（类似于奥卡姆剃刀）有助于减少过拟合，提高模型的泛化能力。

设想你是工厂的生产主管，你有一些芯片在两次测试中的测试结果。对于这两次测试，你想决定芯片是否要被接受或抛弃。为了帮助你做出艰难的决定，你拥有过去芯片的测试数据集，从其中你可以构建一个逻辑回归模型。

老规矩，我们先从数据的可视化开始：

path =  './data/ex2data2.txt'
data2 = pd.read_csv(path, header=None, names=['Test 1', 'Test 2', 'Accepted'])
data2.head()

绘制样本的散点图，区分接受或者被抛弃的样本：

positive = data2[data2['Accepted'].isin([1])]
negative = data2[data2['Accepted'].isin([0])]

fig, ax = plt.subplots(figsize=(12,8))
ax.scatter(positive['Test 1'], positive['Test 2'], s=50, c='b', marker='o', label='Accepted')
ax.scatter(negative['Test 1'], negative['Test 2'], s=50, c='r', marker='x', label='Rejected')
ax.legend()
ax.set_xlabel('Test 1 Score')
ax.set_ylabel('Test 2 Score')
plt.show()

这个数据看起来比前一次的复杂得多。特别地，你会注意到我们找不到线性决策边界来良好的分开两类数据。一个方法是用像逻辑回归这样的线性技术来构造从原始特征的多项式中得到的特征。让我们通过创建一组多项式特征开始入手：

degree = 5
x1 = data2['Test 1']
x2 = data2['Test 2']

data2.insert(3, 'Ones', 1)

for i in range(1, degree):
    for j in range(0, i):
        data2['F' + str(i) + str(j)] = np.power(x1, i-j) * np.power(x2, j)

data2.drop('Test 1', axis=1, inplace=True)
data2.drop('Test 2', axis=1, inplace=True)

data2.head()

现在，我们需要修改上面的代价和梯度函数，并添加正则化项即可。首先是代价函数：

# 正则化代价函数
def cost_with_reg(theta, X, y, learning_rate):
    theta = np.matrix(theta)
    X = np.matrix(X)
    y = np.matrix(y)
    y0 = np.multiply(-y, np.log(fn_sigmoid(X * theta.T)))
    y1 = np.multiply(-(1 - y), np.log(1 - fn_sigmoid(X * theta.T)))
    reg = (learning_rate / (2 * len(X))) * np.sum(np.power(theta[:,1:theta.shape[1]], 2))
    return np.sum(y0 + y1) / len(X) + reg

请注意等式中的 reg 项。还需要注意到另外的一个参数 学习率 。这是一种超参数，用来控制正则化项。现在我们需要设计正则化梯度函数：

如果我们要使用梯度下降法令这个代价函数最小化，因为我们未对 θ 0 \theta_0 θ0 进行正则化，所以梯度下降算法将分如下两种情形：
repeat until convergence { \text{repeat until convergence \{} repeat until convergence {
θ 0 : = θ 0 − α 1 m ∑ i = 1 m ( h θ ( x ( i ) ) − y ( i ) ) x 0 ( i ) , \theta_0 := \theta_0 - \alpha \frac{1}{m} \sum_{i=1}^m \bigg(h_{\theta}(x^{(i)}) - y^{(i)}\bigg)x_{0}^{(i)}, θ0:=θ0−αm1i=1∑m(hθ(x(i))−y(i))x0(i),
θ j : = θ j − α [ 1 m ∑ i = 1 m ( h θ ( x ( i ) ) − y ( i ) ) x j ( i ) + λ m θ j ] , \theta_j := \theta_j - \alpha \bigg[\frac{1}{m} \sum_{i=1}^m \bigg(h_{\theta}(x^{(i)}) - y^{(i)}\bigg)x_{j}^{(i)} + \frac{\lambda}{m} \theta_j \bigg], θj:=θj−α[m1i=1∑m(hθ(x(i))−y(i))xj(i)+mλθj],
for j = 1 , 2 , . . . , n } \text{for $j = 1, 2, ..., n$ \}} for j=1,2,...,n }

带正则化项的梯度下降函数：

def gradient_with_reg(theta, X, y, learning_rate):
    theta = np.matrix(theta)
    X = np.matrix(X)
    y = np.matrix(y)

    parameters = int(theta.ravel().shape[1])
    grad = np.zeros(parameters)

    error = fn_sigmoid(X * theta.T) - y

    for i in range(parameters):
        term = np.multiply(error, X[:,i])

        if (i == 0):
            grad[i] = np.sum(term) / len(X)
        else:
            grad[i] = (np.sum(term) / len(X)) + ((learning_rate / len(X)) * theta[:,i])

    return grad

找出特征值和目标值，并初始化参数 θ \theta θ ：

# 找特征值X和目标值y
cols = data2.shape[1]
X2 = data2.iloc[:,1:cols]
y2 = data2.iloc[:,0:1]

# 转换成numpy数组
X2 = np.array(X2.values)
y2 = np.array(y2.values)
theta2 = np.zeros(11)

给学习率初始化一个合理值，如果有必要的话（即如果惩罚太强或不够强）,我们可以之后再对这个值进行微调：

learning_rate = 1

接下来我们尝试调用新的默认值全为 0 的 θ \theta θ 的正则化函数，以确保其计算工作是正常的：

cost_with_reg(theta2, X2, y2, learning_rate)

# result
0.6931471805599454  

gradient_with_reg(theta2, X2, y2, learning_rate)

# result
array([0.00847458, 0.01878809, 0.05034464, 0.01150133, 0.01835599,
       0.00732393, 0.00819244, 0.03934862, 0.00223924, 0.01286005,
       0.00309594])

现在我们可以使用和上面相同的优化函数来计算优化后的结果：

result2 = opt.fmin_tnc(func=cost_with_reg, x0=theta2, fprime=gradient_with_reg, args=(X2, y2, learning_rate))
result2

# result
(array([ 0.53010248,  0.29075567, -1.60725764, -0.5821382 ,  0.01781027,
        -0.21329508, -0.40024142, -1.37144139,  0.02264303, -0.9503358 ,
         0.0344085 ]), 22, 1)

最后，同样我们可以使用上面的预测函数来查看带正则化的逻辑回归模型在训练数据上的准确度：

predictions = predict(result2[0], X2)
correct = [1 if a==b else 0 for (a, b) in zip(predictions, y2)]
accuracy = (sum(map(int, correct)) % len(correct))
print ('accuracy = {0}%'.format(accuracy))  

# result
accuracy = 78%

虽然我们实现了这些算法，值得注意的是，我们还可以使用 Python 的高级库 scikit-learn 来解决这个问题：

# 使用scikit-learn来建立逻辑回归模型
from sklearn import linear_model
model = linear_model.LogisticRegression(penalty='l2', C=1.0)
model.fit(X2, y2.ravel())

# result
LogisticRegression(C=1.0, class_weight=None, dual=False, fit_intercept=True,
          intercept_scaling=1, max_iter=100, multi_class='ovr', n_jobs=1,
          penalty='l2', random_state=None, solver='liblinear', tol=0.0001,
          verbose=0, warm_start=False)

计算准确率：

model.score(X2, y2)  

# result
0.6610169491525424

这个准确度和我们刚刚实现的差了很多，不过请记住，这个结果可以使用默认参数来计算。我们可能需要做一些参数的调整来获得和我们之前相同的结果。

本文内容由网友自发贡献，版权归原作者所有，本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容，请联系:hwhale#tublm.com(使用前将#替换为@)

机器学习

监督学习