机器学习：线性回归（下）
简介在上一篇文章《机器学习：线性回归（上）》中讨论了二维数据下的线性回归及求解方法，本节中我们将进一步的将其推广至高维情形。
章节安排
[*]背景介绍
[*]最小二乘法
[*]梯度下降法
[*]程序实现
一、背景介绍
1.1 超平面\(L\)的定义
<hr>定义在\(D\)维空间中的超平面\(L\)的方程为：
\[\begin{align*}L:\text w^T\text x+b=0 \tag{1.1}\end{align*}\]
其中：\(\text w^T=\)为不同维度的系数或权重，\(\text x^T=\)为数据样本的特征向量。
在该定义中，超平面\(L\)是由是由法向量\(w\)和偏置项\(b\)决定的。具体来说，超平面\(L\)将\(D\)维空间划分为两个半空间，一个半空间满足\(\text w^T\text x+b>0\)，另一个半空间满足\(\text w^T\text x+b<0\)
，式\((1.1)\)称为矩阵表示法，也可以用标量表示法表示为：
\[\begin{align*}L:\sum_{i=1}^{D}w_ix_i+b=w_1x_1+w_2x_2+\cdots+w_Dx_D+b=0\tag{1.2}\end{align*}\]
在一些情况下，也会将偏置项\(b\)引入向量中，该方法分别对权重\(w\)和特征值\(x\)做增广：
\[\begin{align*}x^T&=\\w^T&=\end{align*}\]
在此基础上，超平面\(L\)的定义可以简化为：
\[\begin{align*}L:\text w^T\text x=0 \tag{1.3}\end{align*}\]
有时也简称
\[\begin{align*}L(\text x)=0 \tag{1.4}\end{align*}\]
示例
为方便读者理解，这里给出一个从二维的直线方程到超平面方程\(L\)的转换
\[\begin{align*}y&=kx+b\\kx-y+b&=0\\\begin{bmatrix}b&k&-1\end{bmatrix}\cdot\begin{bmatrix}1\\x\\y\end{bmatrix}&=0\end{align*}\]
1.2 高维线性回归
<hr>在高维线性回归任务中，采样数据的形式为\(S=\{\text X,\text y\}\)，其中\(X\)称为采样数据，为\(N\times D\)的矩阵，\(y\)称为标签数据，更具体的有：
\[\text X^T=[\text x_0,\text x_1, \dots, \text x_N], \text x_i=, \text x_i \in \mathbb{R}^D\]
\[\text y^T =\]
在高维数据的回归任务中，我们的目标是找到一个权重\(\text w\)，使得其能够对特征数据\(\text X\)给出预测\(\hat{\text y}\)
\[\hat{\text y}=\text X \text w\]
其中：\(\text w^T=\)是大小为\(D*1\)的向量。
同时，我们可以定义均方根误差(MSE)如下：
\[\begin{align*}\text{MSE}=\big \| \text y-\text X\text w\big\|_2^2\end{align*}\]
其中\(\|\cdot\|_2\)为二范数，或欧几里得距离。
线性回归的目标为，最小化损失，下面我们将从最小二乘法和梯度下降法两个角度实现线性回归。
二、最小二乘法
<hr>最小二乘法（Least Squares Method）是一种广泛使用的线性回归问题的求解方法，其核心思想是，均方根误差MSE关于权重\(w\)的偏导为0时所求得的\(w\)为最优解，故对MSE化简如下：
\[\begin{align*}\text{MSE}&=\big \| \text y-\text X\text w\big\|_2^2\\&=(\text y-\text X\text w)^T(\text y-\text X\text w)\\&=\text y^T\text y-\text w^T\text X^T \text y-\text y\text X\text w+\text w^T \text X^T \text X \text w\\\end{align*}\]
由于\(\text w^T\text X^T \text y\)和\(\text y\text X\text w\)是标量，其数值相等，故有：
\[\begin{align*}\text{MSE}&=\text y^T\text y-2\text w^T\text X^T \text y+\text w^T \text X^T \text X \text w\end{align*}\]
求\(\text {MSE}\)关于\(\text w\)的偏导得：
\[\begin{align*}\frac{\partial\text{MSE}}{\partial\text w}&=-2\text X^T\text y+2 \text X^T \text X \text w\end{align*}\]
另偏导等于\(0\)得：
\[\begin{align*}\text X^T\text y&= \text X^T \text X \text w \tag{2.1}\end{align*}\]
该方程称为正规方程（Normal Equation），解该方程可得：
\[\text w =(\text X^T\text X)^{-1}\text X^T \text y\]
2.1 最小二乘法缺点
以下是最小二乘法的主要缺点：
矩阵逆计算的复杂性
最小二乘法的解析解需要计算矩阵\(\text X^T \text X\) 的逆矩阵：
\[\text w = (\text X^T \text X)^{-1} \text X^T \text y \tag{2.2}\]
在高维情况下（即特征数量\(d\)较大），计算\(\text X^T \text X\) 的逆矩阵的计算复杂度很高，甚至可能不可行。具体来说：
[*]计算\(\text X^T \text X\)的时间复杂度为\(O(n d^2)\)，其中\(n\)是样本数量，\(d\)是特征数量。
[*]计算矩阵逆的时间复杂度为\(O(d^3)\)。
因此，当\(d\)很大时，计算逆矩阵的代价非常高。
矩阵不可逆问题
在高维情况下，特征数量\(d\)可能大于样本数量\(n\)，此时矩阵\(\text X^T \text X\)可能是不可逆的（即奇异矩阵），这意味着无法直接计算其逆矩阵。此外，即使矩阵可逆，也可能因为浮点数精度问题导致计算结果不稳定。
对异常值敏感
最小二乘法对异常值非常敏感。因为最小二乘法最小化的是平方误差，所以异常值会对模型的拟合产生较大的影响。这可能导致模型的泛化能力下降。
不适用于稀疏数据
对于稀疏数据（即特征矩阵中有大量零元素），最小二乘法的计算效率较低。稀疏数据通常更适合使用稀疏矩阵的优化方法，如 Lasso 或 Ridge 回归。
过拟合问题
如果没有正则化，最小二乘法容易过拟合，尤其是在特征数量远大于样本数量的情况下。过拟合会导致模型在训练集上表现很好，但在测试集上表现很差。
总结
尽管最小二乘法在许多情况下是一个简单有效的线性回归求解方法，但它也存在一些明显的缺点，特别是在高维数据和复杂情况下。为了克服这些缺点，可以考虑使用其他优化方法，如梯度下降、岭回归（Ridge Regression）、Lasso 回归等，这些方法在计算效率、对异常值的鲁棒性和防止过拟合方面有更好的表现。
三、梯度下降法
<hr>梯度下降法是一种常用的优化算法。通过迭代更新模型的参数，使得均方误差逐步减小，最终达到最优解。
对于单个样本\(\{\text x_i, y_i\}\)，其损失函数定义为：
\
求其关于权重的偏导得：
\[\begin{align*}\frac{\partial}{\partial \text w}J(\text w)&=\frac{\partial}{\partial \text w}(y-\text x_i\text w)^2\\&=2(y-\text x\text w)\text x\tag{3.1}\end{align*}\]
故有参数修正公式如下：
\[\begin{align*}\text w:=\text w -\lambda\cdot \frac{\partial J}{\partial \text w} \tag{3.2}\end{align*}\]
四、程序实现
4.1 生成测试数据
<hr>程序流程：
[*]定义特征维数feature_num及点个数point_num。
[*]定义权重向量w，特征数据X，标签数据y
[*]生成随机数，填充w和X
[*]定义误差向量error，并用随机数填充
[*]计算y
#include <iostream>#include <vector>#include <Eigen/Dense>// Multiple linear regression data generationnamespace MLR { void gen(Eigen::VectorXd& w, Eigen::MatrixXd& X, Eigen::VectorXd& y) { if (w.rows() != X.cols()) { throw std::invalid_argument("Dimension mismatch: The number of rows in w must equal the number of columns in X."); } if (X.rows() != y.rows()) { throw std::invalid_argument("Dimension mismatch: The number of rows in X must equal the number of rows in y."); } w.setRandom(); X.setRandom(); Eigen::VectorXd error(y.rows()); error.setRandom(); error *= 0.02; y = X * w + error; return; }}int main() { const size_t point_num = 10; const size_t feature_num = 7; Eigen::VectorXd w(feature_num); Eigen::MatrixXd X(point_num, feature_num); Eigen::VectorXd y(point_num); MLR::gen(w, X, y); std::cout << "y =\n" << y << "\n"; return 0;}4.2 最小二乘法实现：
<hr>程序流程：
[*]构建向量wp用以存储计算结果
[*]采用公式\((2.2)\)计算权重wp
[*]输出w-wp以观察计算误差
Eigen库中求逆、求转置都需要以矩阵为主体，例如: M.inverse()和M.transpose()。
取名wp是因为Weight prediction的首字母。
void LSM(Eigen::VectorXd& w, Eigen::MatrixXd& X, Eigen::VectorXd& y) { if (w.rows() != X.cols()) { throw std::invalid_argument("Dimension mismatch: The number of rows in w must equal the number of columns in X."); } if (X.rows() != y.rows()) { throw std::invalid_argument("Dimension mismatch: The number of rows in X must equal the number of rows in y."); } w = (X.transpose() * X).inverse() * X.transpose() * y;}int main() { // ... Eigen::VectorXd wp(feature_num); LSM(wp, X, y); std::cout << "w_error =\n" << w-wp << "\n"; return 0;}下图为程序输出结果，由该图可以看出，最小二乘法的估计较为准确。
4.3 梯度下降法实现
<hr>程序流程：
[*]构建向量wp，并初始化为随机权重。
[*]每一个数据样本x，依据公式\((3.2)\)更新一次权重。（GD_step函数功能）
[*]重复步骤2，100次。
[*]输出w-wp以观察计算误差
注意事项：
在该算法中，我们将样本的个数改为100个，即：feature_num = 100
学习率过高会导致发散，详细参考上一篇文章：《机器学习：线性回归（上）》
下式子作用是将矩阵X的第idx行读取为列向量
Eigen::VectorXd x = X.row(idx);
这与我们的使用直觉不符，实际上应为行向量。为避免出错，在后续计算中应使用x.transpose()而非直接使用x。
有一种方法可以规避该问题，即使用点积（内积）进行计算。在代码中给出了相关的示例（注释部分）
void GD_step(Eigen::VectorXd& w, Eigen::MatrixXd& X, Eigen::VectorXd& y, const double& lambda) { if (w.rows() != X.cols()) { throw std::invalid_argument("Dimension mismatch: The number of rows in w must equal the number of columns in X."); } if (X.rows() != y.rows()) { throw std::invalid_argument("Dimension mismatch: The number of rows in X must equal the number of rows in y."); } for (size_t idx = 0; idx < X.rows(); ++idx) { Eigen::VectorXd x = X.row(idx); // 使用点积 // Eigen::VectorXd gradient = 2 * (y(idx) - x.dot(w)) * x; // 因为 y-x*w是标量，且输出结果为VectorXd，因此最后的transpose是可去的。 // Eigen::VectorXd gradient = 2 * (y(idx) - x.transpose() * w) * x.transpose(); Eigen::VectorXd gradient = 2 * (y(idx) - x.transpose() * w) * x; w += lambda * gradient; }}int main() { const size_t point_num = 100; // ... Eigen::VectorXd wp(feature_num); wp.setRandom(); // 生成初始值 double lambda = 2e-3; for (int _ = 0; _ < 100; ++_) { GD_step(wp, X, y, lambda); } std::cout << "w_error =\n" << w - wp << "\n"; return 0;}下图为程序输出结果，由该图可以看出，梯度下降法的估计较为准确。
页:
[1]