ベイジアン・フィルタ

1.What's new

1.1 Lagrangeの運動方程式とは

1.一般化座標

一般化座標「$q_1,q_2,...q_n$」で、デカルト座標を「$x_1,x_2,...x_n$」を以下のように記述する。

$\Large \displaystyle x_i=x_i(q_1,q_2,...q_n)$

上式を時間「$t$」で微分すると、以下のようになる。

$\Large \displaystyle \dot{x_i}=\sum_{j=1}^n \frac{\partial x_i}{\partial q_j}\dot{q_j} \ \ \ \ (i=1,2,..,n)$

2.共役な一般化運動量

一般化座標「$q_1,q_2,...q_n$」を記述すつ場合に、これらのそれぞれに共役な一般化運動量「$p_1,p_2,...p_n$」を以下のように定義する。

$\Large \displaystyle p_i=\frac{\partial T}{\partial \dot{q_i}}$

3.Lagrangeの運動方程式

以下のような量Lを定義して「ラグランジアン」と呼ぶ。

$\Large \displaystyle L=T-U$

以下のような量Lを定義して「ラグランジアン」と呼ぶ。

$\Large \displaystyle p_i=\frac{\partial L}{\partial \dot{q_i}}$

$\Large \displaystyle F_i=\frac{\partial L}{\partial q_i}$

運動量の時間微分が、力だから、以下の「Lagrange」の運動方程式を得る。

$\Large \displaystyle \frac{d}{dt}\left (\frac{\partial L}{\partial \dot{q_i}}\right)=\frac{\partial L}{\partial q_i}$

1.2 Hamiltonianとは

2.確率ロボティクス

2.1 ベイジアン・フィルタとは

1.確率システム制御

一部しかセンシングできない信号
モデル化できない偶然事象
制御系の変化

→

確率論の導入で制御系の合理的な駆動を可能に！

「ロボットの知覚と行動に関する不確実性」に特に注目する比較的新しいロボット工学の取り組みの一つ．
確率論を用いて「不確実性」を陽に表現する．
ある推定に対して，結果を「最も妥当な一つの推定値」で表すのではなく確率分布として推定情報が表される．
→「あいまいさ」や「信頼度」を定量的に表すことができる．
「不確実性」を考慮することによって，行動を慎重にしたり，不確実性を減らすようにロボットを行動させることができる．
ロボットが不確実さに直面した時に「きわどい行動」を回避させることができる．

2.2 確率ロボティクスとは

1.確率ロボティクスの三本柱

マルコフ過程
ベイズ推定
最適制御（価値反復）

2.マルコフ過程

確率ロボティクスでは、状態遷移が確率的に扱われる。
ある状態「$x$」から行動「$u$」を取った時、次の状態が一意に定まらない。
今、どの状態にあるかも既知ではなく、どの状態にあるかは確率密度関数「$bel(x)$」で表現される。
「$bel$」は、「belief(信念)」を意味し、信念の持ち主は、ロボット自身である。
時刻「$t-1$」から行動「$u_t$」を取って、時刻「$t$」になった時、どの状態にあるかは、以下の式で計算される。

$\Large \displaystyle bel(x_t)=\eta\int_Xp(x_t|u_t,x_{t-1}bel(x_{t-1}))dx_{x-1}$

「$\eta$」は「$\int_Xbel(x_t)=1$」にするためのパラメータである。
多くの場合、ロボットが状態遷移（移動）をすると、「$bel(x)$」は曖昧になって行く。
ところで、「x」は単なる記号であり、ユークリッド空間等、距離がてきる空間に存在している訳ではない。
そのため、上の式は、集合「$X$」内の状態「$x$」同士を行動「$u$」によって、結び付ける役割も果たす。
実際に問題を解く時には、Xに位相や座標を入れたり、離散化したりと、より具体的なものになる。
大事な事は、現代制御や有限なマルコフ連鎖よりも抽象度の高いものから式は出発する事である。
また、モンテカルロ法等、ノンパラメトリックな計算手法にも都合が良い表現になっている。

3.ベイズ推定

ロボットの得た観測「$y$」と状態「$x$」を結びつけるため、ベイズの公式が導入される。
ベイズの公式は、正規化のための変数「$\eta$」を使うと、以下のように表現できる。

$\Large \displaystyle bel(x|y)=\eta p(y|x)bel(x)$

この式は、状態についての知識が「$bel(x)$」である時、新たな観測「$y$」があったとする。
この時、もし「$x$」の場合に観測「$y$」が得られる度合い（尤度）を知っていたならば「$bel(x)$」を「$bel(x|y)$」に更新可能であると言う意味である。
情報の独立性等の条件の下、「$bel(x|y)$」で「$bel(x)$」を上書きしてよい。
状態遷移の場合とは逆に、観測が「$bel(x)$」の分布を縮小する。
注意が必要なのは、尤度関数「$p(y|x)$」の扱いである。
ベイズ統計学では、「$p(y|x)$」の定義を「$y$」の観測者（ロボットやその作成者）に委ねている。
これを理解していないと、「確率が間違っている」など不毛な議論に時間を費やす事になる。

4.最適制御（価値反復）

https://www.jstage.jst.go.jp/article/jrsj/36/2/36_36_140/_pdfhttps://www.jstage.jst.go.jp/article/jrsj/36/2/36_36_140/_pdf

3個目の式は、ロボットにミッションを授けるためのものである。
以下の式は、何か行動や制御に目標があり、その目標に対して状態「x」に、どれだけの価値があるかを定義する式である。

$\Large \displaystyle V(x_t)=\gamma \max_{u_t}[r(x_{t-1},u_t) + \int V(x_{t-1})p(x^{'}|u_t, x_{t-1})dx^{'}]$

ここでは上式を「価値反復」の式と呼ぶ。
「価値反復」とは、上式を用いて繰り返し演算で「$V(x)$」を求めるアルゴリズムの一つである。
右式の計算結果を左辺に代入する事を繰り返すと、「$V(x)$」を求める事が出来る。
分かりやすい例として、移動ロボットが目的地（正確には目的の状態）に最短時間で行くと言う問題を考える。
この場合、「$\gamma=1, r(x, u)=-1$☓(uの実行に必要な時間), $V(x)=0$($\forall x\in$目的の状態)」として価値反復で「$V(x)$」を求めると、「$V(x)$」は、状態「$x$」から目的地までの時間に「-1」を掛けた値となる。
求めた「$V(x)$」で、左辺の価値が最大となる行動「$u$」を求めると、ロボットが行動すべき行動となる。

3.ベルマン方程式

「強化学習に出てくるベルマン方程式を理解しよう」から引用。

3.1 ベルマン方程式の概要

細かい話をする前に、ベルマン方程式がどういうものであるのか、そのゴールを示しておきましょう。
ベルマン方程式とは状態「x」を外部からの入力「u」で制御できると考えているときに、ある種の評価「J」の下で「u」を色々変えてみて、いざ評価「J」を最も良くするような「u」が見つかったときに成り立っているべき方程式を表します。
ここで重要なポイントが2つあります。
まず、我々は勝手に状態「x」なるものが「u」で制御できるのだと考えている点です。
例えば、ある回路に流れる電流を適当な印加電圧によって自在に制御できると、勝手に考えているのです。
ではいざその回路の構成を見てみると、実は配線がちゃんとできていなくて、印加電圧をどんな風に変えてみても電流を操作することができなかったとしましょう。
そうしたならば当然、制御できるわけもないのですから、ベルマン方程式だのなんだの言っても意味がありません。

次に、ベルマン方程式はある評価が最も良くなったときに成り立っている方程式であるという点です。
私達は状態「x」を制御するような外部からの入力「u」を求めたいと考えているわけですが、当然どんな「u」が最適であるのかを知りません。
しかし、「もしも最適な「u」が見つかった暁には、その「u」はベルマン方程式なる式を満たしているはずである」ということだけを知っているのです。
つまり、「u」を色々変えてみて最適なものを見つけようとしているその過程において（すなわち今保持している最適でない「u」という制御則について）、ベルマン方程式が成り立っているということではないということです。

これは言い換えれば、今持っている「u」がベルマン方程式を満たしているか、それが最適なものであるのかそうでないのかを確認できるということです。
そうして、いざベルマン方程式を満たしていない「u」を使ってしまっているのだとすれば、どのようにしてそれを修正して、最適なものを見つけていけばよいのかの手段があれば良いということになります。
これが強化学習におけるベルマン方程式の役割です。

3.2 最適制御と評価関数

1.最適制御

ここで一旦寄り道をさせてください。
状態「x」を外部からの入力「u」によって思い通りに操作してやりたい、という問題を制御問題と呼び、工学システムが複雑化していった産業革命あたりで大きく発展してきました。
その中でも、「x」を制御できるような「u」が色々考えられる場合に、その中で最も良さそうな「u」を選定する制御問題を、最適制御問題と呼びます。

ここに制御問題を簡単に見ておきます。
状態「x」に外部入力「u」を与えたときの時間発展が下記の微分方程式に従うような場合に、

$\Large \displaystyle \dot{x}=Ax+Bu$

この微分方程式の解、すなわち「x」の時間波形を「u」を色々変えてみて自在に操りたいというのが制御問題です。
もしも入力「u」を外部から与えなかった場合には、状態「x」は微分方程式「$\dot{x}=Ax$」という単純な式に従って勝手に時間発展していきます（これは、行列指数関数を使えば簡単に一般解が求まります）。
放っておいたら勝手な振る舞いをする「x」を「u」を使ってなんとかしたいということです。
このときの「なんとかしたい」というのは、「x」の収束先を決めてやりたいというだけの場合もありますし、各時刻でとある値を取ってほしいという場合もあります。
前者ならば、指定した値にどれだけ精確にどれだけ速く収束するかが評価対象になりますし、後者の場合は各時刻全体に渡ってその値を評価する必要があるかもしれません。

2.評価関数

ここで、自分が持ち合わせている「u」の良さを評価するために評価関数「J」なるものを考えます。
例えば、典型的には以下のような評価関数が使われます。

$\Large \displaystyle J=\int_0^\infty((x_{ref}−x)^TQ(x_{ref}−x)+u^TRu)dt$

ここで「$x_{ref}$」というのは「x」が各時刻にどういう値を取っておいてほしいかを表す目標値です。
このような評価関数では、状態「x」と自分が定めた目標値との二乗誤差（みたいなもの）を第一項が表現しており、入力「u」の値の大きさを第二項が表しており、式全体では、それらの時間全体での総和を表しています。

その心は「各時刻で状態は目標値に近い値を取っているべきであり、かつ、そのような状態を実現するに当たって外部からの入力はできるだけ小さな値の方が良い」というものになります。

なぜこんな評価関数が典型となっているのかというと、多くの工学的な応用においては「u」という外部入力は、電流であったり電圧であったり、外部から作用される力であったり、何かしらのエネルギーを要するものであるケースが多いからです。
ですから「u」は小さいほうが良いということが組み込まれます。

また単純に二乗ではなく、「Q,R」を間に挟んだ二次形式になっているのは、状態ベクトルや入力ベクトルのどの成分に重きを置くのか（あるいはその相関に重きをおくのか）を表現するためです。

もちろんこの評価関数は自分が達成したい制御に合わせて設計されるべきであり、この形式に限る必要はありません。
入力や状態を受け取り、スカラー値で良し悪しを評価するような汎関数の形式をしていればひとまず良いということです。
ですので、以降は

$\Large \displaystyle J=\int_0^\infty l(x,u,t)dt$

ということにしておきます。（積分の中身は、なるべく小さいほうが嬉しいような値という意味を込めて、loss の頭文字「l」を使いました）

3.価値関数

ここで、上記の評価関数を別の視点で見てみることにします。
どういう視点であるかというと、実は制御が開始してから既に時刻が「t」だけ経過しているようなケースです。
過去にはもう戻れないので、時刻「t」以降についてのみの評価関数を考えます。
それは単純に下記のような形になるでしょう。

$\Large \displaystyle J_t=\int_t^\infty l(x,u,t)dt$

単に積分の開始が「$t$」からになっているだけです。
既に時刻が「$t$」だけ経過したあとに、これからできる限りの努力をしようと考えた場合には上記の評価関数を考えるのが自然です（もちろん積分の中身については、問題に応じて別のものでも構いません）。

さて、既に時刻が経過しているので、本来の評価関数（積分範囲が「0→∞」）に関しての最適な「$u$」は求められないかもしれませんが、これ以降の時刻についてのみ考えた「$J_t$」を最小化するようなことを目指そうということです。
このときに、「$J_t$」の最小値がu を変えてみることで見つかったとしましょう。
そのような最小値を
$\Large \displaystyle V_t=\min_u J_t$
と書き直し、これを価値関数と呼ぶことにします。
言い換えると価値関数とは、時刻「t」から開始した評価関数「$J_t$」の最小値そのもの、ということになります。
ここで、「$t=0$」としてみます。
要するに時刻は実は一切経過していないという場合を考えるのです。
すると
$\Large \displaystyle V_0=\min_u J_0=\min_u J$
となっており、時刻「t」における価値関数はもともとの評価関数の最小値そのものであるということになります。

今、混乱を避けるために、価値関数「$V_t$」が何を引数にとる関数であるのかを明示することにします。

$\Large \displaystyle V(x(t),t)=\min_u J(x(t),u(t),t)$

となっています。
価値関数は既に「$u$」によって評価関数を最適化したあとなので、もはや「$u$」の関数にはなっていないことに注意してください。
これはすなわち、厳密な価値関数というものは、最適な「$u$」が求まった後にわかるということ。
これが冒頭に述べた、ベルマン方程式は最適な入力が既に得られてとして、その最適入力が満たしているべき方程式である、という話につながる。
そして「$V(x,0)$」のケースを考えれば、これは時刻 t から ∞ までの評価値になるということです。

4.ベルマンの最適性原理

さて、価値関数「$V(x(t),t$)」とは、時刻「$t$」以降の評価関数の（「$u$」を調整した末の）最小値となっている関数でした。
今、この関数を次のように書き直します。
それは「時刻$t$から$t+dt$までの僅かな時間の価値関数と時刻$t+dt$以降の価値関数の和」というものです。

$\Large \displaystyle V(x(t),t)=\min_u J(x(t),u(t),t)$
$\Large \displaystyle =\min_u \int_t^\infty l(x,u,t)dt$
$\Large \displaystyle =\min_u \int_t^{t+dt} l(x,u,t)dt+\min_u \int_{t+dt}^\infty l(x,u,t)dt$
$\Large \displaystyle =\min_u \int_t^{t+dt} l(x,u,t)dt+V(x(t+dt),t+dt)$

ここで最後の式変形は、時刻の開始地点が「$t$」から「$t+dt$」に変わっているだけであることに着目したものです。
これは一言でいえば、積分区間を途中で区切って分けただけなのですが、これはベルマンの最適性原理と呼ばれる非常に重要な因縁に基づいています。
全体の時刻における最適な関数 x(t)があったとすれば、それをとある部分的な時刻で切り取った場合にも最適な関数となっていると言う訳です。

3.3 ベルマン方程式

1.価値関数の離散化

価値関数に関して我々は

$\Large \displaystyle V(x(t),t)=\min_u \int_{t+dt}^t l(x,u,t)dt+V(x(t+dt),t+dt)$

という重要な表現を獲得しました。
今時刻は連続的に変化するものとしていますが、これをあるサンプル時間のもとで「t,t+1,t+2」と数えられる数に変形したとします。
すると上記の式は

$\Large \displaystyle V(x(t),t)=\min_u l(x(t),u(t),t)+V(x(t+1),t+1)$

と近似できます。
ここで、「$l$」という時間関数において、離散化されたために第一項の積分が評価できない。
そこで大雑把に、区間 [$t,t+dt$]の適当な時刻における「$l$」の値を代表して採用している。

以降は、離散化されていることを意識するために「$x(t)$」などを「$x_t$」と表記します。
また、時刻「$t$」が少々煩いので、「$x_t$」と「$V$」と「$l$」に関しては時刻「$t$」を引数に取ることを省略する。
（これは、状態「$x$」らが時刻を引数に取っていることを見れば分かるので問題はないはずだ）
すなわち、ここまでで以降の価値関数を

$\Large \displaystyle V(x_t)=\min_u l(x_t,u_t)+V(x_{t+1})$

と書くことにした。

2.状態の時間発展再訪

さて、制御問題に立ち返ると、「$x$」という状態は、

$\Large \displaystyle \dot{x}=Ax+Bu$

という時間発展をするのでした。
今、これは上記のような線形システムである必要はありません。

$\Large \displaystyle \dot{x}=f(x,u)$

という非線形な何らかの関数で表されてもいいです。
そして離散化されていることを考えると、適当なサンプル時刻を使って

$\Large \displaystyle x_{t+1}=g(x_t,u_t)$

となっていることでしょう。
また、この離散時間の非線形システムは確率システムになっていてもよく、適当な確率分布によって

$\Large \displaystyle x_{t+1}\sim p(x_{t+1}∣x_t,u_t)$

となっているかもしれません。
ここで、「$p$」は現在の状態「$x_t$」と現在の入力「$u_t$」によって次の時刻の状態「$x_{t+1}$」がどのような値を取るのかを表す確率分布であり、状態遷移確率と呼ばれます。
強化学習問題では通常、このような確率システムを使うことが多いので、以降は状態の時間発展が確率システムであるとします。

3.ベルマン方程式

強化学習におけるベルマン（最適）方程式は状態「x」で行動「u」と書いたときには

$\Large \displaystyle V(x_t)=\max_u R(x_t,u_t)+\gamma \max_u \sum_{x_t}p(x_{t+1}∣x,u)V(x_{t+1})$

と書かれます（項を分けているのは意図的です）。
ここで「$\gamma$」は割引率と呼ばれます（こいつの意義は後に説明します）。
また「$R$」は状態と行動の組に対して与えられるスカラー報酬関数です。

さて、これまで連続システムの制御問題として見てきた価値関数（の離散化）を振り返ります。

$\Large \displaystyle V(x_t)=\min_u l(x_t,u_t)+V(x_{t+1})$

というものでした。
ベルマン方程式と非常に似ている形をしています。

まず、違いを見ると、「$min$」演算子と「$max$」演算子の違いがあります。
これは単に符号を入れ替えれば問題としては等価になるので表面上の違いしかありません。
最適制御では、入力が大きくなってほしくない、とか、状態が目標値から離れてほしくない、などの「値が小さくなってほしい」ということを評価に盛り込むため「$min$」演算が用いられます。
一方強化学習では、こうなったら嬉しいという報酬を与えるために「$max$」演算子が使われています。

次の違いは第二項についている「$\gamma$」という割引率です。
これも最適制御と強化学習で「悪さを評価する」か「良さを評価するか」が異なることに起因していると考えています。
まず最適制御では悪さを「正の値」で評価します。
すなわち、理想的な最適のときには価値関数が 0 を取るのです。
一方で強化学習では良さを評価するために、問題設定で下手をしたら無限大に値が発散してしまいます。
そこで、時刻無限遠方まで総和を取ったときに、値が発散しないように「$0<\gamma<1$」を掛けておくのです。

また、この「$\gamma$」という値を導入することを肯定的に見る場合は、「現在の状態を良くすることに集中し、未来に得られる評価は予想に過ぎないので、少し重みを小さくしておく」と考えた度合いだと捉えられます。
実際、強化学習では確率システムを扱うケースが多く、未来の評価値は現在から見たら（計算ができないという不都合とは関係なく）不明であり、確率的にしか見積もれないので、このような考え方をするのは合理的です。

実際、第二項の総和の式は、状態遷移確率によって、次の状態が確率的に定まるため、次の時刻における価値関数の期待値しか見れないことを表しています。
（確率システムでなければ、時間発展の微分方程式にしたがって次の状態は決まっているので、このような計算は不要）

また、第二項において「$max$」演算子が登場しているのは、「価値関数というのは「$u$」に関して評価関数を最適化した場合に得られる値を返す関数」であることを思い出しましょう。
状態遷移確率に現れる「$p$」の条件「$u_t$」については、その次の時刻の価値の期待値が u に関して最大化されるように選ばれているべきであることを表しています。

4.まとめ

まとめてしまえば、強化学習で出てくるベルマン（最適）方程式とは、最適制御での価値関数を離散化し、確率システムへ拡張したものであると言えます。
評価関数をベルマンの最適性原理で、時間による分割を行い

$\Large \displaystyle V(x_t)=\min_u l(x_t,u_t)+V(x_{t+1})$

を得ます。
評価の仕方は、 min でも max でも問題に応じて自分が設計するものであるので、

$\Large \displaystyle V(x_t)=\max_u l(x_t,u_t)+V(x_{t+1})$

と書いても良いはずです。
そして、この時刻に関する再帰式で定義される価値が、正の値で和を取られ続けた場合に発散しないように「$\gamma$」という定数を再帰の部分に乗算しておきます。

$\Large \displaystyle V(x_t)=\max_u l(x_t,u_t)+\gamma V(x_{t+1})$

更に、今扱っている状態 x というのは確率的に値が決まるシステムであり、

$\Large \displaystyle x_{t+1}\sim p(x_{t+1}∣x_t,u_t)$

と与えられるものとします。
そうすると、次の時刻の価値関数は決定論的には決まっておらず、次に取りうる xt+1 すべてを考慮した値（すなわち期待値） ∑xt+1p(xt+1∣x,u)V(xt+1) としておく必要があり、

$\Large \displaystyle V(x_t)=\max_u l(x_t,u_t)+\gamma \sum_{x_{t+1}}p(x_{t+1}∣x,u)V(x_{t+1})$

と書けます。
そして、価値関数というのはもともと時刻全体に渡る状態を評価する評価関数が u に対して最適化されたときの値を返す関数であったことを思い出します。
右辺第二項、状態遷移確率によって次の時刻の価値期待値を計算する際、現在の状態 u は計算結果である期待値を最大化しているものが選ばれているべきであり、

$\Large \displaystyle V(x_t)=\max_u R(x_t,u_t)+\gamma \max_u \sum_{x_{t+1}}p(x_{t+1}∣x,u)V(x_{t+1})$

となります。

3.4 最後に

今回は最適制御の価値関数から強化学習のベルマン方程式を見てきました。
そのために、強化学習で使われる「状態 s」や「行動 a 」をそれぞれ「x」と「u」で書いてきました。
若干見辛さはあったかもしれませんが、（最適制御の式に馴染みがある人にとっては？）、それが確率システムであるということを使って自然に拡張できたと思います。

が、本当は、連続時間の価値関数をテイラー近似してしまえば、価値関数に関する偏微分方程式「ハミルトン・ヤコビ・ベルマン方程式」と呼ばれる連続時間版のベルマン方程式が得られます。
いや、実際にはハミルトン・ヤコビ・ベルマン方程式の離散時間版が、強化学習でのベルマン最適方程式になっており、本来はその対応を見る方が良さそうでしたが、都合により割愛しました。
（そのため、価値関数の離散化はえらくご都合主義的な変形になっている）
この記事を書くに当たっては下記を大いに参考にしました。
https://www.jstage.jst.go.jp/article/jrsj/36/2/36_36_140/_pdf

4.最適制御の数理

「最適制御の数理」から引用。

4.1 正準方程式とHJB方程式

まずは最適制御理論の基礎となる式を導出しましょう。

1.問題設定

制御対象のシステムの状態をベクトル「$x$」で、制御入力をベクトル「$u$」で表現します。
このシステムが次の1階微分方程式に従うものとします。これを状態方程式と呼びます。

$\Large \displaystyle \frac{dy}{dx}=f(x,u,t)$

制御入力とは制御する側がある程度自在に決定できる変数のことで、目的を達成するような制御入力がどのようなものか考えるものが制御理論です。
最適制御では、以下の評価値を最小化するように時刻「$t_0$」から「$t_f$」までの制御入力「$u$」を設計します。
評価値は時間の関数「$x(t),u(t)$」に対し実数値を返す関数(汎関数)として定義します。

$\Large \displaystyle J[x(t),u(t)]=\varphi(x(t_f)) + \int_{t_0}^{t_f}L(x(t),u(t),t)$

右辺第一項は終端コストと呼ばれ、最後の瞬間における状態を評価します。
右辺第二項の被積分項はステージコストと呼ばれ、ある瞬間における状態と制御入力を評価します。
例えば、状態と目標値との差分(誤差)や制御入力の大きさ(物理的な意味としては消費電力や燃料の消費量など)を評価します。

2.最適性条件

「$J[x(t),u(t)]$」の最適性の必要条件を考えます。
システムは状態方程式「$\dot{x}=f(x,u,t)$」を満たすという拘束があるので、ラグランジュ未定乗数法を用いて考えます。
この拘束は「$t\in [t_0,t_f]$」で常に満たされるべきなので、ラグランジュ未乗定数を「$t$」の関数「$\lambda(t)$」として、

$\Large \displaystyle J[x(t),u(t)]=\varphi(x(t_f)) + \int_{t_0}^{t_f}L(x(t),u(t),t)$

$\Large \displaystyle =\varphi(x(t_f)) + \int_{t_0}^{t_f}L(x(t),u(t),t)$

5.マルコフ連鎖モンテカルロ法

「【MCMC法シリーズその1】MCMC法の概要とメトロポリス法」から引用。

5.1 MCMC法の概要

統計分析では、パラメータの事後分布を推定する際に、MCMC法がよく用いられます。
事後分布は以下のようなベイズの公式で表される確率分布です。

$\Large \displaystyle P(S|X)=\frac{P(X|S)P(S)}{\int P(X|S)P(S)dS}$

ここで「S」は原因（システムの状態など）で、「X」は結果（観測値など）を表します。
上式の左辺は、結果（観測値など）から原因（システムの状態など）を推定する確率（信念）を示しています。
実際のシステムでは、上式の右項の分母の積分計算が、困難な場合が多いようです。
そのため、この計算を行う際に、MCMC法はとても有効な手法となります。
MCMC法は「マルコフ連鎖」と「モンテカルロ法」の2つの概念の組み合わせた単語です。
したがって、「マルコフ連鎖」と「モンテカルロ法」を切り分けて説明していきたいと思います。

5.2 マルコフ連鎖とは

マルコフ連鎖は、「未来の状態は、現在の状態だけで決まり、過去の情報は影響しないという考え方」を意味します。
マルコフ連鎖を数式で表現すると以下のようになります。

$\Large \displaystyle P(X_{t+1}|X_t,..,X_1)=P(X_{t+1}|X_t)$

この式の意味を簡単に述べると、ある時刻「$t+1$」の状態「$X_{t+1}$」は、その直前の時刻の状態「$X_{t}$」だけによって決まり、「$X_{t}$」より前の状態は影響しないということです。
例えば、天気の統計モデルを考える際、ある日の天気を決めるのはその前日の天気だけと割り切り、2日前以前の天気の影響は無視します。

この仮定を導入することで考慮すべき変数が減り、モデルを組み立てやすくできるというメリットがあります。
また、この式のことを「遷移確率」と呼ぶこともあります。

5.3 MCMC法の一例 : メトロポリス法

1.メトロポリス法の概要

これまで説明してきた「モンテカルロ法」（説明は省略）と「マルコフ連鎖」を組み合わせたのが「マルコフ連鎖モンテカルロ (MCMC) 法」です。
すなわち「直前の状態だけを考慮しながら、乱数を用いて近似計算を行う」という方法として理解できます。

このMCMC法の具体的なアルゴリズム例として、「メトロポリス法」を紹介します。
この手法はMCMC法の中でも、最も基本的な方法です。
ある確率密度関数 f(x) に対して、以下の流れで「$f(x)$」の点列を算出します。

初期値「$\theta (t=0)$」を決めます。
次の地点の候補点「$\theta (a)$」をランダムサンプリングします。
「$f(\theta (t)) ＜ f(\theta (a))$」の場合はパラメータを更新し「$\theta(t+1) ＝ \theta(a)$」とします。

対して、「$f(\theta(t)) ＞ f(\theta(a))$」の場合は以下のルール「a」で決めます。

以下の「$r$」を計算し、確率「$r$」でパラメータを更新「$\theta(t+1) ＝ \theta(a)$」します。

$\Large \displaystyle r=\frac{f(\theta(a))}{f(\theta(t))}$

パラメータが更新されなかった場合は、「$\theta(t+1) ＝ \theta(t)$」とします。
上記の2.と3.の操作を一定回数繰り返します。

この方法は「山登り」に例えることができます。
山を確率密度関数 f(x) とした場合、山頂に向かってひたすら登り、山頂付近ではなるべく留まろうとする登山者の動きとしてイメージできます。

メトロポリス法のサンプリングのイメージ図

確率密度関数の高い位置になるべく長くいることで、確率が高い部分を重点的にサンプリングすることに繋がります。

2.まとめ

本記事ではMCMC法の基礎的な概念についてまとめ、基本的な手法の1つであるメトロポリス法について紹介しました。
MCMC法の大まかな計算の流れは以下の通りです。

3.MCMC法の大まかな計算の流れ

パラメータを「$\theta$」とし、その初期値「$\theta(t=0)$」を定めます。
ある規則に従って、パラメータの候補値「$\theta(a)$」を得ます。
現在のパラメータを更新するか「$\theta(t+1)=\theta(a)$」、もしくは更新しないか「$\theta(t+1)=\theta(t)$」を確率的に決定します。
2.と3.の操作を一定回数繰り返します。

この 2 と 3 の部分の違いによって、様々な手法が提案されています。
今回はメトロポリス法だけ紹介しましたが、次回の記事ではハミルトニアンモンテカルロ (Hamiltonian Monte Carlo：HMC) 法という、より高度な手法について紹介します。

5.4 MCMC法の一例 : Hamiltonian Monte Carlo（HMC）法

1.Hamiltonian Monte Carlo（HMC）

HMCもMetropolis法の一種であり、同様の手続きを踏むことになるが、名前の通りHamiltonianの概念が登場する

2.Leapfrog Method

3.Langevin Equation

6.Chat GPT

「ChatGPT APIとは？料金やできること、キーの取得方法を解説！」から引用。

6.1 ChatGPT APIの始め方

1.ChatGPT API keyを取得する

ChatGPT API keyは、下記の5ステップで取得できます。

公式サイトにログイン
「Personal」メニュー内「View API keys」をクリック
「Create new secret key」をクリック
キーに名前を付ける（任意）
発行されたキーをコピーする

2.ステップ1：OpenAI APIの「公式サイト」にログインする

または、2.ステップ1：OpenAI APIの「Welcome to the OpenAI platform」にログインする

3.ステップ2：「Personal」をクリックし、メニュー内「View API keys」をクリックする

4.ステップ3：「Create new secret key」をクリックする

5.ステップ4：キーに名前を付け、「Create secret key」をクリックする

4.2 OpenAI APIの試用期間後の課金方法

「【ChatGPT】OpenAI APIの試用期間後の課金方法！移行手順を解説」より

1.試用期間かどうかの確認方法

まず、「https://platform.openai.com/account/usage」を開きます。

試用期間中の場合：「EXPIRES」に未来の日付が入っている場合は無料版として利用できます。
試用期間切れの場合：「Expired（期限切れ）」と表示されている場合は、試用期間が切れています。

2.OpenAI APIの料金体系について

体験版について

OpenAIのAPIは、デフォルトでは体験版となっています。体験版には、無料クレジットが18ドル分付与されています。
ただし、無料クレジットには有効期限があり、期限を過ぎると無料クレジットが残っていたとしてもそれは無効になります。期間終了後は、体験版のままでは利用できなくなります。

有料版について

OpenAI APIの有料版は従量課金制です。固定で月額いくらかかるという形式ではなく、使ったら使った分だけ料金が発生するシステムです。
OpenAI APIへの課金は、Web版のChatGPT（ChatGPT Plus）とは異なります。Web版で課金していたとしても、OpenAI APIへは別途課金する必要があります。こちらは次章にて解説します。