[MATLAB]ニューラル ネットワーク回帰モデルの学習(dlnetwork)

ニューラルネットワークを用いた回帰曲線の学習をMATLABで行いたい．
コマンドfitrnetの利用するためのサンプルは次の通り．

この方法ではネットワークのカスタマイズができない．
オブジェクトdlnetworkを利用して，学習も自作する必要がある．

sample_dlnet_dim1.m

clc
close all
clearvars
rng("default");

%% 入力層の次元
dim_in=1;

%% 出力層の次元
dim_out=1;

%% 各層の大きさ
layer_sizes=[dim_in 10 20 10 dim_out];

%% 学習データの個数
num_train=500;

%% ミニバッチ環境：auto / cpu / gpu
mbq_env='auto';

%% ミニバッチサイズ：
% 小規模モデル（MNISTなど）：32, 64, 128
% 中規模CNN（ResNet-18など）：64, 128, 256
% 大規模モデル（ResNet-50、ViTなど）：256, 512, 1024（GPUに応じて）
% NLP（BERTなど）：8, 16, 32（大きい入力のため小さめ）
batch_size=64;

%% 最大エポック数100〜500など
epoch_max=500;

%% ウィンドウ幅（moving window size）
% 1〜2エポック：非常に敏感（ノイズにも反応しやすい）→小規模データや早く終わらせたい実験に使う
% 3〜5エポック：標準的、過剰反応しにくく安定した判断ができる→多くの中〜大規模タスクで推奨される
% 6〜10エポック以上：安定性はあるが、停滞の検出が遅れることも→ロスが大きく揺れるタスクや時間に余裕があるときに使う
win=10;

%% 学習の終了条件：loss<=loss_border
loss_border=1e-6;

%% 学習の終了条件：1-loss_previous/loss_present<=loss_ratio_border（0.1％ぐらい）
loss_ratio_border=1e-3;

%% 学習率の初期値
learn_rate_init=0.01;

%% Warm-UP
% 小規模モデル：1〜3エポック
% 中〜大規模モデル（ResNet, Transformer等）：5〜10エポック
% BERT・ViTなどの大規模事前学習モデル：10エポック
learn_warmup=5;

%% 学習データの生成
x_train=4*rand(num_train,dim_in)-1;
y_train=func_dim1(x_train);
fig1=figure('Position',[0 0 1000 700]);
func_fig_data(fig1,x_train,y_train,'学習データ：y=f(x)');

%% ネットワークの定義
layers=[ featureInputLayer(layer_sizes(1),'Normalization','zscore','Name','input') ];
for i=2:(length(layer_sizes)-1)
    layers(end+1)=fullyConnectedLayer(layer_sizes(i),'Name',"fc_"+(i-1));
    layers(end+1)=reluLayer('Name',"relu_"+(i-1));
end
layers(end+1)=fullyConnectedLayer(layer_sizes(end),'Name',"fc_out");
net=dlnetwork(layers);
net=initialize(net);
fig_net=figure('Position',[1200 0 1000 700]);
func_fig_net(fig_net,net,'ネットワークの定義');

%% 学習データ用のミニバッチキュー
ds_x_train=arrayDatastore(x_train,'IterationDimension',1);
ds_y_train=arrayDatastore(y_train,'IterationDimension',1);
ds_train=combine(ds_x_train,ds_y_train);
mbq_train=minibatchqueue(ds_train, ...
    'MiniBatchSize',batch_size, ...
    'MiniBatchFormat',{'BC','BC'}, ...
    'OutputAsDlarray',true, ...
    'OutputEnvironment',mbq_env);

%% Adamで学習
tic;
disp('training...')
epoch=1;
done=false;
loss_ratio=nan;
loss_train=zeros(epoch_max,1);
figLoss=figure('Position',[0 800 1000 700]);
func_fig_loss(figLoss,epoch_max,loss_train);
grad_ave_train=[];
grad_sqrt_ave_train=[];
while ~done && epoch<=epoch_max
    %% 学習率スケジューラ：エポックごとに減衰させる
    % Cosine Annealing：cosineで急現象
    % Warmup：最初ゆっくり上げてから減衰．⼤規模モデルで主流（BERT, ViT で常⽤）
    if epoch<learn_warmup
        learn_rate=learn_rate_init*epoch/learn_warmup;
    else
        learn_rate=learn_rate_init*0.5*(1+cos(pi*(epoch-learn_warmup)/(epoch_max-learn_warmup)));
    end
    %% 学習データの頭出しとシャッフル
    reset(mbq_train);
    shuffle(mbq_train);
    %% 学習データを全て用いて１エポックの計算
    itr=0;
    loss_train(epoch)=0;
    while hasdata(mbq_train)
        [batch_x,batch_y]=next(mbq_train);
        [loss,grad]=dlfeval(@loss_mse_grad,net,batch_x,batch_y);
        [net,grad_ave_train,grad_sqrt_ave_train]=adamupdate(net,grad,grad_ave_train,grad_sqrt_ave_train,epoch,learn_rate);
        loss_train(epoch)=loss_train(epoch)+double(gather(extractdata(loss)));
        itr=itr+1;
    end
    %% １反復当たりの損失
    loss_train(epoch)=loss_train(epoch)/itr;
    %% 増加したら終了
    if epoch-2*win-1>=1
        loss_ratio=sum(loss_train((epoch-win-1):epoch))/sum(loss_train((epoch-2*win-1):(epoch-win)));
        if ~isnan(loss_ratio) && 1-loss_ratio<=loss_ratio_border
            done=true;
        end
    end
    %% 訓練損失が閾値より十分小さくなったら終了
    if loss_train(epoch)<loss_border
        done=true;
    end
    %% 損失をグラフに表示
    if mod(epoch-1,10)==0
        func_fig_loss(figLoss,epoch_max,loss_train);
    end
    %% 表示
    if done || epoch==1 || epoch==epoch_max || mod(epoch,10)==0
        fprintf('[%04d/%04d] learn_rate=%.03e, loss_train=%.1e%s%.1e, 1-loss_ratio=%.1e%s%.1e\n', ...
            epoch,epoch_max,learn_rate, ...
            loss_train(epoch),tfv(loss_train(epoch)<loss_border,"<",">"),loss_border,...
            1-loss_ratio,tfv(1-loss_ratio<=loss_ratio_border,"<",">"),loss_ratio_border);
    end
    %% 次のエポック
    epoch=epoch+1;
end
epoch=epoch-1;
loss_train=loss_train(1:epoch);
t=toc;

%% テストデータを予測で生成
x_test=(-1:0.01:3)';
y_pred=predict(net,x_test);

%% データの表示
fig3=figure('Position',[1200 800 1000 700]);
func_fig_data2(fig3,x_train,y_train,x_test,y_pred,'学習データと予測データの比較','学習データ','予測データ');

%% 真値との誤差
y_true=func_dim1(x_test);
error=mean((y_pred-y_true).^2);
fprintf('epochs=%d, time=%.1f, loss_train=%.1e, error=%.1e\n',epoch,t,loss_train(epoch),error);

%%========ローカル関数===========

%% 学習データの生成用
function y=func_dim1(x)
y=x.*(x-1).*(x-2);
end

%% モデル損失関数：平均２乗誤差MSE
function [loss,grad]=loss_mse_grad(net,x,y)
y_pred=forward(net,x);
loss=mean((y_pred-y).^2,'all');
grad=dlgradient(loss,net.Learnables);
end

function func_fig_data(fig,x,y,name)
figure(fig);
plot(x,y,'o');
grid on;
xlabel('x');
ylabel('y');
title(name);
drawnow;
end

function func_fig_net(fig,net,name)
figure(fig);
plot(net);
title(name);
drawnow;
end

function func_fig_data2(fig,x1,y1,x2,y2,name,l1,l2)
figure(fig);
plot(x1,y1,'o',x2,y2,'.');
grid on;
xlabel('x');
ylabel('y');
legend(l1,l2);
title(name);
drawnow;
end

function func_fig_loss(fig,epoch_max,loss_train)
figure(fig);
semilogy(1:epoch_max,loss_train,'-');
grid on;
xlabel('Epoch');
ylabel('Loss');
title('学習中の損失MSEの推移');
drawnow;
end