No.004 训练模型

本章主要讲解不同的训练模型的方法。

线性回归

是假设函数，使用模型参数。存在n个特征，m个样本。对于线性模型的训练，设置模型参数，直到最拟合训练集。目标函数是最小化MSE。因此如下成本函数

最小化MSE，得到了

存在以下方法求解线性回归参数

import numpy as np
X = 2 * np.random.rand(100, 1)
# 这里假设y = 4 + 2x
y = 4 + 2 * X + np.random.rand(100, 1)
# X是全1和X组成。
X_b = np.c_[np.ones((100, 1)), X]
# 方法①：使用公式
theta_best = np.linalg.inv(X_b.T.dot(X_b)).dot(X_b.T).dot(y)
# 方法②：调用LinearRegression
from sklearn.linear_model import LinearRegression
lin_reg = LinearRegression()
lin_reg.fit(X, y)
print(lin_reg.intercept_, lin_reg.coef_)
# 方法③：调用最小二乘
theta_best_svd, residuals, rank, s = np.linalg.lstsq(X_b, y, rcond=1e-6)
print(theta_best_svd)
# 方法④：伪逆
print(np.linalg.pinv(X_b).dot(y))

• 伪逆是使用了奇异值分解（SVD）的标准矩阵分解方法，将训练集矩阵X分解为。伪逆的计算公式为。

  计算复杂度

• 是一个(n+1)*(n+1)的矩阵，因此矩阵求逆的复杂度达到了~，如果特征数量翻倍，那么复杂度将会乘以倍。

• 但是通过SVD的方法，复杂度将为。当特征加倍，复杂度提升了4倍。

梯度下降

梯度下降的中心思想：迭代调整参数，是成本函数最小化。首先使用随机的，然后逐步改善，每一步（步长很重要）尝试降低一点成本函数，直到算法收敛得到一个最小值。步长取决于超参数的学习率；如果学习率太低了，需要大量迭代才能够收敛，则会很耗时；但如果学习率太高了，则很容易发散。

批量梯度下降

对参数计算偏导数

对全部的θ计算偏导数，得到了梯度向量包含了全部的成本函数的偏导数

在计算梯度下降的每一步，都是基于完整的训练集X。因此速度会比较慢。但是梯度下降算法将会随着特征数量的上升而提升。
梯度向量乘以学习率，得到了下坡的步长

代码如下

eta = 0.1
n_iterations = 1000
m = 100
theta = np.random.rand(2, 1)
for i in range(n_iterations):
    gradients = 2 / m * X_b.T.dot(X_b.dot(theta) - y)
    theta -= eta * gradients

结果依然是array([[4.45329892], [1.96504494]])。

随机梯度下降

批量梯度下降的主要问题是：需要使用整个训练集计算每一步梯度，因此训练集很大，算法会特别慢。与之相反的是随机梯度下降，每一步在训练集中随机选择一个实例，仅仅基于该单个实例计算梯度。这种方法更快，并且每次迭代只需要操作少量的数据。
由于算法的随机性质，比批量梯度下降更加不规则。成本函数是不断上上下下，但从整体来看是缓缓下降。但即便达到了最小是，依然会反弹，最终会非常接近最小值。
当成本函数非常不规则，随机梯度下降有助于帮助算法挑出局部最小值，从而对于寻找全局最小值更加具有优势。
因此随机性的好处在于可以套利局部最优解，缺点是永远定位不到最小值——解决之道，逐步降低学习率。开始的步长比较大（有助于快速开展，并且逃离局部最小值），然后越来越小，让算法尽量靠近全局最小值，该算法称为模拟退火（simulated annealing）。

n_epochs = 50
t0, t1 = 5, 50
m = 100
def learning_schedule(t):
    return t0 / (t + t1)
theta = np.random.randn(2, 1)
for epoch in range(n_epochs):
    for i in range(m):
        random_index = np.random.randint(m)
        xi = X_b[random_index : random_index + 1]
        yi = y[random_index : random_index + 1]
        gradients = 2 * xi.T.dot(xi.dot(theta) - yi)
        eta = learning_schedule(epoch * m + i)
        theta = theta - eta * gradients

使用随机梯度下降，训练实例需要确保IID分布，从而确保平均而言将参数蜡像全局最优解。如果使用Scikit-Learn中的随机梯度下降，进行线性回归，可以使用SGDRegresor。该类默认优化平方误差成本函数。以下代码最多可以运行1000次轮次，或者直到下一轮次期间损失下降低于0.001为止（max_iter=1000，tol=1e-3）。使用默认的学习调度以0.1开始。不使用任何正则化。

from sklearn.linear_model import SGDRegressor
sgd_reg = SGDRegressor(max_iter=1000, tol=1e-3, penalty=None, eta0=0.1)
sgd_reg.fit(X, y.ravel())
print(sgd_reg.intercept_, sgd_reg.coef_)

小批量梯度下降

在每一步中，不是根据完整的训练集，也不是仅仅基于某一个实例计算梯度，小批量梯度下降称为小型批量的随机实例集上计算梯度。

多项式回归

例如，一个简单的二次方程生成非线性数据没有办法通过直线拟合。因此需要通过Poly-nominal-Features转换训练数据，让训练集中的每个特征的平方（二次多项式）添加为新特征

m = 100
X = 6 * np.random.rand(m, 1) - 3
X = np.array(sorted(list(X)))
y = 0.5 * X ** 2 + X + 2 * np.random.randn(m, 1)
# 添加新特征X**2
from sklearn.preprocessing import PolynomialFeatures
poly_features = PolynomialFeatures(degree=2, include_bias=False)
X_poly = poly_features.fit_transform(X)
print(X[0], X_poly[0])
# 分别为[-1.0904787] 和 [-1.0904787   1.18914379]
# X_poly现在包含了X原始特征+平方项特征。从而能够将LinearRegression模型拟合到扩展之后的训练输数据中。
# 拟合与预测
from sklearn.linear_model import LinearRegression
lin_reg = LinearRegression()
lin_reg.fit(X_poly, y)
print(lin_reg.intercept_, lin_reg.coef_)
# 结果为截距项=0.10470686，系数=0.94094301 0.32822587
y_pred = lin_reg.predict(X_poly)
# 画图
plt.plot(X, y_pred, 'r-', label='Pred')
plt.plot(X, y, 'b.')
plt.xlabel('X')
plt.ylabel('y')
plt.legend()
plt.show()

Poly-Nomial-Features能够将特征的所有组合添加到给定的多项式阶数中。例如degree=3、存在两个特征a和b的时候，存在。因此要小心参数爆炸。

学习曲线

如何判断欠拟合还是过拟合？如何确定模型复杂性。之前通过了交叉检验估计模型的泛化能力。还有一种方法是通过观察学习曲线：绘制的是在训练集和验证机上关于训练集大小的性能函数。

from sklearn.metrics import mean_absolute_error
from sklearn.model_selection import train_test_split
def plot_training_curves(model, X, y):
    # 拆分训练集和验证集
    X_train, X_val, y_train, y_val = train_test_split(X, y, test_size=0.2)
    # 训练误差和验证误差
    train_errors, val_errors = [], []
    # 随着训练样本量的逐渐增大
    for m in range(1, len(X_train)):
        model.fit(X_train[:m], y_train[:m])
        y_train_pred = model.predict(X_train[:m])
        y_val_pred = model.predict(X_val)
        train_errors.append(mean_absolute_error(y_train[:m], y_train_pred))
        val_errors.append(mean_absolute_error(y_val, y_val_pred))
    plt.plot(np.sqrt(train_errors), 'r-+', linewidth=2, label='train')
    plt.plot(np.sqrt(val_errors), 'b-+', linewidth=3, label='val')
    plt.legend()
    plt.xlabel('Training Set Size')
    plt.ylabel('RMSE')
    plt.show()
# X 和 y 来自上一列的二次方方程
plot_training_curves(LinearRegression(), X, y)

• 当训练集只有很少的几个实例的时候，模型能够很好地拟合；但是随着实例的增加，模型不能够完美地拟合训练数据，因为模型并不是线性、而是二次方，因此训练误差会一直上升，直到平稳。
• 而验证集上的误差一开始会很高，因为训练数据太少了，无法泛化；之后随着训练集的增大，验证误差逐渐下降。但是同样因为线性模型并不能够拟合，因此之后的误差保持平稳。
• 两条曲线均表现出欠拟合的状态。

如果使用10阶多项式进行拟合

from sklearn.pipeline import Pipeline
polynomial_reg = Pipeline([
    ('poly_features', PolynomialFeatures(degree=10, include_bias=False)),
    ('lin_reg', LinearRegression())
])
plot_training_curves(polynomial_reg, X, y)

偏差/方差权衡

模型的泛化误差能够表示为三个非常不同的误差之和：

• 偏差：泛化误差的原因在于错误的假设，例如假设数据是线性，但是实际上是二次；高偏差模型最有可能欠拟合训练集。
• 方差：由于模型对训练数据的细微变化过于敏感；具有很多自由度的模型（例如高阶多项式魔性）可能具有较高的方差，因此可能过拟合训练数据。
• 不可避免的误差：数据本身的噪声导致，减少这部分误差的唯一方法是清理数据。

正则化（Regularized）线性模型

减少过拟合的方法是对模型进行正则化（约束模型）：拥有的自由度越少，那么过拟合数据的难度就越大，正则化模型的简单方法是减少多项式的次数。正则化多项式模型的一种简单方法是减少多项式的次数。
对于线性模型，正则化通常是约束模型的权重实现。

岭回归Ridge

通过将等于的正则项调价到成本函数中，不仅拟合数据，而且使模型的权重尽可能小。超参数控制了模型的正则化程度。如果，那么岭回归=线性回归；如果很大，那么所有权重最终都会接近0，结果是经过数据均值的平均线。成本函数如下

注意：下标是从1开始，因此并没有对正则化。如果将定义为特征权重的向量，那么正则项等于表示为权重向量的范数。对于梯度下降，只需将田间道MSE梯度向量即可。
在执行岭回归之前，需要对数据进行缩放（例如，StandardScaler），因为对于输入特征的缩放很敏感。大部分正则化都需要如此。
我们能够通过计算闭合形式（closed-form）的方程，或者执行梯度下降，执行岭回归。闭合形式的岭回归

from sklearn.linear_model import Ridge
m = 100
X = 6 * np.random.rand(m, 1) - 3
X = np.array(sorted(list(X)))
y = 0.5 * X ** 2 + X + 2 * np.random.randn(m, 1)
# 讨论呢各个角落Cholesky求解
ridge_reg = Ridge(alpha=1, solver='cholesky')
ridge_reg.fit(X, y)
print(ridge_reg.predict([[1.5]]))
# 1.55071465
# 使用随机梯度下降法
sgd_reg = SGDRegressor(penalty='l2')
sgd_reg.fit(X, y.ravel())
print(sgd_reg.predict([[1.5]]))

Lasso回归

最小绝对收缩和选择算子回归（Least Absolute Shrinkage and Selection Operator Regression，Lasso）。和岭回归一样，通过向成本函数添加了正则项，但是添加的权重向量是范数。公式如下

Lasso回归的一个重要特点是他倾向于消除不重要的特征的权重，将他们设置为0。为了避免使用Lasso 在梯度下降时在最优解附近反弹，需要逐步降低训练期间的学习率。

from sklearn.linear_model import Lasso
lasso_reg = Lasso(alpha=0.1)
lasso_reg.fit(X, y)
print(lasso_reg.predict([[1.5]]))
from sklearn.linear_model import SGDRegressor
sgd_reg = SGDRegressor(penalty='l1')
sgd_reg.fit(X, y.ravel())
print(sgd_reg.predict([[1.5]]))

弹性网络

弹性网络介于Ridge和Lasso之间，通过控制混合比例r。当r=0的时候，弹性网络等于岭回归；当r=1的时候，弹性网络等于Lasso。成本函数如下

from sklearn.linear_model import ElasticNet
# l1_ratio即为混合比
elastic_net = ElasticNet(alpha=0.1, l1_ratio=0.5)
elastic_net.fit(X, y)
elastic_net.predict([[1.5]])

提前停止

对于梯度下降的迭代类学习算法，存在一个不同的正则化方法，即当验证误差达到最小值的时候停止训练——提前停止法。
经过一轮轮的训练，训练集上的预测误差（RMSE）会不断下降，验证集上的误差也会下降；但是达到某个点，验证误差会停止下降，开始上升，说明模型开始过拟合。因此，通过早起停止，能够在验证误差得到最小值就立刻停止。

from sklearn.base import clone
from sklearn.linear_model import Ridge
from sklearn.preprocessing import StandardScaler
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import PolynomialFeatures
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_absolute_error
X_train, X_val, y_train, y_val = train_test_split(X, y, test_size=0.2)
poly_scaler = Pipeline([
    ('Poly_features', PolynomialFeatures(degree=90, include_bias=False)),
    ('std_scaler', StandardScaler())
])
X_train_poly_scaled = poly_scaler.fit_transform(X_train)
X_val_poly_scaled = poly_scaler.transform(X_val)
sgd_reg = SGDRegressor(max_iter=1, tol=-np.infty, warm_start=True,
                       penalty=None, learning_rate='constant', eta0=0.0005)
min_val_error = float('inf')
best_epoch = None
best_model = None
for epoch in range(1000):
    sgd_reg.fit(X_train_poly_scaled, y_train)
    y_val_predict = sgd_reg.predict(X_val_poly_scaled)
    val_error = mean_absolute_error(y_val, y_val_predict)
    if val_error < min_val_error:
        min_val_error = val_error
        best_epoch = epoch
        best_model = clone(sgd_reg)

逻辑回归

回归算法也能够用于分类。逻辑回归（Logit）被广泛用于估算一个实例属于某个特定类别的概率。如果与预估概率超过了50%，那么该实例属于正类，否则是负类。从而是一个二元分类器。

估计概率

逻辑回归的估计概率

逻辑函数是一个sigmoid函数，即为S形函数，输出介于0~1之间的数字。

一旦逻辑回归模型估算出例属于正类的概率，就能够预测标签

训练和成本函数

逻辑回归模型的训练目标是为了设置参数向量，是模型对正类实例做出高概率估计（），对负类实例做出低概率估计（），例如单个训练实例的成本函数

当概率接近0的时候，将会接近，因此模型将正类实例的概率估计为0，将会导致成本很高；同样，如果将负例实例的概率估计为1，那么成本将会接近。
对于整个训练集的成本函数是所有训练实例的平均成本，能够通过损失函数表示

偏导数

from sklearn import datasets
iris = datasets.load_iris()
print(list(iris.keys()))
X = iris['data'][:, 3:]
y = (iris['target'] == 2).astype(np.int)
from sklearn.linear_model import LogisticRegression
log_reg = LogisticRegression()
log_reg.fit(X, y)
X_new = np.linspace(0, 3, 1000).reshape(-1, 1)
y_proba = log_reg.predict_proba(X_new)
plt.plot(X_new, y_proba[:, 1], "g-", label="Iris virginica")
plt.plot(X_new, y_proba[:, 0], "b--", label="Not Iris virginica")
plt.legend()
plt.show()
# + more Matplotlib code to make the image look pretty

Softmax回归

多元逻辑回归，给定一个实例，首先计算出每个类的分数，然后对这些分数应用softmax函数，估计每个类的概率。
第类的分数——每一类都有自己的指定参数向量。在得到了实例在每个类的分数之后，通过函数计算实例属于类别的概率。该函数计算每个分数的指数，然后对其归一化，除以所有指数的总和。分数通常称为对数或者对数奇数。函数如下

其中，一共有类，是一个向量、包含了实例的每个类别的分数，是实例属于类别的估计概率。
回归分类预测公式如下

argmax返回的是势函数最大化的变量值。在该等式中，返回的是是估计概率最大化的k值。
训练目标是得到一个能够对目标类做出高概率估算的模型（使其他类概率相应降低）。通过成本函数（交叉熵函数）最小化该目标

其中是类别k中第i个实例的目标概率，等于0或者1。
两个概率分布和之间的交叉熵定义为，梯度向量如下

X = iris["data"][:, (2, 3)] # petal length, petal width
y = iris["target"]
softmax_reg = LogisticRegression(multi_class="multinomial",solver="lbfgs", C=10)
softmax_reg.fit(X, y)
print(softmax_reg.predict([[5, 2]]))
# 结果为[2]