Python机器学习-吴恩达课后作业-多元分类-数字识别ex3

1. 数字识别分析

分析对象

自动手写数字识别在今天被广泛使用——从识别信封上的邮政编码(邮政编码)到识别银行支票上写的金额。在本练习中，您将使用逻辑回归和神经网络来识别手写数字(从0到9)。您将扩展之前的逻辑回归实现，并将其应用于one-vs-all分类。

• 数据以.mat格式储存，mat格式是matlab的数据存储格式，按照矩阵保存，与numpy数据格式兼容，适合于各种数学运算，因此主要使用numpy进行运算。
• 数据中有5000个训练样例，其中每个训练样例是一个20×20像素灰度图像的数字，每个像素由一个浮点数表示，该浮点数表示该位置的灰度强度。每个20×20像素的网格被展开成一个400维的向量。这些每个训练样例都变成数据矩阵X中的一行。这就得到了一个5000×400矩阵X，其中每一行都是手写数字图像的训练样例。训练集的第二部分是一个包含训练集标签的5000维向量y，“0”的数字标记为“10”，而“1”到“9”的数字按自然顺序标记为“1”到“9”。

计算程序

import numpy as np
from scipy.io import loadmat
from scipy.optimize import minimize


def load_data(path):
    # 加载数据
    data = loadmat(path)
    return data


def deal_data(data):
    # 处理数据
    # data["X"].shape = (5000, 400)
    # data["y"].shape = (5000, 1)
    X = np.matrix(data["X"])
    # 检查是否已经包含一列全为 1，如果没有，则添加
    if X[:, 0].sum() != X.shape[0]:
        X = np.c_[np.ones((X.shape[0], 1)), X]
    y = np.matrix(data['y'])
    return X, y


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


def gradient(theta, X, y, learningRate):
    # 计算梯度
    theta = np.matrix(theta)
    error = sigmoid(X * theta.T) - y
    grad = ((X.T * error) / len(X)).T + ((learningRate / len(X)) * theta)
    # 截距梯度没有正则化
    grad[0, 0] = np.sum(np.multiply(error, X[:, 0])) / len(X)
    # 将 grad 展开为一维数组，其会将每一行进行拼接，最终的格式为（：，1）
    return np.array(grad).ravel()


def cost(theta, X, y, lamada):
    # 计算代价函数
    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)))
    # Ensure that 'first' and 'second' are 1D arrays
    first = np.squeeze(np.asarray(first))
    second = np.squeeze(np.asarray(second))
    reg = (lamada / (2 * len(X))) * np.sum(np.square(theta[:, 1:]))
    return np.sum(first - second) / len(X) + reg


# 使用训练完毕的分类器预测每个图像的标签
# 构建分类器
def ones_vs_all(X, y, num_labels, lamada):
    rows = X.shape[0]
    params = X.shape[1]
    # k X (n + 1)数组中每个分类器的参数
    all_theta = np.zeros((num_labels, params))
    # 标签是1-indexed而不是0-indexed，遍历每一种分类器
    for i in range(1, num_labels + 1):
        # 初始化当前分类器的 theta 数组
        theta = np.zeros(params)
        # 多元分类中，采用的是一对多的思想，除了当前分类，其他 y 值全部设置为 0
        y_index = np.array([1 if label == i else 0 for label in y])
        y_index = np.reshape(y_index, (rows, 1))
        # 最小化目标函数
        fmin = minimize(fun=cost, x0=theta, args=(X, y_index, lamada), method='TNC', jac=gradient)
        # 将当前分类器的 theta 数组整合到总的 theta 数组中
        # 通过 fmin.x 可以获得对应的 theta
        all_theta[i - 1, :] = fmin.x
    return all_theta


def predict_all(X, all_theta):
    # 计算每个类在每个训练实例上的类概率
    h = sigmoid(X * all_theta.T)
    # 以最大概率创建索引数组
    h_argmax = np.argmax(h, axis=1)
    # 因为我们的数组是0-indexed的，所以我们需要为真标签预测添加一个 indexed
    # 通过对每一行 X 参数进行预测，对应结果的 h 中每一行为分别 10 个类的预测结果，找到预测的最高值所在列数，即认为该行参数对应的是类别为 列数+1
    h_argmax = h_argmax + 1
    return h_argmax


if __name__ == '__main__':
    path = '.\\data\\ex3data1.mat'
    data = load_data(path)
    X, y = deal_data(data)
    print(X[:, 0])
    print(y)
    all_theta = ones_vs_all(X, y, 10, 1)
    y_pred = predict_all(X, all_theta)
    correct = [1 if a == b else 0 for (a, b) in zip(y_pred, data['y'])]
    accuracy = (sum(map(int, correct)) / float(len(correct)))
    print('accuracy = {0}%'.format(accuracy * 100))

运行结果

accuracy = 94.46%

说明

• return np.array(grad).ravel() 这一行的作用是将梯度矩阵展平为一维数组。具体解释如下：

• grad：在梯度下降步骤中计算的梯度矩阵。这是一个二维数组。
• np.array(grad)：将梯度矩阵转换为NumPy数组。
• .ravel()：将二维数组展平为一维数组（一个扁平的向量）。

在逻辑回归实现的上下文中，这个展平后的向量作为梯度被返回，用于优化算法。优化算法期望一个一维数组作为梯度，而 ravel() 是一个方便的方法来实现这种展平。

• ravel()用法举例：

ravel() 是 NumPy 中用于将多维数组展平为一维数组的函数。以下是一个简单的例子：

import numpy as np

# 创建一个二维数组
arr_2d = np.array([[1, 2, 3], [4, 5, 6]])

# 使用 ravel() 将二维数组展平为一维数组
arr_1d = arr_2d.ravel()

print("Original 2D array:")
print(arr_2d)
print("\nFlattened 1D array:")
print(arr_1d)

输出结果：

Original 2D array:
[[1 2 3]
 [4 5 6]]

Flattened 1D array:
[1 2 3 4 5 6]

 在这个例子中，arr_2d 是一个二维数组，使用 ravel() 后得到了一个一维数组 arr_1d。展平后的数组按照原始数组的行优先顺序排列。

• y_i = np.array([1 if label == i else 0 for label in y]) 举例说明：

y_i = np.array([1 if label == i else 0 for label in y]) 这行代码是用于创建一个与原始标签 y 大小相同的新数组 y_i。这个新数组是一个二进制数组，其中元素为 1 表示原始标签等于 i，元素为 0 表示原始标签不等于 i。

假设 y 是一个包含以下标签的数组：[2, 0, 1, 1, 2, 0]，而我们想要创建一个新数组 y_i，其中 i 的值为 1。那么这个新数组的规则是，如果原始标签等于 i（这里是 1），则新数组对应位置的值为 1，否则为 0。

import numpy as np

# 假设原始标签数组 y
y = np.array([2, 0, 1, 1, 2, 0])

# 设置 i 的值
i = 1

# 使用列表推导式创建新数组 y_i
y_i = np.array([1 if label == i else 0 for label in y])

# 输出结果
print(y_i)

输出结果：

[0 0 1 1 0 0]

• 对于第一个元素 2，因为它不等于 i（1），所以对应的 y_i 元素为 0。
• 对于第二个元素 0，因为它等于 i（1），所以对应的 y_i 元素为 1。
• 依此类推，得到新数组 y_i。

• fmin = minimize(fun=cost, x0=theta, args=(X, y, learning_rate), method='TNC', jac=gradient)

这行代码使用 SciPy 库中的 minimize 函数进行优化。

• fun: 代表要最小化的目标函数，这里是 cost 函数。
• x0: 是优化的初始参数值，这里是 theta。
• args: 是传递给目标函数的其他参数，这里是 (X, y, learning_rate)。
• method: 是指定优化算法的方法，这里是 'TNC'，表示 Truncated Newton Conjugate-Gradient。
• jac: 是目标函数的梯度，这里是 gradient 函数。

假设我们有一个简单的目标函数和一个梯度函数：

import numpy as np
from scipy.optimize import minimize

# 目标函数
def simple_cost(theta):
    return (theta[0] - 2) ** 2 + (theta[1] - 3) ** 2

# 梯度函数
def simple_gradient(theta):
    return np.array([2 * (theta[0] - 2), 2 * (theta[1] - 3)])

# 初始参数
theta_init = np.array([0, 0])

# 使用 minimize 进行优化
result = minimize(fun=simple_cost, x0=theta_init, method='TNC', jac=simple_gradient)

# 输出优化结果
print(result)

这里的 simple_cost 是一个简单的二次目标函数，simple_gradient 是其梯度。我们使用 minimize 来找到最小值，初始参数是 [0, 0]。运行结果会包含优化后的参数值以及其他有关优化过程的信息。

输出结果为：

message: Local minimum reached (|pg| ~= 0)
 success: True
  status: 0
     fun: 0.0
       x: [ 2.000e+00  3.000e+00]
     nit: 1
     jac: [ 0.000e+00  0.000e+00]
    nfev: 3

其实就是找到和（2，3）距离最近的点，当然是点本身。