练习1和练习3

---

前言

CG平台上完成

---

一、线性回归（练习1）

1、实现简单示例函数

import numpy as np
a = np.eye(5)
a

输出如下：
array(
[[1., 0., 0., 0., 0.],
[0., 1., 0., 0., 0.],
[0., 0., 1., 0., 0.],
[0., 0., 0., 1., 0.],
[0., 0., 0., 0., 1.]]
)

2、单变量线性回归

在该部分练习中，将实现单变量线性回归并用来预测餐车的利润。
假设你是一家餐厅的领导，正在考虑在不同的城市开设新的分店。该连锁店已经在不同的城市有了餐车，并且你能够获得每个城市的人口和利润数据。
现在需要使用这些数据来帮助你选择下一个被扩展的城市。
文件ex1data1.txt包含线性回归问题的数据集。第一列数据对应城市人口，第二列数据对应那座城市的餐车的利润。利润为负时表示亏损。

2.1 绘制数据

# 引入所需要的库文件
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import os

%matplotlib inline

# 数据存储路径
path = 'ex1data1.txt'

# 读入相应的数据文件
data = pd.read_csv(path, header=None,names=['Population','Profit'])

#查看数据的前五条
data.head(5)

输出如下：

接下来需要实现数据可视化的代码

data.plot(kind='scatter', x='Population', y='Profit',c='red', figsize=(12,8))

输出如下：

2.2 梯度下降

在上一部分的练习中，我们已经将所需要用到的数据加载至变量data中，并为其列分别进行命名。
接下来，我们在数据中添加了一个维度来拟合截距项 𝜃0。并将初始参数值设为0，学习率 𝛼设为0.01。
代码实现

#在列索引为0处添加数据列，该列值均为1
data.insert(0, 'Ones', 1)

#获取数据列数
cols = data.shape[1]

#对变量X和y进行初始化,并将其数据类型转换为矩阵
X = data.iloc[:,0:cols-1]
y = data.iloc[:,cols-1:cols]
X = np.matrix(X.values)
y = np.matrix(y.values)

#学习率、迭代次数的初始化
alpha = 0.01
iterations = 1500

计算成本

def computeCost(X, y, theta):
    inner = np.power(((X * theta.T) - y), 2)
    return np.sum(inner) / (2 * len(X))

theta = np.matrix(np.array([0,0]))
computeCost(X, y, theta)

梯度下降

def gradientDescent(X, y, theta, alpha, iters):
    temp = np.matrix(np.zeros(theta.shape))
    parameters = int(theta.ravel().shape[1])
    cost = np.zeros(iters)
    
    for i in range(iters):
        error = (X * theta.T) - y
        
        for j in range(parameters):
            term = np.multiply(error, X[:,j])
            temp[0,j] = theta[0,j] - ((alpha / len(X)) * np.sum(term))
            
        theta = temp
        cost[i] = computeCost(X, y, theta)
        
    return theta, cost

g, cost = gradientDescent(X, y, theta, alpha, iterations)
g
# 计算最终的参数所得到的成本值
computeCost(X, y, g)

#对拟合曲线进行绘制
x = np.linspace(data.Population.min(), data.Population.max(), 100)
f = g[0, 0] + (g[0, 1] * x)

fig, ax = plt.subplots(figsize=(12,8))
ax.plot(x, f, 'b', label='Prediction')
ax.scatter(data.Population, data.Profit, c='red',label='Traning Data')
ax.legend(loc=2)
ax.set_xlabel('Population')
ax.set_ylabel('Profit')
ax.set_title('Predicted Profit vs. Population Size')

输出如下：

2.3 可视化成本函数

为了更好地理解成本函数的迭代计算，将每一步计算的cost值进行记录并绘制。

fig, ax = plt.subplots(figsize=(12,8))
ax.plot(np.arange(iterations), cost, 'r')
ax.set_xlabel('Iterations')
ax.set_ylabel('Cost')
ax.set_title('Error vs. Training Epoch')

输出如下：

二、多分类问题（练习3）

1、数据集

2、数据可视化

图像在矩阵X中表示为400维向量（其中有5,000个）。 400维“特征”是原始20 x 20图像中每个像素的灰度强度。类标签在向量y中作为表示图像中数字的数字类。
接下来，我们需要加载数据。

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.io import loadmat

data = loadmat('ex3data1.mat')
data

并使用shape内置函数检查数据矩阵 𝑋 , 𝑦的形状：

data[‘X’].shape, data[‘y’].shape

输出：((5000, 400), (5000, 1))

3、逻辑回归的向量化

在该部分练习中，你需要将逻辑回归的实现修改为完全向量化（即没有可替换的 𝑓𝑜𝑟循环）。这是因为向量化代码除了简洁外，还能够利用线性代数优化，并且通常比迭代代码快得多。但是，如果从练习2中看到我们的代价函数已经完全向量化实现了，所以我们可以在这里重复使用相同的实现。

3.1 代价函数的向量化

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

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)

3.2 梯度的向量化

def gradient(theta, X, y):
    theta = np.matrix(theta)
    X = np.matrix(X)
    y = np.matrix(y)
    
    parameters = int(theta.ravel().shape[1])
    error = sigmoid(X * theta.T) - y
    
    grad = ((X.T * error) / len(X)).T
    
    return grad

3.3 正则化逻辑回归的向量化

接下来，需要编写代码实现正则化逻辑回归算法的代价函数和梯度的向量化实现。

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

def gradientReg(theta, X, y, learningRate):
    theta = np.matrix(theta)
    X = np.matrix(X)
    y = np.matrix(y)
    
    parameters = int(theta.ravel().shape[1])
    error = sigmoid(X * theta.T) - y
    
    grad = ((X.T * error) / len(X)).T + ((learningRate / len(X)) * theta)
    
    # intercept gradient is not regularized
    grad[0, 0] = np.sum(np.multiply(error, X[:,0])) / len(X)
    
    return np.array(grad).ravel()

4、多分类-分类器

现在我们已经定义了代价函数和梯度函数，现在我们需要构建一个分类器。对于手写字识别，我们有10个可能的类（0-9），但逻辑回归是一个二分类问题。在本练习中，你的任务是实现一对一全分类方法，其中具有 𝑘 个不同类的标签就有 𝑘个分类器，每个分类器在“类别 𝑖”和“不是 𝑖”之间决定。我们将把分类器训练包含在一个函数中，该函数计算10个分类器中的每个分类器的最终权重，并将权重返回为 [𝑘,𝑛+1] 的数组，其中 𝑛是参数数量。

from scipy.optimize import minimize

def one_vs_all(X, y, num_labels, learning_rate):
    rows = X.shape[0]
    params = X.shape[1]
    
    # k个分类器的参数，形状为(k,n+1)
    all_theta = np.zeros((num_labels, params + 1))
    
    # 插入值为1的列，用于计算截距项
    X = np.insert(X, 0, values=np.ones(rows), axis=1)
    
    # 将分类标签转换为0-1标识
    for i in range(1, num_labels + 1):
        theta = np.zeros(params + 1)
        y_i = np.array([1 if label == i else 0 for label in y])
        y_i = np.reshape(y_i, (rows, 1))
        
        # 使用minimize函数最小化代价函数
        fmin = minimize(fun=costReg, x0=theta, args=(X, y_i, learning_rate), method='TNC', jac=gradientReg)
        all_theta[i-1,:] = fmin.x
    
    return all_theta

检查下需要初始化的变量，以及变量的形状：

rows = data['X'].shape[0]
params = data['X'].shape[1]

all_theta = np.zeros((10, params + 1))

X = np.insert(data['X'], 0, values=np.ones(rows), axis=1)

theta = np.zeros(params + 1)

y_0 = np.array([1 if label == 0 else 0 for label in data['y']])
y_0 = np.reshape(y_0, (rows, 1))

X.shape, y_0.shape, theta.shape, all_theta.shape

输出：((5000, 401), (5000, 1), (401,), (10, 401))
其中， 𝑡ℎ𝑒𝑡𝑎是一维数组，因此当它被转换为计算梯度代码中的矩阵时，它变成了形状为 (1,401) 的矩阵。同时，我们需要检查 𝑦 中的类标签。

np.unique(data[‘y’])#看下有几类标签

输出：array([ 1, 2, 3, 4, 5, 6, 7, 8, 9, 10], dtype=uint8)
接下来，为了确保训练函数正确运行，运行以下代码查看是否得到合理的输出。

all_theta = one_vs_all(data[‘X’], data[‘y’], 10, 1)
all_theta

输出如下：

5、使用分类器进行预测

def predict_all(X, all_theta):
    rows = X.shape[0]
    params = X.shape[1]
    num_labels = all_theta.shape[0]
    
    # 与之前一样，需要插入一列确保矩阵形状
    X = np.insert(X, 0, values=np.ones(rows), axis=1)
    
    # 将其转换为矩阵
    X = np.matrix(X)
    all_theta = np.matrix(all_theta)
    
    # 计算每个训练样本所属每个类别的概率
    h = sigmoid(X * all_theta.T)
    
    # 创建具有最大概率的索引数组
    h_argmax = np.argmax(h, axis=1)
    
    # 因为我们的数组是零索引的，所以我们需要为真正的标签预测+1
    h_argmax = h_argmax + 1
    
    return h_argmax
    
	y_pred = predict_all(data['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))

---

总结

本次练习让我学习了线性回归和多分类问题的原理，对机器学习这一课程有了很深的理解。