机器学习练习三：逻辑回归

最新推荐文章于 2021-11-21 23:34:39 发布

lc_1123

最新推荐文章于 2021-11-21 23:34:39 发布

阅读量1.5k

点赞数

分类专栏： Python 机器学习文章标签： python 机器学习

Python 同时被 2 个专栏收录

4 篇文章 0 订阅

订阅专栏

机器学习

3 篇文章 1 订阅

订阅专栏

    <div class="markdown_views"><hr>

转载出处：http://blog.csdn.net/and_w/article/details/52836716

作者：John Wittenauer
翻译：GreatX
源：Machine Learning Exercises In Python, Part 3

这篇文章是一系列 Andrew Ng 在 Coursera 上的机器学习课程的练习的一部分。这篇文章的原始代码，练习文本，数据文件可从这里获得。

Part 1 简单线性回归（Simple Linear Regression）
Part 2 多元线性回归（Multivariate Linear Regression）
Part 3 逻辑回归（Logistic Regression）
Part 4 多元逻辑回归（Multivariate Logistic Regression）
Part 5 神经网络（Neural Networks）
Part 6 支持向量机（Support Vector Machines）
Part 7 K-均值聚类与主成分分析（K-Means Clustering & PCA）
Part 8 异常检测与推荐（Anomaly Detection & Recommendation）

在这个系列的第二部分用梯度下降完成了多元线性回归的实现并且将它应用于房价数据集。在这篇文章中我们将目标从预测连续值（回归）转到将结果分为两个或多个离散桶（ bucket）（分类）并将它应用于学生招生问题。假设你是大学部门的管理者，你想要基于两次考试的结果来确定每个申请人的入学机会。你有之前申请人的历史数据你可以用之作训练集。对于每个训练样本，你有该申请人的两次考试的分数和最终录取决定。为了完此目的，我们将使用逻辑回归建立一个基于考试分数的分类模型（classification model ）以估计录取的可能性。

逻辑回归

你或许会疑惑为什么将“回归”算法用在分类问题上。虽然这个名字似乎指示着什么，但是逻辑回归实际上是分类算法。我猜它之所以这样命名是因为在它的学习方法上和线性回归相似，但是代价（cost）和梯度（gradient）函数的表达不同。特别的，逻辑回归使用 S型函数（sigmoid）或 “logistic” 激活函数（activation function ）而不是线性回归的连续输出。当我们深入到实现中去，我们会了解到更多。
在开始之前，让我们导入并检查将要处理的数据集。

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

import os  
path = os.getcwd() + '\data\ex2data1.txt'  
data = pd.read_csv(path, header=None, names=['Exam 1', 'Exam 2', 'Admitted'])  
data.head()  
 
 1
2
3
4
5
6
7
8
9
 
 
  
  
 
 
 
 1
2
3
4
5
6
7
8
9

	Exam 1	Exam 2	Admitted
0	34.623660	78.024693	0
1	30.286711	43.894998	0
2	35.847409	72.902198	0
3	60.182599	86.308552	1
4	79.032736	75.344376	1

在数据集中有两个连续独立变量—— Exam 1 和 Exam 2。我们预测目标是“录取（Admitted）”标签，它二值的。1 表示被录取 0 表示没有被录取。让我们看看两次分数的散点图并用颜色编码来可视化样本的正类和负类（positive or negative）（录取情况）。

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')  
 
 1
2
3
4
5
6
7
8
9
 
 
  
  
 
 
 
 1
2
3
4
5
6
7
8
9

这里写图片描述

从图中我们可以看出有一个近似线性决策边界（ linear decision boundary）。它有点弧度所以我们不能用一条直线将所有样本正确分类，但我们能得到一个相当接近的。现在我们需要实现逻辑回归，我们可以训练一个模型来找出最优决策边界并作类预测。第一步实现 sigmoid 函数。

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

这个函数是逻辑回归输出的激活函数。它将连续输入映射到 0-1 之间。这一值可以被解释为类概率（class probability），或者是输入样本被分为正类的可能性。用这个可能性和阈值一起，我们可以获得一个离散标签预测。它有助于可视化函数输出，来看看它真正在做什么。

nums = np.arange(-10, 10, step=1)

fig, ax = plt.subplots(figsize=(12,8))  
ax.plot(nums, sigmoid(nums), 'r')
 
 1
2
3
4
 
 
  
  
 
 
 
 1
2
3
4

这里写图片描述

我们下一步是写一个代价函数。代价函数在训练数据中给定一组模型参数评估模型的性能。下面是逻辑回归的代价函数。

def cost(theta, X, y):  
    theta = np.matrix(theta)
    X = np.matrix(X)
    y = np.matrix(y)
    first = np.multiply(-y, np.log(sigmoid(X * theta.T)))
    second = np.multiply((1 - y), np.log(1 - sigmoid(X * theta.T)))
    return np.sum(first - second) / (len(X))
 
 1
2
3
4
5
6
7
 
 
  
  
 
 
 
 1
2
3
4
5
6
7

注意：我们将输出降到单个标量值，该标量是“误差”的总和，而“误差”被量化为一个函数，该函数是模型分配的类概率和样本的真标签（true label）之差。这个实现是完全矢量化的——它将模型对整个数据集的预测的计算集于一个表达式sigmoid(X * theta.T)。如果这里的数学不能理解，请参考我在上面给出的链接中的练习文本有详细地解释。

我们可以看看代价函数能否正常工作，但我们先需要设置一番。

# add a ones column - this makes the matrix multiplication work out easier
data.insert(0, 'Ones', 1)

# set X (training data) and y (target variable)
cols = data.shape[1]  
X = data.iloc[:,0:cols-1]  
y = data.iloc[:,cols-1:cols]

# convert to numpy arrays and initalize the parameter array theta
X = np.array(X.values)  
y = np.array(y.values)  
theta = np.zeros(3)  
 
 1
2
3
4
5
6
7
8
9
10
11
12
 
 
  
  
 
 
 
 1
2
3
4
5
6
7
8
9
10
11
12

我喜欢检查数据结构的 shape，我总是让自己能够感知数据。这项技术在实现矩阵乘法时非常有用。

X.shape, theta.shape, y.shape  
 
 1
 
 
  
  
 
 
 
 1

((100L, 3L), (3L,), (100L, 1L))
 
 1
 
 
  
  
 
 
 
 1

现在让我们来计算模型参数为零时初始解的 cost，在这里用 θ 表示模型参数。

cost(theta, X, y) 
 
 1
 
 
  
  
 
 
 
 1

0.69314718055994529
 
 1
 
 
  
  
 
 
 
 1

现在我们有了一个可以正确的代价函数，下一步就是写一个函数计算模型参数的梯度，搞清楚在训练数据中参数们对模型的输出是怎样改善的。回想一下梯度下降函数，我们并不是仅仅随意地给出几个参数值，然后看看谁表现地最好。在每一次训练迭代中我们更新参数保证让它们在同一方向移动以降低训练误差（例如 cost ）。我们之所以能这么做是因为代价函数可微。这里所涉及的微积分方程推导远远超出了本文讨论的范围，完整的方程请看练习文本。下面是该函数。

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 = 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
 
 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
 
 
  
  
 
 
 
 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15

注意在这里我们并不执行梯度下降——我们计算单独的梯度。在这个练习中，一个叫 fminunc 的 Octave 函数是用来优化给定函数的参数的以计算代价和梯度的。既然我们在用 Python ，我们可以使用 SciPy 的优化 API 来实现相同功能。

import scipy.optimize as opt  
result = opt.fmin_tnc(func=cost, x0=theta, fprime=gradient, args=(X, y))  
cost(result[0], X, y)  
 
 1
2
3
 
 
  
  
 
 
 
 1
2
3

0.20357134412164668
 
 1
 
 
  
  
 
 
 
 1

现在我们的数据集有了优化的模型参数。接下来我们需要写一个函数用学到的参数 θ 的预测。我们可以使用这个函数对我们的分类器精度进行评估。

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

theta_min = np.matrix(result[0])  
predictions = predict(theta_min, X)  
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)  
 
 1
2
3
4
5
6
7
8
9
 
 
  
  
 
 
 
 1
2
3
4
5
6
7
8
9

accuracy = 89%
 
 1
 
 
  
  
 
 
 
 1

此时我们的逻辑回归分类器有 89% 的几率正确预测一个学生能否被录取。还不错！记住这虽然是数据集的精度，但我们并没有保持一个 hold-out 集也没有使用交叉验证（cross-validation）来得到真正的精确度近似值所以这个数字可能高于它的真实性能（这个主题在下一个练习中讨论）。

逻辑回归正则化（Regularized Logistic Regression）

现在我们有一个可行的逻辑回归实现，我们将通过添加正则化来改善该算法。在代价函数中正则化是一个项，它使算法更喜欢“简单的”模型（在这种情况下，模型将缩减系数）。该理论认为，这有助于减少过拟合（overfitting），提高模型的泛化（generalize）能力。我们将应用逻辑回归的正则化版本到一个稍微更加有挑战性的问题上。假设你是一个工厂的产品经理，你有一些芯片两种不同检测方法的检测结果。从两种检测中，你将决定哪些芯片可用。为了帮助你作出决定，你有通过的芯片的检测结果数据集，从其中你可以建立逻辑回归模型。

让我们通过数据可视化开始。

path = os.getcwd() + '\data\ex2data2.txt'  
data2 = pd.read_csv(path, header=None, names=['Test 1', 'Test 2', 'Accepted'])

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')  
 
 1
2
3
4
5
6
7
8
9
10
11
12
 
 
  
  
 
 
 
 1
2
3
4
5
6
7
8
9
10
11
12

这里写图片描述

这数据看起来比之前的样本要复杂一些。特别的，你将注意到在这个数据并没表现非常好的线性决策边界。解决这个问题的一种方法是用一种像逻辑回归之类的线性技术从原始特征多项式中分离出构造特征。我们可以试着创造一些多项式特征反馈到分类器中。

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()  
 
 1
2
3
4
5
6
7
8
9
10
11
12
13
14
 
 
  
  
 
 
 
 1
2
3
4
5
6
7
8
9
10
11
12
13
14

	Accepted	Ones	F10	F20	F21	F30	F31	F32
0	1	1	0.051267	0.002628	0.035864	0.000135	0.001839	0.025089
1	1	1	-0.092742	0.008601	-0.063523	-0.000798	0.005891	-0.043509
2	1	1	-0.213710	0.045672	-0.147941	-0.009761	0.031616	-0.102412
3	1	1	-0.375000	0.140625	-0.188321	-0.052734	0.070620	-0.094573
4	1	1	-0.513250	0.263426	-0.238990	-0.135203	0.122661	-0.111283

现在我们需要修改代价和梯度函数使之包含正则化项。在每种情况下，正则化项加在先前的计算上。下面是更新后的代价函数。

def costReg(theta, X, y, learningRate):  
    theta = np.matrix(theta)
    X = np.matrix(X)
    y = np.matrix(y)
    first = np.multiply(-y, np.log(sigmoid(X * theta.T)))
    second = np.multiply((1 - y), np.log(1 - sigmoid(X * theta.T)))
    reg = (learningRate / 2 * len(X)) * np.sum(np.power(theta[:,1:theta.shape[1]], 2))
    return np.sum(first - second) / (len(X)) + reg
 
 1
2
3
4
5
6
7
8
 
 
  
  
 
 
 
 1
2
3
4
5
6
7
8

注意，我们添加了一个叫 reg 的新变量，它是参数值的函数。随着参数的增长，加在代价函数上的惩罚（ penalization）也在增加。也请注意，我们给函数加了一个新的“学习率”参数。这也是正则化项的一部分。学习率给我们一个新的超参数（hyper-parameter）让我们能够调整正则化在代价函数中的比重。

接下来我们将添加正则化到梯度函数。

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

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

    error = 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)) + ((learningRate / len(X)) * theta[:,i])

    return grad
 
 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
 
 
  
  
 
 
 
 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19

就和代价函数一样，正则化项加在原计算上。然而，不像代价函数，我们包含逻辑（ logic）以确保第一个参数不被正则化。这个决定后的直觉是第一个参数考虑到了模型的偏差或截距所以不应被惩罚。

可以像之前那样对新函数测试一下。

# set X and y (remember from above that we moved the label to column 0)
cols = data2.shape[1]  
X2 = data2.iloc[:,1:cols]  
y2 = data2.iloc[:,0:1]

# convert to numpy arrays and initalize the parameter array theta
X2 = np.array(X2.values)  
y2 = np.array(y2.values)  
theta2 = np.zeros(11)

learningRate = 1

costReg(theta2, X2, y2, learningRate)  
 
 1
2
3
4
5
6
7
8
9
10
11
12
13
 
 
  
  
 
 
 
 1
2
3
4
5
6
7
8
9
10
11
12
13

0.6931471805599454
 
 1
 
 
  
  
 
 
 
 1

也可以重用之前的优化代码来找出最优模型参数。

result2 = opt.fmin_tnc(func=costReg, x0=theta2, fprime=gradientReg, args=(X2, y2, learningRate))  
result2  
 
 1
2
 
 
  
  
 
 
 
 1
2

(array([ 0.35872309, -3.22200653, 18.97106363, -4.25297831, 18.23053189, 20.36386672, 8.94114455, -43.77439015, -17.93440473, -50.75071857, -2.84162964]), 110, 1)
 
 1
 
 
  
  
 
 
 
 1

最后，可以应用之前相同的方法来创建对训练数据的标签预测并评估模型性能。

theta_min = np.matrix(result2[0])  
predictions = predict(theta_min, X2)  
correct = [1 if ((a == 1 and b == 1) or (a == 0 and b == 0)) else 0 for (a, b) in zip(predictions, y2)]  
accuracy = (sum(map(int, correct)) % len(correct))  
print 'accuracy = {0}%'.format(accuracy)  
 
 1
2
3
4
5
 
 
  
  
 
 
 
 1
2
3
4
5