机器学习训练算法十一(模型训练算法-Matlab实验)

在机器学习中经常用梯度下降法、牛顿下降法、高斯牛顿法、列文伯格-马夸尔特法来训练模型，为了方便理解，本文严谨且详细的阐述了这几个算法的数学推导原理，并且用matlab程序加以实验验证。

1、测试数据

为了方便研究最小二乘法问题,现提供如下车辆时间和行驶距离的观测数据用于讨论和分析。

2、数学模型

令：车辆的初速度$v_0$ 、加速度$a$、时间$t$、距离 $\hat{y}$ ，第$i$组观测值为$(t_i,{\hat{y}}_i)$，时间$t$与距离$\hat{y}$满足的数学模型如下：
$$\hat{y}=f(v_0,a,t)=\frac{1}{2}a{t}^2+v_{0}t(公式59)$$$$\text{min}F(a,v_0)=\frac{1}{2}\sum _{i=1}^{10}(f(v_0,a,t_i)-\widehat{y_i}{)}^2(公式60)$$
观察公式 60 发现， $\widehat{y_i}$ 和 $t_i$ 是 常 数 , $a$ 和 $v_0$ 是变量，可知该最小二乘问题是函数$F$关于 $a$ 和 $v$ 的最小极值问题。

3、梯度下降法

3.1、数学原理推导

链接地址

3.2、Matlab程序

% ====================================================
% 梯度下降算法
% ====================================================
f = @(args,t)(0.5*args(1)*power(t,2)+args(2)*t); % 定义数学模型
l = @(args,t,y)((f(args,t)-y)); % 定义残差函数
o = @(lv)(0.5.*power(lv,2)); % 定义目标函数
% ======================偏导函数========================
delta = 0.0001;
oj1 = @(args,t,y)((o(l(args+[delta 0]',t,y))-o(l(args,t,y)))./delta); % 关于目标函数第一个参数的一阶偏导函数
oj2 = @(args,t,y)((o(l(args+[0 delta]',t,y))-o(l(args,t,y)))./delta); % 关于目标函数第二个参数的一阶偏导函数
% ======================初始数据========================
t = [0 1 2 3 4 5 6 7 8 9 10]; % 加速度实验的观测值,时间(单位:秒)
y = [0 11.52 26.2 43.5 64.12 87.57 114.12 143.5 176.3 211.5 250.12]; % 加速度实验的观测值,距离(单位:米)
x = [1;1]; % 初始化参数[args(1);args(2)]
goal = 0.0001; % 目标误差
epochs = 10000; % 训练次数
ov = 0; % 初始化目标值
for k=1:epochs
    % ======================误差计算========================
    tl = l(x,t,y); % 计算残差值
    tmp_ov = sum(o(tl)); % 计算最新目标值
    tmp_error = abs(tmp_ov-ov); % 目标值的变化
    disp([num2str(k),'->目标值:[ ',num2str(tmp_ov),']误差:[',num2str(tmp_error),']']);
    ov = tmp_ov; % 更新目标值
    if tmp_error<goal
        disp(['训练次数: ',num2str(k),'次...']);
        break;
    end
    % ======================执行训练========================
    J=[sum(oj1(x,t,y))
       sum(oj2(x,t,y))]; % 雅克比矩阵
    % ======================更新参数========================
    x=x-0.0001*J;
end
% ========================================================
% 打印训练结果
% ========================================================
disp(['训练后的参数: ',mat2str(x)]);
plot(t,y,'o',min(t):0.001:max(t),f(x,min(t):0.001:max(t)),'-');

3.3、执行日志

2301->目标值:[ 0.066603]误差:[0.00010549]
2302->目标值:[ 0.066498]误差:[0.000105]
2303->目标值:[ 0.066393]误差:[0.00010451]
2304->目标值:[ 0.066289]误差:[0.00010403]
2305->目标值:[ 0.066186]误差:[0.00010355]
2306->目标值:[ 0.066083]误差:[0.00010307]
2307->目标值:[ 0.06598]误差:[0.0001026]
2308->目标值:[ 0.065878]误差:[0.00010212]
2309->目标值:[ 0.065776]误差:[0.00010165]
2310->目标值:[ 0.065675]误差:[0.00010118]
2311->目标值:[ 0.065575]误差:[0.00010071]
2312->目标值:[ 0.065474]误差:[0.00010025]
2313->目标值:[ 0.065375]误差:[9.9784e-05]
训练次数: 2313次...
训练后的参数: [3.00520620930594;9.99137523198286]

3.4、训练结果

4、牛顿下降法

4.1、数学原理推导

链接地址

4.2、Matlab程序

% ====================================================
% Newton算法
% ====================================================
f = @(args,t)(0.5*args(1)*power(t,2)+args(2)*t); % 定义数学模型
l = @(args,t,y)((f(args,t)-y)); % 定义残差函数
o = @(lv)(0.5.*power(lv,2)); % 定义目标函数
% ======================偏导函数========================
delta = 0.0001;
oj1 = @(args,t,y)((o(l(args+[delta 0]',t,y))-o(l(args,t,y)))./delta); % 关于目标函数第一个参数的一阶偏导函数
oj2 = @(args,t,y)((o(l(args+[0 delta]',t,y))-o(l(args,t,y)))./delta); % 关于目标函数第二个参数的一阶偏导函数
oh00 = @(args,t,y)((oj1(args+[delta 0]',t,y)-oj1(args,t,y))/delta); % 关于目标函数的(第一个参数,第一个参数)的二阶偏导函数
oh01 = @(args,t,y)((oj1(args+[0 delta]',t,y)-oj1(args,t,y))/delta); % 关于目标函数的(第一个参数,第二个参数)的二阶偏导函数
oh10 = @(args,t,y)((oj2(args+[delta 0]',t,y)-oj2(args,t,y))/delta); % 关于目标函数的(第二个参数,第一个参数)的二阶偏导函数
oh11 = @(args,t,y)((oj2(args+[0 delta]',t,y)-oj2(args,t,y))/delta); % 关于目标函数的(第二个参数,第二个参数)的二阶偏导函数
% ======================初始数据========================
t = [0 1 2 3 4 5 6 7 8 9 10]; % 加速度实验的观测值,时间(单位:秒)
y = [0 11.52 26.2 43.5 64.12 87.57 114.12 143.5 176.3 211.5 250.12]; % 加速度实验的观测值,距离(单位:米)
x = [1;1]; % 初始化参数[args(1);args(2)]
goal = 0.0001; % 目标误差
epochs = 10000; % 训练次数
ov = 0; % 初始化目标值
for k=1:epochs
    % ======================误差计算========================
    tl = l(x,t,y); % 计算残差值
    tmp_ov = sum(o(tl)); % 计算最新目标值
    tmp_error = abs(tmp_ov-ov); % 目标值的变化
    disp([num2str(k),'->目标值:[ ',num2str(tmp_ov),']误差:[',num2str(tmp_error),']']);
    ov = tmp_ov; % 更新目标值
    if tmp_error<goal
        disp(['训练次数: ',num2str(k),'次...']);
        break;
    end
    % ======================执行训练========================
    H=[sum(oh00(x,t,y)),sum(oh01(x,t,y))
       sum(oh10(x,t,y)),sum(oh11(x,t,y))]; % 黑塞矩阵
    J=[sum(oj1(x,t,y))
       sum(oj2(x,t,y))]; % 雅克比矩阵
    % ======================更新参数========================
    x=x-inv(H)*J;
end
% ========================================================
% 打印训练结果
% ========================================================
disp(['训练后的参数: ',mat2str(x)]);
plot(t,y,'o',min(t):0.001:max(t),f(x,min(t):0.001:max(t)),'-');

4.3、执行日志

1->目标值:[ 55574.2293]误差:[55574.2293]
2->目标值:[ 0.044544]误差:[55574.1847]
3->目标值:[ 0.04453]误差:[1.4246e-05]
训练次数: 3次...
训练后的参数: [2.9945866146396;10.0356844035264]

4.4、训练结果

5、高斯牛顿法

5.1、数学原理推导

链接地址

5.2、Matlab程序

% ====================================================
% Gauss-Newton算法
% ====================================================
f = @(args,t)(0.5*args(1)*power(t,2)+args(2)*t); % 定义数学模型
l = @(args,t,y)((f(args,t)-y)); % 定义残差函数
o = @(lv)(0.5.*power(lv,2)); % 定义目标函数
% ======================偏导函数========================
delta = 0.0001;
j1 = @(args,t,y)((l(args+[delta 0]',t,y)-l(args,t,y))./delta); % 关于残差函数第一个参数的一阶偏导函数
j2 = @(args,t,y)((l(args+[0 delta]',t,y)-l(args,t,y))./delta); % 关于残差函数第二个参数的一阶偏导函数
% ======================初始数据========================
t = [0 1 2 3 4 5 6 7 8 9 10]; % 加速度实验的观测值,时间(单位:秒)
y = [0 11.52 26.2 43.5 64.12 87.57 114.12 143.5 176.3 211.5 250.12]; % 加速度实验的观测值,距离(单位:米)
x = [1;1]; % 初始化参数[args(1);args(2)]
goal = 0.0001; % 目标误差
epochs = 10000; % 训练次数
ov = 0; % 初始化目标值
for k=1:epochs
    % ======================误差计算========================
    tl = l(x,t,y); % 计算残差值
    tmp_ov = sum(o(tl)); % 计算最新目标值
    tmp_error = abs(tmp_ov-ov); % 目标值的变化
    disp([num2str(k),'->目标值:[ ',num2str(tmp_ov),']误差:[',num2str(tmp_error),']']);
    ov = tmp_ov; % 更新目标值
    if tmp_error<goal
        disp(['训练次数: ',num2str(k),'次...']);
        break;
    end
    % ======================执行训练========================
    tmp_j1 = j1(x,t,y); % 关于第一个参数的一阶偏导函数的值
    tmp_j2 = j2(x,t,y); % 关于第二个参数的一阶偏导函数的值
    H=[sum(tmp_j1.*tmp_j1),sum(tmp_j1.*tmp_j2)
       sum(tmp_j2.*tmp_j1),sum(tmp_j2.*tmp_j2)]; % 近似黑塞矩阵
    J=[sum(tl.*tmp_j1)
       sum(tl.*tmp_j2)]; % 近似雅克比矩阵
    % ======================更新参数========================
    x=x-inv(H)*J;
end
% ========================================================
% 打印训练结果
% ========================================================
disp(['训练后的参数: ',mat2str(x)]);
plot(t,y,'o',min(t):0.001:max(t),f(x,min(t):0.001:max(t)),'-');

5.3、执行日志

1->目标值:[ 55574.2293]误差:[55574.2293]
2->目标值:[ 0.044455]误差:[55574.1848]
3->目标值:[ 0.044455]误差:[1.4919e-15]
训练次数: 3次...
训练后的参数: [2.9952026286965;10.0333143483027]

5.4、训练结果

6、列文伯格-马夸尔特法(LM 法)

6.1、数学原理推导

链接地址

6.2、Matlab程序

% ====================================================
% LM算法
% ====================================================
f = @(args,t)(0.5*args(1)*power(t,2)+args(2)*t); % 定义数学模型
l = @(args,t,y)((f(args,t)-y)); % 定义残差函数
o = @(lv)(0.5.*power(lv,2)); % 定义目标函数
% ======================偏导函数========================
delta = 0.0001;
j1 = @(args,t,y)((l(args+[delta 0]',t,y)-l(args,t,y))./delta); % 关于残差函数第一个参数的一阶偏导函数
j2 = @(args,t,y)((l(args+[0 delta]',t,y)-l(args,t,y))./delta); % 关于残差函数第二个参数的一阶偏导函数
% ======================初始数据========================
t = [0 1 2 3 4 5 6 7 8 9 10]; % 加速度实验的观测值,时间(单位:秒)
y = [0 11.52 26.2 43.5 64.12 87.57 114.12 143.5 176.3 211.5 250.12]; % 加速度实验的观测值,距离(单位:米)
x = [1;1]; % 初始化参数[args(1);args(2)]
goal = 0.0001; % 目标误差
epochs = 10000; % 训练次数
ov = 0; % 初始化目标值
for k=1:epochs
    % ======================误差计算========================
    tl = l(x,t,y); % 计算残差值
    tmp_ov = sum(o(tl)); % 计算最新目标值
    tmp_error = abs(tmp_ov-ov); % 目标值的变化
    disp([num2str(k),'->目标值:[ ',num2str(tmp_ov),']误差:[',num2str(tmp_error),']']);
    ov = tmp_ov; % 更新目标值
    if tmp_error<goal
        disp(['训练次数: ',num2str(k),'次...']);
        break;
    end
    % ======================执行训练========================
    tmp_j1 = j1(x,t,y); % 关于第一个参数的一阶偏导函数的值
    tmp_j2 = j2(x,t,y); % 关于第二个参数的一阶偏导函数的值
    H=[sum(tmp_j1.*tmp_j1),sum(tmp_j1.*tmp_j2)
       sum(tmp_j2.*tmp_j1),sum(tmp_j2.*tmp_j2)]; % 近似黑塞矩阵
    J=[sum(tl.*tmp_j1)
       sum(tl.*tmp_j2)]; % 近似雅克比矩阵
    % ======================计算值λ========================
    if k==1
        % 初始化值λ
        lambda = 0.000001*max(diag(H));
        v=2;
    else
        % 更新值λ
        tmp_delta = delta*ones(size(J));
        beta = (tmp_ov-o(l(x+tmp_delta,t,y)))/(tmp_ov-(tmp_ov+J'*tmp_delta+0.5.*tmp_delta'*H*tmp_delta));
        if beta>0
            lambda = lambda * max([1/3;1-power(2*beta-1,3)]);
            v = 2;
        else
            lambda = lambda*v;
            v = 2*v;
        end
    end
    disp(['-->值λ: ',mat2str(lambda)]);
    % ======================更新参数========================
     x=x-inv(H+lambda*eye(size(H)))*J;
end
% ========================================================
% 打印训练结果
% ========================================================
disp(['训练后的参数: ',mat2str(x)]);
plot(t,y,'o',min(t):0.001:max(t),f(x,min(t):0.001:max(t)),'-');

6.3、执行日志

1->目标值:[ 55574.2293]误差:[55574.2293]
-->值λ: 0.00633324999999457
2->目标值:[ 0.044517]误差:[55574.1847]
-->值λ: 0.0126664999999891
3->目标值:[ 0.044455]误差:[6.1711e-05]
训练次数: 3次...
训练后的参数: [2.99520293572588;10.033313067223]

6.4、训练结果

7、方法比较