【深度之眼吴恩达机器学习第四期】笔记（七）

模型优化

准备数据：

import numpy as np
import scipy.io as sio
import scipy.optimize as opt
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns

def load_data():
    # 因为ex5data的d['X']的shape=(12, 1)，而且pandas不能用2维数据构造dataframe,所以在返回前ravel一下，把它变成一维的
    d = sio.loadmat('ex5data1.mat')
    return map(np.ravel, [d['X'], d['y'], d['Xval'], d['yval'], d['Xtest'], d['ytest']])

X, y, Xval, yval, Xtest, ytest = load_data()
# 展示数据
df = pd.DataFrame({'water_level':X, 'flow':y})
sns.lmplot('water_level', 'flow', data=df, fit_reg=False, height=7)
plt.show()

# 三个数据集都改回二维，并加上一列x0
X, Xval, Xtest = [np.insert(x.reshape(x.shape[0], 1), 0, np.ones(x.shape[0]), axis=1) for x in (X, Xval, Xtest)]

展示数据：

损失函数

def cost(theta, X, y):
# INPUT：参数值theta，数据X,标签y
# OUTPUT：当前参数值下代价函数

    # STEP1：获取样本个数
    m = X.shape[0]
    
    # STEP2：计算代价函数
    error = X@theta.T-y
    cost = np.sum(np.power(error,2))/(2*m)
    return cost

theta = np.ones(X.shape[1])
# 303.9515255535976
cost(theta, X, y)

梯度函数

def gradient(theta, X, y):
# INPUT：参数值theta，数据X,标签y
# OUTPUT：当前参数值下梯度

    # STEP1：获取样本个数
    m = X.shape[0]
    
    # STEP2：计算梯度
    grad = X.T@(X@theta.T-y)/m
    return grad

# array([-15.30301567, 598.16741084])
gradient(theta, X, y)

正则化梯度与代价函数

def regularized_gradient(theta, X, y, l=1):
# INPUT：参数值theta，数据X,标签y
# OUTPUT：当前参数值下梯度

    # STEP1：获取样本个数
    m = X.shape[0]
    
    # STEP2：计算正则化梯度
    regularized_term = l/m*theta
    regularized_term[0] = 0  # 不用对theta0正则化
    return gradient(theta, X, y) + regularized_term

# array([-15.30301567, 598.25074417])
regularized_gradient(theta, X, y)

def regularized_cost(theta, X, y, l=1):
    m = X.shape[0]
    # theta[1:]表示不用对theta0正则化
    regularized_term = (l / (2 * m)) * np.power(theta[1:], 2).sum()
    return cost(theta, X, y) + regularized_term

学习曲线

def linear_regression_np(X, y, l=1):
# INPUT：数据X,标签y，正则化参数l
# OUTPUT：当前参数值下使用梯度下降得到的结果

    # STEP1：初始化参数
    theta = np.ones(X.shape[1])
    
    # STEP2：调用优化算法拟合参数
    res = opt.minimize(fun=regularized_cost,
                       x0=theta,
                       args=(X,y,l),
                       method='TNC',
                       jac=regularized_gradient,
                       options={'disp': True})
    return res

theta = np.ones(X.shape[0])
final_theta = linear_regression_np(X, y, l=0).get('x')

b = final_theta[0]
m = final_theta[1]

# 原始数据的点
plt.scatter(X[:,1], y, label="Training data")
# 预测的直线
plt.plot(X[:, 1], X[:, 1]*m + b, label="Prediction")
plt.legend(loc=2)
plt.show()

使用不同数量的训练数据集：

training_cost, cv_cost = [], []

# STEP1：获取样本个数，遍历每个样本
m = X.shape[0]
for i in range(1, m+1):
    # STEP2：计算当前样本的代价
    # 在lambda=0的时候，取前i个数据进行训练
    res = linear_regression_np(X[:i, :], y[:i], l=0)
    # 计算梯度下降得到参数theta的损失，不要加上正则项
    tc = cost(res.get('x'),X[:i, :], y[:i])
    cv = cost(res.get('x'),Xval, yval)
    
    # STEP3：把计算结果存储至预先定义的数组training_cost, cv_cost中
    training_cost.append(tc)
    cv_cost.append(cv)

# 画出训练损失和验证损失随数据量变化而变化的图像
plt.plot(np.arange(1, m+1), training_cost, label='training cost')
plt.plot(np.arange(1, m+1), cv_cost, label='cv cost')
plt.legend(loc=1)
plt.show()

这个图在数据量小的时候两个损失相差很大，然后损失间隔迅速减少，最后增加数据量间隔也没有太大的变化，这是欠拟合的情况。这种情况还有一个特征，就是最后的损失相对较高，这里的y在[-5,30]，每个数据平均有√(25*2)≈7的误差失其实也挺高的了。

准备高次项数据：

def prepare_poly_data(*args, power):
    # args: 会保持输入参数X, Xval等的顺序
    def prepare(x):
        # 特征映射，变为高次项的组合
        df = poly_features(x, power=power)
        # 归一化处理
        ndarr = normalize_feature(df).values
        # 添加偏置项x0
        return np.insert(ndarr, 0, np.ones(ndarr.shape[0]), axis=1)
    return [prepare(x) for x in args]

def poly_features(x, power, as_ndarray=False):  #特征映射
    data = {'f{}'.format(i): np.power(x, i) for i in range(1, power + 1)}
    df = pd.DataFrame(data)
    return df.as_matrix() if as_ndarray else df
    
def normalize_feature(df):
    # 归一化
    return df.apply(lambda column: (column - column.mean()) / column.std())

X, y, Xval, yval, Xtest, ytest = load_data()
poly_features(X, power=3)

原来的x→(x,x²,x³)

x→x,……,x^8，然后归一化为[-1,1]，最后输出前三行

X_poly, Xval_poly, Xtest_poly= prepare_poly_data(X, Xval, Xtest, power=8)
X_poly[:3, :]

学习曲线函数：

def plot_learning_curve(X, y, Xval, yval, l=0):
# INPUT：训练数据集X,y，交叉验证集Xval，yval，正则化参数l
# OUTPUT：当前参数值下的学习曲线

    # STEP1：初始化参数，获取样本个数，开始遍历
    training_cost, cv_cost = [], []
    m = X.shape[0]
    for i in range(1, m + 1):
        # STEP2：调用之前写好的拟合数据函数进行数据拟合
        res = linear_regression_np(X[:i,:],y[:i],l)
        
        # STEP3：计算样本代价
        tc = cost(res.get('x'),X[:i,:],y[:i])
        cv = cost(res.get('x'),Xval,yval)
        
        # STEP3：把计算结果存储至预先定义的数组training_cost, cv_cost中
        training_cost.append(tc)
        cv_cost.append(cv)
    plt.plot(np.arange(1, m + 1), training_cost, label='training cost')
    plt.plot(np.arange(1, m + 1), cv_cost, label='cv cost')
    plt.legend(loc=1)

# lambda=0
plot_learning_curve(X_poly, y, Xval_poly, yval, l=0)
plt.show()
# lambda=1
plot_learning_curve(X_poly, y, Xval_poly, yval, l=1)
plt.show()
# lambda=100
plot_learning_curve(X_poly, y, Xval_poly, yval, l=100)
plt.show()

lambda=0的时候，训练损失一直很小，而且随着训练数据量的增多，训练损失和验证损失间的间隔依旧很大，这是过拟合的表现。

lambda=1的学习曲线，验证集损失相对较低了。

lambda=100的时候，间隔迅速减小，最后损失依旧很高，所以这是欠拟合了。

选择最优的超参数lambda

# 候选的超参数lambda
l_candidate = [0, 0.001, 0.003, 0.01, 0.03, 0.1, 0.3, 1, 3, 10]
training_cost, cv_cost = [], []
for l in l_candidate:
    res = linear_regression_np(X_poly, y, l)    
    tc = cost(res.x, X_poly, y)
    cv = cost(res.x, Xval_poly, yval)    
    training_cost.append(tc)
    cv_cost.append(cv)

plt.plot(l_candidate, training_cost, label='training')
plt.plot(l_candidate, cv_cost, label='cross validation')
plt.legend(loc=2)
plt.xlabel('lambda')
plt.ylabel('cost')
plt.show()

损失随lambda变化而变化的图像，我们应该选择验证集损失最小的哪一个超参数，这里是lambda=1。

选择验证集损失最小的超参数，这里是lambda=1，计算准确率。

best_l=l_candidate[np.argmin(cv_cost)]
theta = linear_regression_np(X_poly, y, best_l).x
# test cost(l=1) = 7.466265914249742
print('test cost(l={}) = {}'.format(best_l, cost(theta, Xtest_poly, ytest)))