逻辑回归

题目1

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

可视化数据

导入需要用到的库

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

# 读取数据
path = "D:\gitcode\Coursera-ML-AndrewNg-Notes-master\Coursera-ML-AndrewNg-Notes-master\code\ex2-logistic regression\ex2data1.txt"
data = pd.read_csv(path, names=["exam1", "exam2", "admitted"])
# print(data.head())

让我们创建两个分数的散点图，并使用颜色编码来可视化，如果样本为1的（被接纳）或0的（未被接纳）。

# 可视化数据
positive = data[data['admitted'].isin([1])]
negative = data[data['admitted'].isin([0])]
fig, ax = plt.subplots(figsize=(6, 5))
ax.scatter(positive['exam1'], positive['exam2'], c='b', label='admitted')
ax.scatter(negative['exam1'], negative['exam2'], c='r', label='not admitted')
plt.xlabel("exam1 score")
plt.ylabel("exam2 score")
plt.legend()
plt.show()

sigmoid 函数

g 代表一个常用的逻辑函数（logistic function）为S形函数（Sigmoid function），公式为： \[g\left( z \right)=\frac{1}{1+{{e}^{-z}}}\]
得到逻辑回归模型的假设函数：
\[{{h}_{\theta }}\left( x \right)=\frac{1}{1+{{e}^{-{{\theta }^{T}}X}}}\]
定义sigmoid函数

# sigmoid函数
def sigmoid(z):
    return 1 / (1 + np.exp(-z))

代价(损失)函数

现在，我们需要编写代价函数来评估结果。
代价函数：
\(J\left( \theta \right)=\frac{1}{m}\sum\limits_{i=1}^{m}{[-{{y}^{(i)}}\log \left( {{h}_{\theta }}\left( {{x}^{(i)}} \right) \right)-\left( 1-{{y}^{(i)}} \right)\log \left( 1-{{h}_{\theta }}\left( {{x}^{(i)}} \right) \right)]}\)
在用代码实现时，需要注意，theta维1维向量，使用numpy进行矩阵相乘时，需要是numpy数组，否则可能会存在问题

# 定义损失函数
def computeCost(theta, X, y):
    '''
    :param theta: 1维向量 n个特征
    :param X: numpy数组 m行n列
    :param y: numpy数组 m行1列
    :return: 损失数值
    '''
    theta = np.reshape(theta, (1, theta.shape[0]), order='C')
    m = X.shape[0]
    first = (-y) * np.log(sigmoid(X @ theta.T))
    second = (1 - y) * np.log(1 - sigmoid(X @ theta.T))
    return np.sum(first - second) / m

梯度函数

这是梯度公式（batch gradient descent）：

\[\frac{\partial J\left( \theta \right)}{\partial {{\theta }_{j}}}=\frac{1}{m}\sum\limits_{i=1}^{m}{({{h}_{\theta }}\left( {{x}^{(i)}} \right)-{{y}^{(i)}})x_{_{j}}^{(i)}} \]

# 计算梯度
def gradient(theta, X, y):
    '''
    :param theta: n维向量 n个特征
    :param X: numpy数组 m行n列
    :param y: numpy数组 m行n列
    :return: n维梯度 n个特征(n,)
    '''
    theta = np.reshape(theta, (1, theta.shape[0]), order='C')
    m = X.shape[0]
    return (((sigmoid(X @ theta.T) - y).T @ X) / m).flatten()

在实现时，使用矩阵乘法进行得到梯度矩阵，最终返回得到一个n维向量，因为在后续使用tnc算法进行优化时，要求optimize的函数返回值为向量。此外，注意，我们实际上没有在这个函数中执行梯度下降，我们仅仅在计算一个梯度步长。在练习中，一个称为“fminunc”的Octave函数是用来优化函数来计算成本和梯度参数。由于我们使用Python，我们可以用SciPy的“optimize”命名空间来做同样的事情。

SciPy's truncated newton（TNC）实现寻找最优参数

我们可以利用scipy.optimize提供的fmin_tnc进行寻找最优参数

# 利用tcn进行最优化拟合
# 方法一直接利用opt.fmin_tnc进行优化
import scipy.optimize as opt
result = opt.fmin_tnc(func=computeCost, x0=theta, fprime=gradient, args=(X, y))
theta = np.reshape(result[0], (1, 3), order='C')

也可以利用opt.minimize进行优化，设置优化方法为TNC

# 方法二利用opt.minimize进行优化,设置方法为TNC
result = opt.minimize(fun=computeCost, x0=theta, jac=gradient, args=(X, y), method='TNC')
theta = np.reshape(result['x'], (1,3), order='C')
print(result)

但是，注意，两种方法得到的结果的类型不同，可以将结果打印出来后取出所需要的最优参数

JZTXT

逻辑回归练习