初心者が書いたMatlab強化学習のまとめ

Matlabで強化学習をしてみた
=============

いきなり大学のセミナーで「はい、じゃ君たちの班の課題は強化学習ね！」と放り出され、6週間でなんとか一応結果を出すことができた。しかしかなりきつかった。

Matlabの強化学習はPythonのそれと異なり、Matlab使用者数が少ないため参考にできる事例が少なく（英語でもかなり少ない）、またMatlabの説明が非常に難解なため、読んでも意味不明で試行錯誤する中で意味をつかんでいくことが多々あり、さらにツールボックスがかなり新しいことも相まって、現時点では新たに取り組むハードルが高い。少しでもこれからMatlabで強化学習していく人たちの助けになればと思い、この度共有することにした。


読む上での注意
-----------
- 筆者はqittaとMarkdownとGitHub初心者（そもそもこれが初めての投稿）なので、読みにくい場合は、あなたの人生に渡って超汎用的にトレーニングされた脳内で適宜補正・補完・変換をしてほしい。
- GitHubで共有したプログラムやレポート、ReadMeには英語とドイツ語が入り乱れているので、適度に翻訳をかけて欲しい（直すのが面倒）。
- Matlabもプロではないので「書き方が汚い！」とかの批判は無しで。読みたくなければそっと閉じてね。




##この記事の目標

Matlabを用いた強化学習には、そのやり方において2種類あることを紹介し、それぞれにおいて参考になる事例のリンクの共有し、また私のコードの共有により参考にできる強化学習の事例をそれぞれ一つずつ増やすことが目標である。
さっさと例だけ拾って参考にしたい人は下のリンクからどうぞ
https://github.com/CreamyFoam/RLforQiita/tree/main


##強化学習(RL)とは
[Wikipedia](https://ja.wikipedia.org/wiki/%E5%BC%B7%E5%8C%96%E5%AD%A6%E7%BF%92)曰く
>強化学習（きょうかがくしゅう、英: reinforcement learning）とは、ある環境内におけるエージェントが、現在の状態を観測し、取るべき行動を決定する問題を扱う機械学習の一種。エージェントは行動を選択することで環境から報酬を得る。強化学習は一連の行動を通じて報酬が最も多く得られるような方策（policy）を学習する。環境はマルコフ決定過程として定式化される。代表的な手法としてTD学習やQ学習が知られている。



RLの詳しい説明はこのブログの焦点ではないので、詳しくは各自で調べて欲しい。


RLをより詳しく知りたい人はMatlabのRLの[強化学習入門Onramp](https://jp.mathworks.com/learn/tutorials/reinforcement-learning-onramp.html)をすると良いだろう。前半のRLの概念の学習は非常によくできている[^1]。なお残念ながら日本語版はまだ存在しない。

[^1]: しかし、後半のNeural Network(NN)の作成の仕方の説明は正直必要なかった。むしろ、defaultのAgentをどうやってSimulinkやmファイルに実装していくかを重点的に解説してほしかった。完全に理解している人が作成したため、無駄に細かく説明されて混乱を招いている感じ。最初に理論から入った人には前半は退屈だし、実践から入った人（野口英世がいうならば「パラシュート型学習」）には後半は不要。


##MatlabでのRLの使い方（2種類）
MatlabではRLToolboxは主に2種類の使い方ができる。

###mファイルでクラスを作成する方法
シミュレーションしたい事柄を記述したハンドルクラスを作成し、それを用いて学習を行う。
こちらではSimulinkを用いることなく、全てmファイル上で作成し実行することができる。Simulinkアレルギーのある方にはうってつけの方法である。
参考になるページはこちら
[Game2024](https://qiita.com/HYCE/items/102b2d58218271cc61eb)
[単振り子](https://qiita.com/motorcontrolman/items/958ad922583b414dca96)

__Agentの作成の仕方などはSimulinkを用いたものと基本的に同じであるため、これらのページを一通り読んでから読み進めて欲しい。(後々説明するのが面倒なので...)__

mファイルを用いたシミュレーションの方法についてはこれらのページだけで十分なので、ここでは説明しない。読むだけでは理解できない場合は、試しに簡単な運動を記述するハンドルクラスを単振り子の例と同様に作成して実際に動かしてみることをお勧めする。私も練習がてら小型飛行機であるセスナ172の二次元運動モデル(垂直方向)を作成したため、参考資料として共有する。

[mファイルの機械学習Cessna172](https://github.com/CreamyFoam/RLforQiita/tree/main/rl2DmFileModel)

注!: なお、mファイルの方はSimulinkの方より先に書いたため、AgentのOption等の使い方が正しくない可能性があります。疑わしい場合はご自身でお調べください。



###Simulinkを用いる方法
こちらではSimulinkのシミュレーションモデルを用いてRLを行うことができる。すでに複雑なSimulinkモデル等が手元にある場合、この方法を用いれば新たにモデルを作成し直す必要がない。
それ以外のメリットとしては、isDoneやReward等の関数を視覚的に配置し管理することが可能であることがあげられる。

Simulinkがよくわからない人はまずSimulink Onrampをすると良いだろう。こちらは日本語版がある。
[Simulnk Onramp](https://jp.mathworks.com/learn/tutorials/simulink-onramp.html)

またこの方法の参考になるページとしてMatlab公式の二足歩行ロボットがある。残念ながらこちらには日本語版はない。
[Train Biped Robot to Walk](https://jp.mathworks.com/help/reinforcement-learning/ug/train-biped-robot-to-walk-using-reinforcement-learning-agents.html)

また以下が私のコードである。

[Simulinkを用いた機械学習3D飛行モデル](https://github.com/CreamyFoam/RLforQiita/tree/main/rl3DSimulinkModel)

このコードを開いて比較しながら読み進めて欲しい。


##Simulinkを用いる方法の実装　
以下が私のSimulinkのRLのモデルの概要である。
EasyModel.slxを開いてみてほしい。

![Picture1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1848682/f672f987-8a48-ee5f-5929-62869f3d3b5f.png)




赤い枠で囲まれたブロックがRLに必須の部分であり、緑の枠はアニメーションの表示や座標変換等のそれぞれの課題に依存する付加的なブロックである。

赤い枠のブロックとして以下の要素がある。

- Simulationモデル
- Agent
- Reward
- isDone

またsimulinkには登場しない必須の要素として、Resetfunctionがある。

__Simulationモデル__には機械学習で動かしたいシミュレーションのモデルが入っている。AgentからのAction(指示)を受け、このモデルは指定した時間の長さ(Sample Time)において受け取ったActionに基づいてシミュレーションを実行し、その状態をシミュレーション時間ステップごとに出力する。例えば、車なら、Agentはハンドルの操作角をSample Time時間一定として指定し、Simulationモデルはそれをもとに車のシミュレーションを行う。もし道路が乾いていたら車は機敏に曲がるだろうし、雨の日なら曲がりきれずにスリップしてしまうかもしれない。Sample Timeが0.1秒で車の時間ステップが0.01秒なら、10ステップ分は同じハンドル操作量が保持されるわけである
そしてSimulationモデルはその時の車の状態(位置、速度ベクトル等)を時間ステップごとに出力する。



__Agent__とは先述したとおり、Simulationモデルの情報(Observation)を元にNNで「考え」Actionを与える存在である。Agentにはその入力と出力の種類の組み合わせで大まかに4種類に分けられる。入力出力ともに離散・連続の2種類があり、離散ではスイッチのオン・オフのように離散した値をとることができ、連続では指定した範囲(例えば[-10,10])の値を自由にとることができる。それら2種類ずつの組み合わせで合計4種類のAgentの種類がある。


__Reward__とは、AgentのActionに基づいて実行されたシミュレーションの状態を評価し、スカラー値として返す関数である。AgentはRewardがより多くなるような方向にActionを修正していく。その修正の方法はAgentの種類による。厳密な仕組みについては各自で調べて欲しい。

__isDone__とは、シミュレーションを継続するか判断する関数である。例えば、車が壁にぶつかった場合はそれ以上のシミュレーションは不要なのでisDoneはtrueを返す。

__Resetfunction__とはisDoneによってターミネートされたシミュレーションの状態を次のシミュレーションのために初期化するものである。一番目のシミュレーションやその他追加のパラメータを設定するためにPreloadファイルをメインのプログラムの先頭で読み込むことが多いが、基本的にはその部分のシミュレーションに関わる部分を書き写せばよい。

以上が各要素の大まかな説明である。



###取り組んだ課題
今回取り組んだ課題は、二つのWP(目標物・経由地)を通過しながら、その間に配置された障害物を避けるというものである。（課題の当初は4つのWPと4つの障害物が複雑に配置されたものを行なったが、どうしても達成できなかったため、課題を簡略化した。(課題設定も自らで行うというセミナー課題だったため）)

![Screenshot 2021-08-05 at 23.12.48.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1848682/909bc892-6b80-aca3-e71f-01260514b581.png)

緑色の球体がWPであり、赤い線の先端が表す飛行機はこの球体の内部を通過することが目標である。黒色の球体が障害物であり、これに衝突することは避けなければならない。この図はシミュレーションの初期状態を表しており、飛行機は高度50mから水平に飛行を開始する。

実際には、与えられた3次元飛行モデルが非常に複雑であり、計算に時間を要するので、簡易化されたSimulationモデルでAgentを学習させ、そのAgentを実際の複雑なSimulationモデルに投入するという手段をとった。Agentが学習したのは簡易モデルにおいてであり、また著作権の関係上複雑なSimulationモデルは公開できないため、今回は簡易モデルのみを説明する。


####Simulationモデルの運動方程式
飛行機の位置を$(X,Y,Z)$で表し、その機首方位を$\chi,\gamma$で表し、機速を一定($VK$)とすれば、その単位時間ステップ($TS$)ごとの運動は以下のように表すことができる。(Simulationブロックをひらけば同じ式を見つけられる）



```math
\begin{pmatrix} X(n) \\ Y(n) \\ Z(n) \end{pmatrix} = \begin{pmatrix} X(n-1) \\ Y(n-1) \\ Z(n-1) \end{pmatrix} + VK \cdot \begin{pmatrix} cos(\gamma(n)) cos(\chi(n)) \\ cos(\gamma(n)) sin(\chi(n)) \\ -sin(\gamma(n)) \end{pmatrix} \cdot Ts
```



```math
\begin{pmatrix} \chi(n) \\ \gamma(n) \end{pmatrix} = \begin{pmatrix} \chi(n-1) \\ \gamma(n-1) \end{pmatrix} + \begin{pmatrix} \dot{\chi}(n) \\ \dot{\gamma}(n) \end{pmatrix} \cdot Ts
```

なお、$n$は$n$番目のタイムステップを表している。ここではこの式の導出については詳しくは述べない。
後述するAgentはActionとして$\dot{\chi},\dot{\gamma}$を出力し、Simulationモデルは受けとったActionをここに代入する。
このような、運動方程式等を離散化してシミュレーションモデルを作る方法については、すでに紹介した[単振り子](https://qiita.com/motorcontrolman/items/958ad922583b414dca96)のページを確認してほしい。



####Agentの作成
Simulink上では簡単にAgentブロックをReinforcement Learning Toolboxから見つけ出し配置することができるが、残念ながらそれだけではAgentの設定は終わっていない。このAgentはWorkspaceに保存される必要があり、そのAgentの作成方法はmファイルを用いたAgentの作成方法と同じである。（故に詳しい説明は上記[単振り子](https://qiita.com/motorcontrolman/items/958ad922583b414dca96)と[Game2024](https://qiita.com/HYCE/items/102b2d58218271cc61eb)と[Train Biped Robot to Walk](https://jp.mathworks.com/help/reinforcement-learning/ug/train-biped-robot-to-walk-using-reinforcement-learning-agents.html)を見ていただきたいのだが、それではこの投稿の意味がなくなってしまうので、以下にAgent作成のプログラムとコメントアウトされた解説を載せる。
このコードはrlEasyModel.m内に含まれる。

<details><summary>Agentの作成のプログラム</summary><div>

```matlab
%% Choose what you want to do
NewAgent = false;%create new agent. if false:The later specified agent will be tained further(reuse of agent).
doTraining = false;%if false:run simulation
maxEpisodes = 2500;%Max number of Episodes. (Each simulation is called "Episode". )

Ts = 0.01; % Simulation step time
Tf = 300; % Simulation end time


%% Load Model
%load the Simulink file
mdl = 'EasyModel';%if Error "Unrecognized function or variable 'mdl'".=> try run this code several times
open_system(mdl)

numObs = 10;%Number of observation values.
obsInfo = rlNumericSpec([numObs 1]);
obsInfo.Name = 'observations';

numAct = 2;
actInfo = rlNumericSpec([numAct 1],'LowerLimit',[-0.1;-0.1],...
                                   'UpperLimit',[0.1 ;0.1]);%[chiDot,gammaDot] determined by the maneuverbility
actInfo.Name = 'actions';

blk = [mdl,'/RLAgent'];%install the agent which is made in this code to the agent block in simulink
env = rlSimulinkEnv(mdl,blk,obsInfo,actInfo);%making environment
env.ResetFcn = @(in) ResetFcn(in);%assign reset function which you made

%% Options
%Reuse of same agent
agentOpts = rlDDPGAgentOptions;%to save buffer of training for reuse of pretrained agent
agentOpts.SaveExperienceBufferWithAgent = true;% Reference https://jp.mathworks.com/help/reinforcement-learning/ref/rlddpgagentoptions.html?searchHighlight=SaveExperienceBufferWithAgent&s_tid=srchtitle
agentOpts.ResetExperienceBufferBeforeTraining = false;
agentOpts.NumStepsToLookAhead = 1;
agentOpts.ExperienceBufferLength  = 1e6;
%agentOpts.MiniBatchSize = 256;
agentOpts.DiscountFactor = 0.99;
agentOpts.SampleTime = 0.2;% Agent sample time. This defines the time gap between two actions. Once a set of action is given, the values will be held for Sampletime seconds 


maxSteps = floor(Tf/Ts);
trainOpts = rlTrainingOptions(...
    'MaxEpisodes',maxEpisodes,...
    'MaxStepsPerEpisode',maxSteps,...
    'ScoreAveragingWindowLength',250,...
    'Verbose',false,...
    'Plots','training-progress',...%or 'none'
    'StopTrainingCriteria','EpisodeCount',...
    'StopTrainingValue',maxEpisodes,...
    'SaveAgentCriteria','EpisodeCount',...
    'SaveAgentValue',maxEpisodes);


trainOpts.UseParallel = false;%requires at least 3 CPU cores. 
trainOpts.ParallelizationOptions.Mode = 'async';
trainOpts.ParallelizationOptions.StepsUntilDataIsSent = 32;
trainOpts.ParallelizationOptions.DataToSendFromWorkers = 'Experiences';
trainOpts.SaveAgentDirectory = "savedAgents";
trainOpts.StopOnError = "off";


%% new agent or not.

if NewAgent == true
 %agent = rlDDPGAgent(obsInfo, actInfo, agentOpts);%default agent Agent type "DDPG"
 agent = createDDPGAgent(numObs,obsInfo,numAct,actInfo,0.3); %use of agent setting of RLWalker which is an example from Mathworks. Reference https://jp.mathworks.com/help/reinforcement-learning/ug/train-biped-robot-to-walk-using-reinforcement-learning-agents.html
 %With the use of agent setting for Walker, the programs
 %"createDDPGAAgent.m" and "createNetworks.m" will be used.
else
 load('savedAgents/WPObstWP/Agent4000.mat','saved_agent'); %use of an already trained agent    %Reference https://jp.mathworks.com/matlabcentral/answers/495436-how-to-train-further-a-previously-trained-agent
 agent = saved_agent;
end

```

</div></details>



####Reward
Rewardとして、CloserBonusとWPBonus、CrashPenaltyを用意した。(Rewardブロックを開いてみてほしい)

CloserBonusは飛行機がWP(目標地点)に近づくほど値が大きくなる関数で、以下のような分布をしている。この関数は微分可能であり、Agentが方向を定めることを手助けする。

![Picture1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1848682/2cc34f66-1790-183f-52ba-2feade506d4d.png)


WPBonusはWPに到達した時に与えられる報酬である。CloserBonusと比べ非常に大きく設定することで、あたりを長く飛び回って小さな報酬をかき集めるより、WPに到達することのほうが良いとAgentに学習させることができる。

CrashPenaltyは飛行機が障害物に衝突した場合にマイナスの報酬を与える。これによって飛行機は障害物の回避運動を学習するようになる。

####isDone
飛行機が墜落した場合や次のWPへのコースから大きく外れた場合、障害物に衝突した場合にisDoneはtrueを出力する。それによって、そのトレーニングエピソードは終了し、次の新たなエピソードが開始される。(isDoneブロックを開いてみてほしい)

####Resetfunction
Resetfunctionを書いた当時、考えるの放棄していた筆者はPreloadの関数の中身を完全にコピペしていた。flightplanとobstacleはPreloadに含まれるので本来はResetfunctionには必要ないものである。(ResetFcn.mを開いてみてほしい)

<details><summary>Resetfunction</summary><div>

```matlab
function in = ResetFcn(in)
flightplan = [2      0,   0;... 
              150,   0, -60;...
              350,   0, -40;...
              zeros(18,3)];%本来必要ない
 
obstacles =  [1,     0,   0,  0;...
              0,     0,   0,  0;
              250,   0, -50,  30;...
              zeros(19,4)]; %本来必要ない
              

%%以下の部分がシミュレーションを初期化する初期値である。
r0 = [0; 0; -50];
Dir0 = [0;0];
VK = 25;
end
```

</div></details>


###すべての要素が用意できたら
すべての要素が用意できたら、(必要に応じてPreload関数を先頭に挿入しつつ）Agentの作成のコードの下に以下のコードを追加する。

```matlab
%% Training? Simulation?

if doTraining    
    % Train the agent.
    trainingStats = train(agent,env,trainOpts);
else
    % Simulation
    % Load a pretrained agent for the selected agent type.
    load('savedAgents/WPObstWP/Agent4000.mat','agent');
    rng(0)%if values are randamised. This specifies the seed of "rand".

    simOptions = rlSimulationOptions('MaxSteps',maxSteps);
    experience = sim(env,agent,simOptions);

end
```
最初に指定したdoTrainingがtrueのときはトレーニングが開始され、それ以外ではすでに学習済みのAgentを用いてシミュレーションを行う。


##RLの結果
適切なHyper parametersを見つけたり、Rewardを調整したりと試行錯誤を(前の複雑なシミュレーションでの学習を含め)2週間くらい続け、最終的に以下の結果を得ることができた。

3000回のエピソード(シミュレーション回数・試行回数)により、定常的に一つ目のWPに到達できるようになり、また障害物をそれなりの確率で避けられるようになった。
![Picture1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1848682/1aede5ad-2e62-c73a-a2e0-558dfdae869a.png)
3000回のエピソードのAgentでのシミュレーションがこちら。
![Picture1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1848682/83f21a9c-639c-1db6-4466-c694097892e0.png)
一つ目の緑色のWPを通過し、黒色の障害物を避けているが、二つ目のWPの場所がわからないようで、あらぬ方向に飛んでいってしまっている。


これまでで一番良い結果だったので、同じAgentを用いてさらに4000回のトレーニングを行った。（若干のRewardの変更も行った)
その結果、トレーニング後半において、大きな報酬の増加があり、最後の100エピソードほどで定常的にWP2に到達できるようになった。はっきりとした報酬の収束はみられないが、これはWP1では通過した瞬間に次のWPであるWP2を行き先として示され一回しか1000ポイントを得るチャンスがないのに対し、WP２ではそれが最終目的地であるために、WP2内を飛行している限り1000ポイントを加算し続けることができ、WP2の緑の球体の通過の仕方で1000ポイント単位で報酬が変動してしまうためである。
![Picture2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1848682/ae79828c-372b-b610-0850-ae8407a69276.png)
そのAgentを用いた結果がこちら。
![Picture1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1848682/c26d9153-81f8-6115-0e62-d018dc0e2ae3.png)
![Picture1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1848682/ec5ab198-b076-5179-2030-0291a8ce818d.png)

AgentはWP1(上記画像において奥の緑）はかするだけで良く、黒い障害物はギリギリでかわし、WP2はできるだけ内部にとどまれるように、その中心を目指して飛行していることがわかる。

###その後
先述したとおり、このAgentは後に非常に複雑なシミュレーションモデルに投入された。
![Screenshot 2021-08-06 at 08.07.25.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1848682/6a66a4e6-0987-0a16-92d3-f02aa3b3f2ba.png)
残念ながら、最初のWP1(画像右の球体)を通過しなかったため、WP2に到達する前にシミュレーションがisDoneによって終了させられている。

Agentの応用がうまくいかなかった原因として2点考えられる。

- 簡易モデルでは速度が一定だったため、速度変化のある複雑なモデルではAgentはどう反応して良いかわからなかったこと。
- 簡易モデルでは角速度をActionとしてSimulationモデルに与えたため、タイムラグなしに機体を操作した。しかし、複雑なモデルでは目標角をSimulationモデルは受け取り、それに対して内部のモデルが適切に方向舵を調整し、機体がそれに反応するというダイナミクスが存在し、目標角の変化に対する実際の機体の機動は伝達関数で表される。そのボード線図をみると、AgentのActionのとりうる周波数領域では位相の遅延とゲインの減少があり、Agentが予想したような機体の機動にならなかったことが推測される。(AgentのActionである角速度からSimulationモデルのインプットである目標角へは特別な関数を用意し変換した。その関数の出来がよくなかったということである。)

##追加の考察
実は、Simulinkを用いたものとmファイルを用いたものではSimulationモデルに与えるActionが大きく異なる。Simulinkのほうでは目標角のみを与え、機体の飛行の制御は内部の飛行モデルに任せていた。つまり、飛行機は飛ぶ方法はわかった上で、その飛行経路のみを勉強すればよかった。一方でmファイルの方ではスロットルや方向舵といった、人間がパイロットとして操縦するような値をActionとして与えた。そのため、Agentはまず機体を安定して飛ばす方法を学習する必要がある。目標角のような入力はhigh levelな制御であり、一方機体を直接制御することはlow levelな制御である。low levelな制御はこれまでの制御理論のほうが速いうえに確実であり（安定性解析ができる）。機械学習はそれよりもhigh levelに位置する制御に用いるほうが良い。実際mファイルでは機体は飛行を学習することができなかった。

また、Simulinkの方ではWPと障害物の位置と大きさをエピソードごとにランダムに変化させる学習をしたかった。今回の学習では課題が常に一定であったため、Agentはこの環境に特化した学習をしていると思われる。ランダム性を入れることで、WPと障害物の位置をObservationから理解し、それぞれに到達する/避けるという「本質」を学習させることができるだろう。(それにはよりよいコンピュータが必要でもある。）

##よくわからなかった人へ
残念ながら、書いている私も書けば書くほどMatlabのDescriptionが(難解ではあるものの)重要だと思うようになった。それぞれのコマンドの役割などについて、初見では全く理解できないこともあると思う（私もWalkerの例を見た時はそうだった)。その場合はわからない部分をMatlabのヘルプページで調べてみて、とりあえずその説明を読んだ後にそのパラメータを変えてみるなどして、その部分の役割などを理解していけば良いと思う。

##終わりに
記事を書き始めるまでは、この記事を読めばとりあえず動くSimulinkを用いたRLのモデルが作れコードが書けるような内容にしようと意気込んでいたが、書き進めるうちに何をどこまで書くべきか悩み、結局自らのコードを共有するだけの記事となってしまった。非常に参考にさせていただいた「単振り子」と「Game2024」方々のような素晴らしい記事ではないが、Matlab強化学習の一つの参考にできる事例として活用してくだされば幸いである。(再記　mファイルの方はSimulinkの方より先に書いたため、AgentのOptionの使い方などが正しくない可能性があります。疑わしい場合はご自身でお調べください。)私の稚拙な文章にも関わらず最後までお読みくださり、ありがとうございました。