確率推論としての最適制御

SICE制御理論若手合宿2024 確率推論としての最適制御 名古屋大学 本田 康平 https://kohonda.github.io/ 1

• https://kohonda.github.io/

内容 1. 確率推論としての最適制御の概要 2. 「最適な」制御分布とは？ 3. 最適制御分布の考察 4. 最適制御分布から解をサンプルする方法 5. 応用・発展 2

内容 1. 確率推論としての最適制御の概要 2. 「最適な」制御分布とは？ 3. 最適制御分布の考察 4. 最適制御分布から解をサンプルする方法 5. 応用・発展 3

「由緒正しい」最適制御 数理モデルで予測 数理最適化 最適制御問題 min 𝑥𝑥0:𝑇𝑇 ,𝑢𝑢0:𝑇𝑇−1 s. t. 𝐽𝐽(𝑥𝑥0:𝑇𝑇 , 𝑢𝑢0:𝑇𝑇−1 ) 𝑥𝑥0 = 𝑥𝑥(𝑡𝑡) 𝑥𝑥𝑡𝑡+1 = 𝑓𝑓(𝑥𝑥𝑡𝑡 , 𝑢𝑢𝑡𝑡 ) 𝑔𝑔 𝑥𝑥𝑡𝑡 , 𝑢𝑢𝑡𝑡 ≤ 0 数理最適化によって最適制御問題を解く ∀ 𝑡𝑡 ∈ 0, … , 𝑇𝑇 • 解ける問題のクラスが制限される (x 非線形 x 解析微分不可) • 数理最適化ソルバが実質ブラックボックス • 込み入ったことをしようとすると魔改造・鬼調整が発生 4

確率推論として最適制御問題を解く 推論時には勾配が不要 (なことが多い) 解の確率的な揺らぎを表現できる ソルバーの実装や並列化が簡単 (なことが多い) Predictive horizon 𝑇𝑇 𝑢𝑢𝑡𝑡+3 𝑢𝑢𝑡𝑡+2 𝑢𝑢𝑡𝑡+1 𝑢𝑢𝑡𝑡 ① 𝑢𝑢𝑡𝑡+3 𝑢𝑢𝑡𝑡+2 𝑢𝑢𝑡𝑡+1 𝑢𝑢𝑡𝑡 ~ 𝜋𝜋 ∗ (𝑢𝑢𝑡𝑡:𝑡𝑡+3 ) 𝑇𝑇 ② 𝑢𝑢𝑡𝑡 最適制御分布を ①どのように表現し、②どのように解をサンプルするか︖ 5

実装は簡単 例 : モデル予測経路積分制御 (MPPI) (G. Williams+ 2018) 1. ガウスノイズを付与した𝐾𝐾個の解をサンプル : 𝑢𝑢𝑡𝑡𝑘𝑘 = 𝑢𝑢� 𝑡𝑡 + 𝜖𝜖, 𝜖𝜖~𝑁𝑁(0, Σ) 𝑘𝑘 𝑘𝑘 = 𝐽𝐽 𝑥𝑥0:𝑇𝑇 , ⋅ , s. t. 𝑥𝑥𝑡𝑡+1 = 𝑓𝑓(𝑥𝑥𝑡𝑡 , 𝑢𝑢𝑡𝑡 ) 2. 各軌道のコスト関数を計算 : 𝐽𝐽𝜏𝜏𝑘𝑘 𝑢𝑢0:𝑇𝑇−1 3. 重みづけ平均による解の更新 : 𝑢𝑢� 𝑡𝑡+1 ← ∑𝑘𝑘 Softmax −𝜆𝜆−1 𝐽𝐽𝜏𝜏𝑘𝑘 𝑢𝑢𝑡𝑡𝑘𝑘 https://github.com/kohonda/proj-svg_mppi https://github.com/kohonda/mppi_playground 6

• https://github.com/kohonda/mppi_playground
• https://github.com/kohonda/proj-svg_mppi

内容 1. 確率推論としての最適制御の概要 2. 「最適な」制御分布とは︖ 3. 最適制御分布の考察 4. 最適制御分布から解をサンプルする方法 5. 応用・発展 7

最適制御分布はどのように表現される? 最適制御問題 min 𝑥𝑥0:𝑇𝑇 ,𝑢𝑢0:𝑇𝑇−1 s. t. 𝑢𝑢𝑡𝑡+3 … 𝑢𝑢𝑡𝑡+2 𝑢𝑢𝑡𝑡+1 𝑢𝑢𝑡𝑡 ~ 𝜋𝜋 ∗ (𝑢𝑢0:𝑇𝑇−1 ) 𝐽𝐽(𝑥𝑥0:𝑇𝑇 , 𝑢𝑢0:𝑇𝑇−1 ) 𝑥𝑥0 = 𝑥𝑥(𝑡𝑡) 𝑥𝑥𝑡𝑡+1 = 𝑓𝑓(𝑥𝑥𝑡𝑡 , 𝑢𝑢𝑡𝑡 ) 今回は無視︕𝑔𝑔 𝑥𝑥 , 𝑢𝑢 𝑡𝑡 𝑡𝑡 ≤0 𝑢𝑢𝑡𝑡 いくつかの仮定を置くと、最適制御分布の確率密度関数は求まる 最適制御分布の確率密度 𝜋𝜋 ∗ 𝑢𝑢0:𝑇𝑇−1 = 𝑍𝑍 −1 exp −𝜆𝜆−1 𝐽𝐽𝜏𝜏 𝑢𝑢0:𝑇𝑇−1 ボルツマン分布 𝑝𝑝(𝑢𝑢0:𝑇𝑇−1 ) 事前分布 ※ 同様の最適分布形状は他の多くのトピックでも現れる (ボルツマンマシン、Direct Preference Optimization) 8

1/3 ∗ 最適制御分布𝝅𝝅 の導出 (S. Levine 2018)(Z. Wang+ 2021)など (1) 最適性変数𝒪𝒪𝒕𝒕 によるグラフィカルモデル化 最適制御問題 min 𝑥𝑥0:𝑇𝑇 ,𝑢𝑢0:𝑇𝑇−1 s. t. 𝐽𝐽(𝑥𝑥0:𝑇𝑇 , 𝑢𝑢0:𝑇𝑇−1 ) 最適(𝒪𝒪 = 1)か否か(𝒪𝒪 = 0)の最適性変数 𝑥𝑥0 = 𝑥𝑥(𝑡𝑡) 𝑥𝑥𝑡𝑡+1 = 𝑓𝑓(𝑥𝑥𝑡𝑡 , 𝑢𝑢𝑡𝑡 ) 𝒪𝒪0 𝒪𝒪1 𝒪𝒪2 𝑢𝑢0 𝑢𝑢1 𝑢𝑢2 𝑥𝑥0 𝑥𝑥1 ベイズの定理より (軌道𝜏𝜏 ≔ [𝑥𝑥0:𝑇𝑇 , 𝑢𝑢0:𝑇𝑇−1 ]) 𝑝𝑝 𝜏𝜏|𝒪𝒪 = 1 ∝ 𝑝𝑝 𝒪𝒪 = 1 𝜏𝜏 × 𝑝𝑝 𝑥𝑥0 � 最適な軌道の分布 最適性分布 𝑇𝑇−1 𝑡𝑡=1 𝑥𝑥2 𝑝𝑝 𝑥𝑥𝑡𝑡+1 𝑥𝑥𝑡𝑡 , 𝑢𝑢𝑡𝑡 𝑝𝑝 𝑢𝑢𝑡𝑡 軌道𝜏𝜏の分布 … (1) 9

2/3 ∗ 最適制御分布𝝅𝝅 の導出 (S. Levine 2018)(Z. Wang+ 2021)など (2) 最適軌道分布𝒑𝒑(𝝉𝝉|𝓞𝓞 = 𝟏𝟏)を𝓞𝓞を使わずに𝝅𝝅(𝝉𝝉)で近似 𝑝𝑝∗ 𝜏𝜏 = min 𝔻𝔻KL 𝜋𝜋 𝜏𝜏 ||𝑝𝑝(𝜏𝜏|𝒪𝒪 = 1) 𝜋𝜋 (KL divergenceの最小化) = min 𝔼𝔼𝜋𝜋(𝜏𝜏) log 𝜋𝜋(𝜏𝜏) − log(𝑝𝑝(𝜏𝜏|𝒪𝒪 = 1) 𝜋𝜋 (定義より) 𝜋𝜋 ∗ 𝑢𝑢0:𝑇𝑇−1 = min 𝔼𝔼𝜋𝜋(𝜏𝜏) 𝔼𝔼𝑝𝑝(𝑥𝑥0:𝑇𝑇 |𝑢𝑢0:𝑇𝑇−1) − log 𝑝𝑝 𝒪𝒪 = 1 𝜏𝜏 最適制御分布 𝜋𝜋 最適性分布としてボルツマン分布を仮定 𝑝𝑝 𝒪𝒪 = 1 𝜏𝜏 ≔ 𝑍𝑍 −1 exp −𝜆𝜆−1 𝐽𝐽𝜏𝜏 (𝑢𝑢0:𝑇𝑇−1 ) (式 (1) & 軌道の分布より) + 𝔻𝔻KL (𝜋𝜋(⋅)||𝑝𝑝(⋅)) 𝐽𝐽𝜏𝜏 𝑢𝑢0:𝑇𝑇−1 ≔ 𝔼𝔼𝑝𝑝 𝑥𝑥0:𝑇𝑇 𝑢𝑢0:𝑇𝑇−1 𝐽𝐽(𝑥𝑥0:𝑇𝑇 , 𝑢𝑢0:𝑇𝑇−1 ) = min 𝔼𝔼𝜋𝜋(𝜏𝜏) 𝜆𝜆−1 𝐽𝐽𝜏𝜏 (𝑢𝑢0:𝑇𝑇−1 ) + 𝔻𝔻KL (𝜋𝜋(𝑢𝑢0:𝑇𝑇−1 )||𝑝𝑝(𝑢𝑢0:𝑇𝑇−1 )) 𝜋𝜋 コスト関数の期待値 (2) 事前分布と乖離に対する正則化項 ※ ELBOの最小化によっても導出可能 10

3/3 ∗ 最適制御分布𝝅𝝅 の導出 (S. Levine 2018)(Z. Wang+ 2021)など (3) ラグランジュ未定乗数法による𝝅𝝅∗ の算出 式 (2) を ∫ 𝜋𝜋 𝑢𝑢0:𝑇𝑇−1 d𝑢𝑢0:𝑇𝑇−1 = 1 の制約の下で最小化する ラグランジュ関数 𝐿𝐿(𝜋𝜋, 𝛼𝛼) = 𝔼𝔼𝜋𝜋 𝜏𝜏 𝜆𝜆−1 𝐽𝐽𝜏𝜏 𝑢𝑢0:𝑇𝑇−1 + 𝔻𝔻KL (𝜋𝜋(𝑢𝑢0:𝑇𝑇−1 )| 𝑝𝑝 𝑢𝑢0:𝑇𝑇−1 ∗ (𝑢𝑢 𝜕𝜕𝐿𝐿 𝜋𝜋 0:𝑇𝑇−1 ) = 𝜆𝜆−1 𝐽𝐽𝜏𝜏 + log − 𝛼𝛼 = 0 𝜕𝜕𝜋𝜋 𝑝𝑝(𝑢𝑢0:𝑇𝑇−1 ) より � 𝜋𝜋 d𝑢𝑢0:𝑇𝑇−1 = exp 𝛼𝛼 − 1 � exp −𝜆𝜆−1 𝐽𝐽𝜏𝜏 d𝑢𝑢0:𝑇𝑇−1 = 1 + 𝛼𝛼 (� 𝜋𝜋 𝑢𝑢0:𝑇𝑇−1 d𝑢𝑢0:𝑇𝑇−1 − 1) 𝜋𝜋 ∗ (𝑢𝑢0:𝑇𝑇−1 ) = exp(𝛼𝛼 − 1)exp(−𝜆𝜆−1 𝐽𝐽𝜏𝜏 𝑢𝑢0:𝑇𝑇−1 ) より Z = exp −𝛼𝛼 + 1 = � exp −𝜆𝜆−1 𝐽𝐽𝜏𝜏 d𝑢𝑢0:𝑇𝑇−1 最適制御分布の確率密度 𝜋𝜋 ∗ 𝑢𝑢0:𝑇𝑇−1 = 𝑍𝑍 −1 exp −𝜆𝜆−1 𝐽𝐽𝜏𝜏 𝑢𝑢0:𝑇𝑇−1 最適性分布 (ボルツマン分布) 𝑝𝑝(𝑢𝑢0:𝑇𝑇−1 ) 事前分布 (3) 11

∗ 最適制御分布𝝅𝝅 の導出 (S. Levine 2018)(Z. Wang+ 2021)など 補足 ボルツマン分布による最適性分布の仮定 最適性分布 (軌道𝜏𝜏が得られたときの最適な確率) ボルツマン分布: 𝑝𝑝 𝒪𝒪 = 1 𝜏𝜏 ≔ 𝑍𝑍 −1 exp −𝜆𝜆−1 𝐽𝐽𝜏𝜏 (𝑢𝑢0:𝑇𝑇−1 ) 正規化定数 コスト関数 温度パラメータ コスト関数 例: ボルツマン分布 軌道 軌道 コスト関数の低い軌道ほど高い確率密度 (=高い最適性) となる ※ ボルツマン分布はコスト𝐽𝐽𝜏𝜏 の情報エントロピー最大化から導かれる 12

内容 1. 確率推論としての最適制御の概要 2. 「最適な」制御分布とは︖ 3. 最適制御分布の考察 4. 最適制御分布から解をサンプルする方法 5. 応用・発展 13

1/4 ∗ 最適制御分布𝝅𝝅 の考察 最適な行動分布 = (最適性分布) × (入力の事前分布) Drunken Spider (H J Kappen 2005) 酔い具合 = 事前分布のノイズレベル ※ 酔っている時に狭い橋を渡ると 池に落ちてしまう 酔っていない時の 最適経路 池 クモ 酔いが酷い時の 最適経路 入力ノイズ (事前分布) によって最適な行動が変化する 14

∗ 最適制御分布𝝅𝝅 の考察 2/4 「対称性の破れ」による意思決定の遅延 (H J Kappen 2005) 分岐点に近づくにつれて𝝅𝝅∗ のモード数が減少する (対称性の破れ) 𝜋𝜋 ∗ 𝜋𝜋 ∗ 𝜋𝜋 ∗ 時刻 対称性が破れることで最適分布のモードが一つに絞られる ※ 決定論的な制御よりも意思決定を先延ばしにする傾向有り 15

∗ 最適制御分布𝝅𝝅 の考察 3/4 温度パラメータ𝝀𝝀による分布形状の変化 コスト関数 𝜋𝜋 ∗ = 𝑍𝑍 −1 exp −𝝀𝝀−𝟏𝟏 𝐽𝐽𝜏𝜏 𝑝𝑝(𝑢𝑢0:𝑇𝑇−1 ) 最適性分布 𝜆𝜆 = 0.1 小 𝜆𝜆 = 1 温度パラメータ𝜆𝜆 𝜆𝜆 = 10 大 温度パラメータを小さくするほどOne-hotな分布形状となる ※ 式 (2) から解釈すると、正則化項の影響が小さくなる 16

4/4 ∗ 最適制御分布𝝅𝝅 の考察 サンプル複雑性の支配要因 𝜋𝜋 ∗ = 𝑍𝑍 −1 exp −𝜆𝜆−1 𝐽𝐽𝜏𝜏 𝑝𝑝(𝑢𝑢0:𝑇𝑇−1 ) 𝝅𝝅∗ を精度良く近似するために必要なサンプル数 (=サンプル複雑性) は温度パラメータ𝝀𝝀とコスト関数𝑱𝑱𝝉𝝉 に依存する ■ (前提) 分布がOne-hotであるほどサンプル複雑性は高まる → 𝜆𝜆が小さく、 𝐽𝐽𝜏𝜏 が大きいほど多くのサンプル数が必要 =事前分布と大きく異なるような最適分布 (式 (2)より) コストが大きな状態 17

内容 1. 確率推論としての最適制御の概要 2. 「最適な」制御分布とは︖ 3. 最適制御分布の考察 4. 最適制御分布から解をサンプルする方法 5. 応用・発展 18

最適制御分布から解をサンプルする 𝑢𝑢𝑡𝑡+1 𝑢𝑢𝑡𝑡 ~ 𝜋𝜋 ∗ (𝑢𝑢𝑡𝑡:𝑡𝑡+3 ) 最適制御分布 サンプル 𝑇𝑇 ロボットに印加 𝑢𝑢𝑡𝑡 𝜋𝜋 ∗ (𝑢𝑢0:𝑇𝑇−1 ) = 𝑍𝑍 −1 exp −𝜆𝜆−1 𝐽𝐽𝜏𝜏 (𝑢𝑢0:𝑇𝑇−1 ) 𝑝𝑝(𝑢𝑢0:𝑇𝑇−1 ) ■問題 : 𝜋𝜋 ∗ は複雑なので厳密な推論は手に負えない (Interactable) = 𝜋𝜋 ∗ から直接サンプルをすることは非現実的 → アプローチ : MCMC、ラプラス近似、変分推論 並列化が容易なので(︖) 現在最も主流 19

変分推論による最適分布の近似 ■ アイデア : 簡単な変分分布𝝅𝝅𝜽𝜽 によって𝝅𝝅∗ を近似する 𝜋𝜋 ∗ 手に負えない︕ min 𝔻𝔻𝐾𝐾𝐾𝐾 (𝜋𝜋 ∗ ||𝜋𝜋𝜃𝜃 ) 𝜃𝜃 𝜋𝜋𝜃𝜃∗ 𝜋𝜋𝜃𝜃 サンプル容易 ※ 𝜋𝜋𝜃𝜃 の選び方は任意 (例 : 多変量ガウス分布、混合ガウス分布) 20

例 : MPPI (G. Williams+ 2018) 𝒕𝒕 ∏ 最適制御分布𝝅𝝅 を多変量ガウス分布𝝅𝝅𝝁𝝁𝜽𝜽 = 𝒕𝒕 𝓝𝓝(𝝁𝝁𝜽𝜽 , 𝚺𝚺𝒕𝒕 )で近似 𝜇𝜇 𝜃𝜃 ∗ ∗ ∗ = min 𝔻𝔻KL 𝜋𝜋 ||𝜋𝜋𝜇𝜇𝜃𝜃 𝜇𝜇𝜃𝜃 = min 𝔼𝔼𝜋𝜋∗ log 𝜋𝜋 ∗ − log 𝜋𝜋𝜇𝜇𝜃𝜃 𝜇𝜇𝜃𝜃 = max 𝔼𝔼𝜋𝜋∗ � log 𝒩𝒩(𝜇𝜇𝜃𝜃𝑡𝑡 , Σ𝑡𝑡 ) = max � 𝜋𝜋 ∗ log 𝒩𝒩 ⋅ d𝑢𝑢0:𝑇𝑇−1 𝜇𝜇𝜃𝜃 𝜕𝜕ℒ𝜃𝜃 ∗ 𝑡𝑡 = 0 ⟺ � 𝜋𝜋 ⋅ 𝜕𝜕𝜇𝜇𝜃𝜃 𝑡𝑡 𝑡𝑡 𝑢𝑢𝑡𝑡 − 𝜇𝜇𝜃𝜃∗ 凸関数! d𝑢𝑢𝑡𝑡 = 0 𝜇𝜇𝜃𝜃 ℒ𝜃𝜃 大域最適解 −1 −1 exp 𝜆𝜆 𝐽𝐽 𝑝𝑝(𝑢𝑢 )d𝑢𝑢 𝑢𝑢 ∫ 𝔼𝔼 exp −𝜆𝜆 𝐽𝐽𝜏𝜏 𝑢𝑢𝑡𝑡 𝑡𝑡 𝜏𝜏 𝑡𝑡 𝑡𝑡 𝑝𝑝 𝑢𝑢𝑡𝑡 𝑡𝑡 = ⟺ 𝜇𝜇𝜃𝜃∗ = −1 𝔼𝔼𝑝𝑝 𝑢𝑢𝑡𝑡 [exp −𝜆𝜆−1 𝐽𝐽𝜏𝜏 ] ∫ exp 𝜆𝜆 𝐽𝐽𝜏𝜏 𝑝𝑝(𝑢𝑢𝑡𝑡 )d𝑢𝑢𝑡𝑡 ※ Σ𝑡𝑡 を同時に推定することも可能だが、大域最適性は失われる (モンテカルロ法により推定可能) 21

内容 1. 確率推論としての最適制御の概要 2. 「最適な」制御分布とは︖ 3. 最適制御分布の考察 4. 最適制御分布から解をサンプルする方法 5. 応用・発展 22

変分推論MPCの分布表現の種類 事前分布と変分分布の分布が解の性質・性能を左右 最適制御分布 : 𝜋𝜋 ∗ (𝑢𝑢0:𝑇𝑇−1 ) = 𝑍𝑍 −1 exp −𝜆𝜆−1 𝐽𝐽𝜏𝜏 (𝑢𝑢0:𝑇𝑇−1 ) 𝒑𝒑(𝒖𝒖𝟎𝟎:𝑻𝑻−𝟏𝟏 ) ∗ ∗ 𝜋𝜋 𝑢𝑢 = min 𝔻𝔻 𝜋𝜋 ||𝝅𝝅𝜽𝜽 変分推論 : 𝜃𝜃 0:𝑇𝑇−1 KL 𝜃𝜃 変分分布 ≒ 解の多様性 (モダリティ) 変分分布 事前分布 事前分布 ≒ 解の探索範囲 • ガウス分布 (MPPI 2018)(CEM 2003) • 前時刻で最適化した分布 (MPPI 2018), etc. • 混合ガウス分布 (Osa 2020)(VI-MPC 2020)(Wang+ 2021) • Normalizing Flow (Flow-MPPI 2022)(Sacks+ 2023) • パーティクル表現(SV-MPC 2020)(DuSt-MPC 2021) • 強化学習方策分布 (TD-MPC 2022) ※ モデルフリーRLの文脈では(MPO 2018)(SAC 2018)などで同じ枠組みが利用されている (Control as Inference)23

KL divergenceの非対称性 非対称性を利用して「行動決定の遅延」を低減 ∗ ∗ 𝜋𝜋 𝑢𝑢 = min 𝔻𝔻 𝜋𝜋 ||𝝅𝝅𝜽𝜽 変分推論 : 𝜃𝜃 0:𝑇𝑇−1 KL 𝜃𝜃 Forward KL div. 確率密度 Forward KL div.最小化では「最適解」が複数モードを覆いかぶさる その結果、行動決定が更に遅延 𝝅𝝅∗ 𝝅𝝅∗𝜽𝜽 Forward 𝝅𝝅∗𝜽𝜽 Reverse → Reverse KL div. 𝔻𝔻𝐊𝐊𝐊𝐊 𝝅𝝅𝜽𝜽 ||𝝅𝝅∗ を利用してmode-seekな解を算出 (Kobayashi+ 2022)(Honda+ 2024) 24

拡散モデル as 最適化ソルバ (C Pan+ 2024) 拡散モデルを用いて最適制御分布からサンプル (学習不要) • デノイジングに用いるスコア関数が最適制御問題から推定可能 • MPPIの更新則によってスコア関数の一部を推定する → ただし、1回のサンプルにMPPIを拡散ステップ数 (100~1000) 回解く必要がある… 最適制御分布からのサンプル デノイジング デノイジング ノイズ𝒩𝒩(0, 𝐼𝐼) デノイジングされる様子 https://zenn.dev/takuya_fukatsu/articles/91c5dbd785cec9 25

• https://zenn.dev/takuya_fukatsu/articles/91c5dbd785cec9

モデルベース強化学習における役割 プランニングパートで探索と活用を担う 表現学習 環境 環境モデル 確率推論MPC • 環境ダイナミクス • 報酬/制約モデル, etc. プランニング エージェント アクション 潜在空間でCEM(PlaNet 2019) / MPPI(TD-MPC 2022) / VI-MPC(Okada+ 2020) 多様なデータを集めるために探索 & 学習した環境モデルを活用 26