• 设为首页
  • 点击收藏
  • 手机版
    手机扫一扫访问
    迪恩网络手机版
  • 关注官方公众号
    微信扫一扫关注
    迪恩网络公众号

[线性模型] 线性回归(Linear Regression)原理及MATLAB实现

原作者: [db:作者] 来自: [db:来源] 收藏 邀请

线性回归解决的问题

“线性回归” 试图学得一个通过属性的线性组合来进行预测的函数,以尽可能准确地预测实值输出标记,一般形式为

\[f(\boldsymbol{x})=\boldsymbol{w}^T\boldsymbol{x}+b \tag 1 \]

其中 \(\boldsymbol{x}\) 表示一组属性,长度为 \(n\) 的列向量. \(\boldsymbol{w}=(w_1;w_2;w_3;...;w_n)\) 表示一组参数,长度为 \(n\) 的列向量,每个 \(w_i\) 可以理解为属性 \(x_i\) 的 “重要程度”. \(b\) 也为参数,可以理解为对 \(\boldsymbol{w}^T\boldsymbol{x}\) 的值进行修正. \(f(\boldsymbol{x})\) 表示通过函数算得的预测值.

比如,对同一个地区的父亲和儿子的身高进行记录,把所有父亲的身高作为 \(\boldsymbol{x}\),所有儿子的身高作为 \(\boldsymbol{y}\)通过一些数学方法,解出一组 \(\boldsymbol{w}\)\(b\) 使得对于每一个输入 \(x_i\),尽可能满足

\[f(x_i) = wx_i+b\approx y_i \tag 2 \]

之后,对于任何一个不在 \(\boldsymbol{x}\) 中的父亲身高,都可以使用函数预测出其对应的儿子身高.

显然,该模型预测的准确程度取决于训练集在坐标轴上的分布,如果大致呈线性的话,使用线性回归可以得到较好的结果. 因此,在一开始,可以输出训练集的坐标图,粗略地看一下,再决定是否使用线性回归训练。

\(\boldsymbol{w}\)\(b\) 的求解过程

常用的解决思路有两种:

  1. 使用最小二乘法进行参数估计,求得参数的闭式解;
  2. 给定一个合适的学习率,使用梯度下降法,求出参数的数值解。

最小二乘法

首先,要定一个优化目标,考虑式 \((2)\),回归任务中常使用均方误差作为损失函数,并使其最小化,可以表示为

\[\begin{align} E_{(w,b)} &= \sum\limits_{i=1}^{m}(f(x_i)-y_i)^2 \\ & = \sum\limits_{i=1}^{m}(y_i-wx_i-b)^2 \tag 3\\ (w^*, b^*) &= \arg\min\limits_{(w,b)}E_{(w,b)} \tag 4\\ \end{align} \]

该公式可以理解为,使所有样本到直线上的欧式距离和最小. 求解 \(w\)\(b\) 使 \(E_{(w,b)}\) 最小化的过程,可以使用概率论中参数估计的方法. 将 \(E_{(w,b)}\)\(w\)\(b\) 分别求导,得:

\[\begin{align} \frac{\partial E_{(w,b)}}{\partial w} &= 2\left( \sum\limits_{i=1}^{m}x_i^2- \sum\limits_{i=1}^{m}(y_i-b)x_i\right) \tag 5\\ \frac{\partial E_{(w,b)}}{\partial b} &= 2\left(mb- \sum\limits_{i=1}^{m}(y_i-wx_i)\right) \tag 6\\ \end{align} \]

此处 \(E_{(w,b)}\) 为关于 \(w\)\(b\) 的凸函数,因此当式 \((5)(6)\) 都为 \(0\) 时,得到参数 \(w\)\(b\) 的最优的闭式解

\[\begin{align} w &= \frac{\sum\limits_{i=1}^{m}y_i(x_i-\bar x)}{\sum\limits_{i=1}^{m}x_i^2 - \frac{1}{m}\left(\sum\limits_{i=1}^{m}x_i\right)^2}\ = \frac{\sum\limits_{i=1}^{m}(x_i-\bar x)(y_i-\bar y)}{\sum\limits_{i=1}^{m}(x_i-\bar x)^2} \tag 7\\ b &= \frac{1}{m}\sum\limits_{i=1}^{m}(y_i-wx_i) \tag 8\\ \end{align} \]

\(\boldsymbol x_d = (x_1-\bar x, x_2-\bar x, ..., x_m-\bar x), \boldsymbol y_d= (y_1-\bar y, y_2-\bar y, ..., y_m-\bar y)\), 式 \((7)\) 可以改写为

\[w = \frac{\boldsymbol x_d^T \boldsymbol y_d}{\boldsymbol x_d^T \boldsymbol x_d} \tag 9\\ \]

这些表示是为了与 ”多元线性回归“ 中 \(\boldsymbol w\) 的解对应. 接下来,我们讨论求解一开始的线性模型,即属性 \(\boldsymbol x\) 为一个列向量,而不是单一的数,大多数情况下线性规划处理的也是这样的问题。

求解 “多元线性回归” 中参数的过程与上述类似,仍旧采用最小二乘法。令 \(\hat{\boldsymbol w} = (\boldsymbol w;b), \boldsymbol y=(y_1;y_2;...;y_m)\),将所有的属性向量构造成一个 \(m(n+1)\) 大小的矩阵 \(\boldsymbol X\), 表示为

\[\boldsymbol X = \begin{pmatrix}        x_{11} & x_{12} & \cdots & x_{1n} &1\\        x_{21} & x_{22} & \cdots & x_{2n} &1\\       \vdots & \vdots & \ddots & \vdots & \vdots\\        x_{m1} & x_{m2} & \cdots & x_{mn} &1 \\        \end{pmatrix} \]

类似式 \((3)(4)\) 此时均方误差可以表示为

\[\begin{align} E_{\hat{\boldsymbol w}} &= (\boldsymbol y - \boldsymbol X\hat{\boldsymbol w})^T(\boldsymbol y - \boldsymbol X\hat{\boldsymbol w}) \tag {11}\\ \hat{\boldsymbol w}^* &= \arg\min\limits_{\hat{\boldsymbol w}^*}E_{\hat{\boldsymbol w}} \tag {10} \end{align} \]

\(E_{\hat{\boldsymbol w}}\)\(\hat{\boldsymbol w}\) 求偏导,得

\[\begin{align} \frac{\partial E_{\hat{\boldsymbol w}}}{\partial \hat{\boldsymbol w}} &= \frac{\partial \boldsymbol y^T\boldsymbol y}{\partial \hat{\boldsymbol w}} +\frac{\partial \boldsymbol y^T\boldsymbol X \hat{\boldsymbol w}}{\partial \hat{\boldsymbol w}} + \frac{\partial \hat{\boldsymbol w}^T \boldsymbol X^T \boldsymbol y}{\partial \hat{\boldsymbol w}} + \frac{\partial \hat{\boldsymbol w}^T \boldsymbol X^T \boldsymbol X \hat{\boldsymbol w}}{\partial \hat{\boldsymbol w}}\\ & = 0-\boldsymbol X^T \boldsymbol y -\boldsymbol X^T \boldsymbol y +(\boldsymbol X^T \boldsymbol X + \boldsymbol X^T \boldsymbol X) \hat{\boldsymbol w}\\ & = 2\boldsymbol X^T (\boldsymbol X\hat{\boldsymbol w}-\boldsymbol y) \tag {11}\\ \end{align} \]

令其等于 \(0\), 当 \(\boldsymbol X^T \boldsymbol X\) 为满值矩阵或正定矩阵时,解出 \(\hat{\boldsymbol w}\)

\[\hat{\boldsymbol w} = (\boldsymbol X^T \boldsymbol X)^{-1} \boldsymbol X^T \boldsymbol y \tag{12} \]

与式 \((9)\) 拥有类似的结构. 当 \(\boldsymbol X^T \boldsymbol X\) 不满秩时,通常引入正则化项,决定其中一组 \(\hat{\boldsymbol w}\) 作为最优解.

梯度下降

对于式 \((11)\) 可以采用一种迭代求解的方法,即梯度下降法. 首先随机给定 \(\hat{\boldsymbol w}\) 的一组初始解,可以是零向量. 然后,让 \(\hat{\boldsymbol w}\) 沿着梯度下降最快的方向递减,当偏导数为 \(0\) 时,得到 \(\hat{\boldsymbol w}\) 的最优解. 其中 \(\hat{\boldsymbol w}\) 的迭代式为

\[\hat{\boldsymbol w} \leftarrow \hat{\boldsymbol w} - \alpha\frac{\partial E_{\hat{\boldsymbol w}}}{\partial \hat{\boldsymbol w}}\tag {13} \]

其中 \(\alpha\) 表示为 \(\hat{\boldsymbol w}\)沿着梯度下降最快方向递减的速率,也即学习率. 将式 \((11)\) 代入式 \((13)\) 并令 \(\alpha \leftarrow 2\alpha\),得

\[\begin{align} \hat{\boldsymbol w} \leftarrow \hat{\boldsymbol w} - \alpha\boldsymbol X^T (\boldsymbol X\hat{\boldsymbol w}-\boldsymbol y)\tag {14}\\ \hat{\boldsymbol w}_j \leftarrow \hat{\boldsymbol w}_j - \alpha\boldsymbol x_j^T (\boldsymbol X\hat{\boldsymbol w}-\boldsymbol y)\tag {15}\\ \end{align} \]

每一次迭代使用训练集的所有数据的算法称为批量梯度下降法,一次迭代的时间复杂度为 \(O(nm)\). 因此,当训练集样本数和为维数都很大时,得出一组解需要大量的时间。针对这样的情况,引入随机梯度下降法.

每次更新 \(\hat{\boldsymbol w}\) 仅使用一个训练样本,然后判断是否收敛. 当有几万条以上的训练样本时,采用随机梯度下降会快很多,其伪代码如下

do until Converge 
    foreach x in X 
        foreach j in w 
            w[j] = w[j] - a*x[j]*(x[i]*w-y)

判断收敛方法一般由两种:

  1. \(\hat{\boldsymbol w}\) 中每个参数的变化均小于一个阈值;
  2. \(E_{\hat{\boldsymbol w}}\) 的变化小于一个阈值;

使用随机梯度下降,虽然速度快,但是收敛性能可能会不太好. 综合批量梯度下降法和随机梯度下降法,可每次迭代随机从训练集中挑选一定数量的样本训练,即小批量梯度下降法.

MATLAB 实现

代码采用批量梯度下降法,使用随机生成在坐标系中线性分布的 \(100\) 条训练样本.

% 生成随机训练样本,大致为 y=0.7x+12;
kb = [0.7, 12];
x = zeros(100,2);
y = zeros(100,1);

for i = 1:100
    x(i,1) = randi(1000,1);
    x(i,2) = 1;
end
for i = 1:100
    y(i) = kb(1)*x(i)+kb(2) + randi(200,1)-100;
end
hold on;
% w: 权重
% X: 训练集属性
% alpha: 学习率
% eps: 阈值
% m: 训练集样本数
% n: 训练集维数
function w = BDG(X, y, alpha, eps)
    [m,n] = size(X);
    X = [X ones(m,1)];
    n = n+1;
    w = zeros(n,1);
    prew = zeros(n,1);
    while(1)
        flag = 0;
        w(:) = prew(:) - alpha*X(:,:)\'*(X(:,:)*prew(:)-y(:));
        for j = 1:n
            if abs(w(j)-prew(j))>eps
                flag = 1;
                break;
            end
       end
       prew = w; 
       if flag==0
           break;
       end
    end
end
% 测试并输出图像
alpha = 0.00000001;
eps = 0.01;
w = BDG(x, y, alpha, eps);
p = zeros(1000,1);
q = zeros(1000,1);
for j = 1:1000
    p(j,1) = j;
    p(j,2) = j;
    q(j) = w(1)*p(j,1)+w(2);
end
plot(x(:,1), y, \'*\');
hold on;
plot(p(:,1),  q, \'-r\');

直线如图所示

参考文献

  1. 周志华《机器学习》
  2. 南瓜书 https://github.com/datawhalechina/pumpkin-book
  3. 线性回归与梯度下降算法 https://www.cnblogs.com/eczhou/p/3951861.html

鲜花

握手

雷人

路过

鸡蛋
该文章已有0人参与评论

请发表评论

全部评论

专题导读
上一篇:
Delphi DbGridEh实现表格没有内容的渐变效果发布时间:2022-07-18
下一篇:
基于matlab的升压斩波电路仿真发布时间:2022-07-18
热门推荐
阅读排行榜

扫描微信二维码

查看手机版网站

随时了解更新最新资讯

139-2527-9053

在线客服(服务时间 9:00~18:00)

在线QQ客服
地址:深圳市南山区西丽大学城创智工业园
电邮:jeky_zhao#qq.com
移动电话:139-2527-9053

Powered by 互联科技 X3.4© 2001-2213 极客世界.|Sitemap