简单RNN的MATLAB实现

近期学习了RNN（循环神经网络），但网上的代码大多都是python编写的，且基本都是用工具箱。于是自己结合网上的代码，用MATLAB进行了编写，大致框架如下，但可能存在问题，希望与读者多交流，后面的激活函数可以选择sigmid/tanh/Ruel. 

% implementation of RNN 
% 以自己编写的函数为例进行计算
clc
clear
close all
%% training dataset generation
dx = 0.15;
x1 = 0:dx:1;   % 输入值
x2 = x1.^2;
y_start = x1./(x2+1);   % 真实输出值
the_dim = length(x1);  % 输入时间序列的长度
%% input variables
alpha = 0.1;       % 学习率
input_dim = 2;     % 输入变量维数
hidden_dim = 16;   % 隐含层大小
output_dim = 1;    % 输出变量维数

%% initialize neural network weights
synapse_0 = 2*rand(input_dim,hidden_dim) - 1;    % 输入层与隐含层权值矩阵：input_dim*hidden_dim
synapse_1 = 2*rand(hidden_dim,output_dim) - 1;   % 隐含层与输出层权值矩阵：hidden_dim*output_dim
synapse_h = 2*rand(hidden_dim,hidden_dim) - 1;   % 前一时刻的隐含层与现在时刻的隐含层全值矩阵：hidden_dim*hidden_dim

  % 储存权值更新
synapse_0_update = zeros(size(synapse_0));   
synapse_1_update = zeros(size(synapse_1));
synapse_h_update = zeros(size(synapse_h));

%% train logic
for j = 1:1000   % 迭代训练
    
    y_out = zeros(size(y_start));    % 初始化空的二进制数组，储存神经网络的预测值
    
    overallError = 0;   % 重置误差值
    
    % layer_2为输出层输出； layer_1为隐含层输出； layer_1x为隐含层输入
    layer_2_deltas = [];      % 记录layer_2的导数值
    layer_1_values = [];      % 记录layer_1的值
    layer_1x_values = [];      % 记录layer_1x的值
    layer_1_values = [layer_1_values; zeros(1, hidden_dim)]; % 储存隐含层的输出值，并提前设置上一时刻的值
    
    % 开始对一个序列进行处理
    for position1 = 1:the_dim       % 输入的时间序列长度
        X = [x1(position1) x2(position1)];   % X 是 input
        y = [y_start(position1)]';           % Y 是label，用来计算最后误差
        
        % 这里是RNN，因此隐含层比较简单
        % X ------------------------> input
        % synapse_0 ----------------> U_i  输入层与隐含层间的权值
        % layer_1_values(end, :) ---> previous hidden layer （S(t-1)）  
        % synapse_h ----------------> W_i
        % layer_1 ------------------> new hidden layer (S(t))
        layer_1x = X*synapse_0 + layer_1_values(end, :)*synapse_h;      % 隐含层输入
        layer_1 = acf(X*synapse_0 + layer_1_values(end, :)*synapse_h);  % 隐含层输出
        
        % layer_1 ------------------> hidden layer (S(t))
        % layer_2 ------------------> 最终的输出结果，其维度应该与 label (Y) 的维度是一致的
        % 这里的 sigmoid 其实就是一个变换，将 hidden layer (size: 1 x 16) 变换为 1 x 1
        % 有写时候，如果输入与输出不匹配的话，使可以使用 softmax 进行变化的
        % output layer (new binary representation)
        layer_2x = layer_1*synapse_1;  % 输出层输入
        layer_2 = acf(layer_1*synapse_1); % 输出层输出
        
        % 计算误差，根据误差进行反向传播
        % layer_2_error ------------> 此次（第 position+1 次的误差）
        % y 是真实结果
        % layer_2 是输出结果
        % layer_2_deltas 输出层的变化结果，使用了反向传播，见那个求导（输出层的输入是 layer_2，那就对输入求导即可，然后乘以误差就可以得到输出的diff）
        layer_2_error = y - layer_2;                        % 误差
        layer_2_deltas = [layer_2_deltas; layer_2_error*acf_output_to_derivative(layer_2x)];  
        
        % 总体的误差（误差有正有负，用绝对值）
        overallError = overallError + abs(layer_2_error); % 累计绝对误差
        
        % 记录神经网络的输出
        y_out(position1) = layer_2;
        
        % 记录下此次的隐含层输出 (h(t))，以便在反向传播时使用
        layer_1_values = [layer_1_values; layer_1];     % 因为多设置了零时刻的隐含层，所以长度比样本时间序列多一
        layer_1x_values = [layer_1x_values; layer_1x];
    end
    
    % 计算隐含层的diff，用于求参数的变化，并用来更新参数，还是每一个timestep来进行计算
    future_layer_1_delta = zeros(1, hidden_dim);  % 提前设置储存量
    
    % 重置权值参数变化量
    synapse_0_update = synapse_0_update * 0;
    synapse_1_update = synapse_1_update * 0;
    synapse_h_update = synapse_h_update * 0;
    % 开始进行反向传播，计算 hidden_layer 的diff，以及参数的 diff
    for position2 = the_dim:-1:1   % 误差反向传播，同时也是按照时间序列从尾传到头
        % 因为是通过输入得到隐含层，因此这里还是需要用到输入的
        % x -> (operation) -> y, x_diff = derivative(x) * y_diff
        % 注意这里从最后开始往前推
        X = [x1(position2) x2(position2)];
        % layer_1 -----------------> 表示隐含层输出 hidden_layer (h(t))
        % prev_layer_1 ------------> (h(t-1))
        layer_1 = layer_1_values(position2+1, :);     % 提取当前时刻隐含层输出
        prev_layer_1 = layer_1_values(position2, :);  % 提取前一时刻隐含层输出
        layer_1x =layer_1x_values(position2, :);
        % layer_2_delta -----------> 就是隐含层的diff
        % hidden_layer_diff,根据这个可以推算输入的diff以及上一个隐含层的diff
        % error at output layer
        layer_2_delta = layer_2_deltas(position2, :);
        % 这个地方的 hidden_layer 来自两个方面，因为 hidden_layer -> next timestep, hidden_layer -> output，
        % 因此其反向传播也是两方面
        % error at hidden layer
        layer_1_delta = (future_layer_1_delta*(synapse_h') + layer_2_delta*(synapse_1')) ...
                        .* acf_output_to_derivative(layer_1x);
        
        % let's update all our weights so we can try again
        synapse_1_update = synapse_1_update + (layer_1')*(layer_2_delta);
        synapse_h_update = synapse_h_update + (prev_layer_1')*(layer_1_delta);
        synapse_0_update = synapse_0_update + (X')*(layer_1_delta);
        
        future_layer_1_delta = layer_1_delta;
    end
    
    synapse_0 = synapse_0 + synapse_0_update * alpha;
    synapse_1 = synapse_1 + synapse_1_update * alpha;
    synapse_h = synapse_h + synapse_h_update * alpha;
    
%     if(mod(j,1000) == 0)
%         err = sprintf('Error:%s\n', num2str(overallError)); fprintf(err);
%         d = bin2dec(num2str(d));
%         pred = sprintf('Pred:%s\n',dec2bin(d,8)); fprintf(pred);
%         Tru = sprintf('True:%s\n', num2str(c)); fprintf(Tru);
%         out = 0;
%         size(c)
%         sep = sprintf('-------------\n'); fprintf(sep);
%     end
    
end
plot(y_out,'--')
hold on
plot(y_start)

%%% 激活函数
% % sigmid
function yy = acf(xx) 
  yy = 1./(1 + exp(-xx)); 
end

function yy = acf_output_to_derivative(xx) 
  yy = acf(xx).*(1-acf(xx)); 
end

% % tanh
% function yy = acf(xx) 
%    yy = (exp(xx)-exp(-xx))./(exp(xx)+exp(-xx));
% end
% 
% function yy = acf_output_to_derivative(xx) 
%     yy = 1-acf(xx).^2;
% end

% % Ruel
% function yy = acf(xx) 
%    yy = max(0,xx);
% end
% 
% function yy = acf_output_to_derivative(xx) 
%     yy = xx./abs(xx);
%     yy = max(0,yy);
% end