PID神经网络原理与MATLAB实现（SISO）

最新推荐文章于 2024-06-11 19:08:37 发布

大强强小强强

最新推荐文章于 2024-06-11 19:08:37 发布

阅读量9.4k

点赞数 26

分类专栏：算法

本文链接：https://blog.csdn.net/weixin_41869763/article/details/112849902

版权

算法专栏收录该内容

8 篇文章 21 订阅

订阅专栏

序

最近想实现一下PID神经网络，但是书籍和博客都令人头疼，主要是卡在误差反向传播的计算过程中。找一篇通俗易懂的文章实在不易，最终，只能自己静下心，仔细琢磨。只要每一步在逻辑上都是合理的，我们有理由相信能够得到正确的结果。抱着这样的心态，由浅入深，来实现一下。

一、网络结构定义

简单起见，假设一个受控系统单输入 $u$ 单输出 $x$ ，使用一个PID神经网络来作为控制器，使得系统输出达到目标值 $z$ 。

网络结构定义

其中， $W_{3\times 2}$ 和 $V_{1\times 3}$ 分别为输入层—隐含层权重矩阵和隐含层—输出层权重矩阵。 $y_i$ 和 $y_o$ 分别是PID计算前后向量。 $z, x, u$ 分别为期望值、实际值和网络输出，都是数。

1.1 前向传播

根据以上结构，容易得出前向传播如下。注意简介起见，第 $(k)$ 次的标注省略。

$\begin{aligned} & y_{i1} = w_{11} z + w_{12} x \\ & y_{i2} = w_{21} z + w_{22} x \\ & y_{i3} = w_{31} z + w_{32} x \tag{1} \end{aligned}$

$\begin{aligned} & y_{o1} = y_{i1} \\ & y_{o2} = y_{o2}(k-1) + y_{i2} \\ & y_{o3} = y_{i3}(k) - y_{i3}(k-1) \tag{2} \end{aligned}$

$v_{11} y_{o1} + v_{12} y_{o2} + v_{13} y_{o3} \tag{3}$

不妨设受控对象为一个类似二次函数的时变系统（用于仿真），则输出及损失函数为

$u^2(k) + 0.1x(k-1) \tag{4}$

$\frac{1}{2} (z-x)^2 \tag{5}$

1.2 反向传播

前向传播计算完成，以下开始反向传播。首先考虑权重 $V$ 对误差 $J$ 的影响，也即计算 $\frac{\partial J}{\partial V}$ ，根据前向传递函数关系，易得

$\begin{aligned} \frac{\partial J}{\partial v_{11}} &= \frac{\partial J}{\partial x} \frac{\partial x}{\partial u} \frac{\partial u}{\partial v_{11}} \\ &= -(z-x) \cdot \frac{\partial x}{\partial u}\cdot y_{o1} \end{aligned} \tag{6}$

由于受控对象的数学模型往往不知道或者偏导难以计算，中间的 $\frac{\partial x}{\partial u}$ 常用下式子代替

$\frac{\partial x}{\partial u} ={\rm sign} \frac{x(k) - x(k-1)}{u(k) - u(k-1)} \tag{7}$

以下 $\frac{\partial x}{\partial u}$ 皆由式(7)得出，不在赘述。

其中 ${\rm sign}(\cdot)$ 是符号函数，输入正数输出 $1$ ，输入负数输出 $- 1$ ，输入 $0$ 输出 $0$ 。使用符号函数的好处是归一化到 $[- 1, 1]$ ，避免分母很接近0导致输出很大的问题。

同理易得

$\frac{\partial J}{\partial v_{12}} = -(z-x) \cdot \frac{\partial x}{\partial u} \cdot y_{o2} \tag{8}$

$\frac{\partial J}{\partial v_{13}} = -(z-x) \cdot \frac{\partial x}{\partial u} \cdot y_{o3} \tag{9}$

接下来，计算真正的难点 $\frac{\partial J}{\partial W}$ 。

$\begin{aligned} \frac{\partial J}{\partial w_{11}} &= \frac{\partial J}{\partial x} \frac{\partial x}{\partial u} \frac{\partial u}{\partial y_{o1}} \frac{\partial y_{o1}}{\partial y_{i1}} \frac{\partial y_{i1}}{\partial w_{11}} \\ &=-(z-x) \cdot\frac{\partial x}{\partial u} \cdot v_{11} \cdot 1 \cdot z \end{aligned} \tag{10}$

同理易得

$\begin{aligned} \frac{\partial J}{\partial w_{12}} &= \frac{\partial J}{\partial x} \frac{\partial x}{\partial u} \frac{\partial u}{\partial y_{o1}} \frac{\partial y_{o1}}{\partial y_{i1}} \frac{\partial y_{i1}}{\partial w_{12}} \\ &=-(z-x) \cdot\frac{\partial x}{\partial u} \cdot v_{11} \cdot 1 \cdot x \end{aligned} \tag{11}$

$\begin{aligned} \frac{\partial J}{\partial w_{21}} &= \frac{\partial J}{\partial x} \frac{\partial x}{\partial u} \frac{\partial u}{\partial y_{o2}} \frac{\partial y_{o2}}{\partial y_{i2}} \frac{\partial y_{i2}}{\partial w_{21}} \\ &=-(z-x) \cdot \frac{\partial x}{\partial u} \cdot v_{12} \cdot \frac{\partial y_{o2}}{\partial y_{i2}} \cdot z \end{aligned} \tag{12}$

其中

$\frac{\partial y_{o2}}{\partial y_{i2}} = \frac{\Delta y_{o2}}{\Delta y_{i2}} = \frac{y_{o2}(k) - y_{o2}(k-1)}{y_{i2}(k) - y_{i2}(k-1)} = \frac{y_{i2}(k) }{y_{i2}(k) - y_{i2}(k-1)}$

和之前一样，使用符号函数计算可得

$\frac{\partial y_{o2}}{\partial y_{i2}} = {\rm sign}\frac{y_{i2}(k) }{y_{i2}(k) - y_{i2}(k-1)} \tag{13}$

同理

$\begin{aligned} \frac{\partial J}{\partial w_{22}} &= \frac{\partial J}{\partial x} \frac{\partial x}{\partial u} \frac{\partial u}{\partial y_{o2}} \frac{\partial y_{o2}}{\partial y_{i2}} \frac{\partial y_{i2}}{\partial w_{22}} \\ &=-(z-x) \cdot \frac{\partial x}{\partial u} \cdot v_{12} \cdot \frac{\partial y_{o2}}{\partial y_{i2}} \cdot x \end{aligned} \tag{14}$

$\begin{aligned} \frac{\partial J}{\partial w_{31}} &= \frac{\partial J}{\partial x} \frac{\partial x}{\partial u} \frac{\partial u}{\partial y_{o3}} \frac{\partial y_{o3}}{\partial y_{i3}} \frac{\partial y_{i3}}{\partial w_{31}} \\ &=-(z-x) \cdot \frac{\partial x}{\partial u} \cdot v_{13} \cdot \frac{\partial y_{o3}}{\partial y_{i3}} \cdot z \end{aligned} \tag{15}$

其中

$\frac{\partial y_{o3}}{\partial y_{i3}} = \frac{\Delta y_{o3}}{\Delta y_{i3}} = \frac{y_{o3}(k) - y_{o3}(k-1)}{y_{i3}(k) - y_{i3}(k-1)} = \frac{y_{o3}(k) - y_{o3}(k-1)}{y_{o3}(k)}$

依旧使用符号函数归一化

$\frac{\partial y_{o3}}{\partial y_{i3}} ={\rm sign} \frac{y_{o3}(k) - y_{o3}(k-1)}{y_{o3}(k)}$

同理

$\begin{aligned} \frac{\partial J}{\partial w_{32}} &= \frac{\partial J}{\partial x} \frac{\partial x}{\partial u} \frac{\partial u}{\partial y_{o3}} \frac{\partial y_{o3}}{\partial y_{i3}} \frac{\partial y_{i3}}{\partial w_{32}} \\ &=-(z-x) \cdot \frac{\partial x}{\partial u} \cdot v_{13} \cdot \frac{\partial y_{o3}}{\partial y_{i3}} \cdot x \end{aligned} \tag{16}$

1.3 权值更新

如果使用传统的梯度下降法，则更新公式为

$\alpha_1 \frac{\partial J}{\partial W} \tag{17}$

$\alpha_2 \frac{\partial J}{\partial V} \tag{18}$

其中， $\alpha_1, \alpha_2$ 为学习率。

有时，也增加动量项，则更新公式为

$\begin{aligned} & W_d(k) = \beta_1 W_d(k-1) + (1-\beta_1)\frac{\partial J}{\partial W}\\ & W = W - \alpha_1 W_d \tag{19} \end{aligned}$

$\begin{aligned} & V_d(k) = \beta_2 V_d(k-1) + (1-\beta_2)\frac{\partial J}{\partial V}\\ & V = V - \alpha_2 V_d \tag{20} \end{aligned}$

其中， $\beta_1, \beta_2$ 为动量因子，一般取 0.9； $W_d, V_d$ 可取初始值为0进行迭代。

1.4 实际问题

以上算法纯属理论，在仿真中可能都会遇到问题。一个典型的问题是除数为0的情况，在这里很容易发生，例如 $u (k) = u (k - 1)$ 就是如此。另一个问题是数值过大，比如除数很小就很容易导致这种现象。因此有必要规避除数为0，也要做好限幅工作。写程序时需要注意这些细节。

二、仿真

现假设受控系统是一个式(4) 描述的系统，设期望值为一个阶梯信号，整个程序如下

%% 参数初始化
% 初始化
rng(1);                         % 设置随机数种子
z = 1;                          % 目标值
x = 0;                          % 实际值初始化
W = 0.3*(rand(3,2) - 0.5);      % 隐含层—输入层权重矩阵
V = 0.3*(rand(3,1) - 0.5);      % 输出层—隐含层权重矩阵

%学习率
alpha1 = 0.06;                  % W权重相关学习率
alpha2 = 0.02;                  % V权重相关学习率

% 动量因子
beta1 = 0.9;
beta2 = 0.9;

% 限幅值
xmax = 1;
ymax = 1;
umax = 1;

% 变量初始化，_1 结尾为上一次的值，_2 结尾表示上两个时刻值
yi3_1 = 0;                      
u_1 = 0;
u_2 = 0;
x_1 = 0;
yi2_1 = 0;
yo3_1 = 0;
dW = 0;                         % W变化量，初始化为0
dV = 0;
yo = zeros(3,1);     

xdata = [];                     % 记录实际值
zdata = [];                     % 记录期望值

N = 2000;                       % 迭代次数
for k=1:N
    %% 设置目标值
    if k <500
        z = 0.5;
    elseif k < 1000
        z = 0.2;
    elseif k < 1500 
        z = 0.8;
    else
        z = 0;
    end
    
    %% 正向传播
    % 隐含层输入
    yi = W * [z; x];                
    
    % 隐含层输出
    yo(1) = yi(1);                  % P  
    yo(2) = yo(2) + yi(2);          % I
    yo(3) = yi(3) - yi3_1;          % D
    yo(find(yo>ymax)) = ymax;       % 限幅
    yo(find(yo<-ymax)) = -ymax;
    
    % 输出层输出
    u = V' * yo;                    
    u(find(u>umax)) = umax;
    u(find(u<-umax)) = -umax;
    
    % 受控对象输出
    x = u^2 + 0.1*x_1;                  
    x(find(x>xmax)) = xmax;         
    x(find(x<-xmax)) = -xmax;
    
    % 损失函数
    J = 1/2 * (z-x)^2;              

    %% 反向传播
    
    % 计算J对V的偏导数 pJ_pV
    pJ_px = -(z-x);
    
    temp = u - u_1;                     % u(k)-u(k-1)
    if temp >= 0 && temp < 1e-20        % 除数限幅在(-∞,-1e-20]∪[1e-20,+∞)
        temp = 1e-20;                   % 防止除数为0
    elseif temp < 0 && temp > -1e-20
        temp = -1e-20;
    end
    px_pu = sign((x - x_1) / temp);     % 使用sign，限幅在{-1,0,1}
    
    pJ_pV = pJ_px * px_pu * yo;         % 误差对V的偏导,yo为3x1向量
    
    % 计算J对W的偏导数 pJ_pW
    pJ_pw11 = pJ_px * px_pu * V(1) * 1 * z;
    pJ_pw12 = pJ_px * px_pu * V(1) * 1 * x;
    
    temp = yi(2) - yi2_1;               % yi(k) - yi(k-1)
    if temp >= 0 && temp < 1e-20        % 除数限幅在(-∞,-1e-20]∪[1e-20,+∞)
        temp = 1e-20;                   % 防止除数为0
    elseif temp < 0 && temp > -1e-20
        temp = -1e-20;
    end
    pJ_pw21 = pJ_px * px_pu * V(2) * (sign(yi(2)/temp)) * z;
    pJ_pw22 = pJ_px * px_pu * V(2) * (sign(yi(2)/temp)) * x;

    temp = yo(3) ;                      % yo(k)
    if temp >= 0 && temp < 1e-20        % 除数限幅在(-∞,-1e-20]∪[1e-20,+∞)
        temp = 1e-20;                   % 防止除数为0
    elseif temp < 0 && temp > -1e-20
        temp = -1e-20;
    end
    pJ_pw31 = pJ_px * px_pu * V(3) * (sign((yo(3) - yo3_1)/temp)) * z;
    pJ_pw32 = pJ_px * px_pu * V(3) * (sign((yo(3) - yo3_1)/temp)) * x;
  
    pJ_pW = [pJ_pw11, pJ_pw12; pJ_pw21, pJ_pw22; pJ_pw31, pJ_pw32]; % 误差对W的偏导
    
    %% 动量法权重更新
    dW = beta1 * dW + (1-beta1) * pJ_pW;        % 动量法计算增量
    dV = beta2 * dV + (1-beta2) * pJ_pV;
    
    W = W - alpha1 * dW;                        % 权重更新
    V = V - alpha2 * dV;
    
    %% 其他数据更新
    u_2 = u_1;
    u_1 = u;                                    % u(k-1)
    x_1 = x;                                    % x(k-1)
    yi2_1 = yi(2);                              % yi2_1 为上一次 yi(2)
    yo3_1 = yo(3);                              % yo3_1 为上一次 yo(3)
    
    xdata = [xdata; x];                         % 记录数据方便绘图
    zdata = [zdata; z];
end

%% 可视化
plot(1:N, zdata, 'LineWidth',1.5);hold on;      % 目标值
plot(1:N, xdata, '.'); hold off                 % 实际值
legend('目标值','实际值')                        % 标注
title('PID神经网络');                           % 标题

运行结果如下

在这里插入图片描述

可见，跟踪效果很不错。不过，和所有神经网络一样，参数要小心谨慎，如果设置不合理，就达不到期望的效果。

四、总结与展望

和PID和神经网络相比，PID神经网络有什么美妙之处呢？其实，它结合了两者的优势，PID使用了反馈控制，好的PID参数使得实际值快速稳定地接近目标值，但是有时候很难调节出一组好的PID参数。神经网络的优势是，通过误差反馈传播调整权值，这个过程是自动完成的。PID参数已经通过权值的方式隐藏在网络结构中，使用反向传播时该参数自动得到优化。

其实，PID神经网络大显身手的地方是多输入多输出的耦合系统的控制，下一次，就用来做一个无人机姿态控制器。

— 完 —

大强强小强强

关注

26
点赞
踩
140

收藏

觉得还不错? 一键收藏
打赏
19
评论
PID神经网络原理与MATLAB实现（SISO）

序最近想实现一下PID神经网络，但是书籍和博客都令人头疼，主要是卡在误差反向传播的计算过程中。找一篇通俗易懂的文章实在不易，最终，只能自己静下心，仔细琢磨。只要每一步在逻辑上都是合理的，我们有理由相信能够得到正确的结果。抱着这样的心态，由浅入深，来实现一下。一、网络结构定义简单起见，假设一个受控系统单输入 uuu 单输出 xxx ，使用一个PID神经网络来作为控制器，使得系统输出达到目标值 zzz。其中，W3×2W_{3\times 2}W3×2和V1×3V_{1\times 3}V1×3.
复制链接

扫一扫