方策勾配法でLQRを解く

LQRの問題はコスト関数の最小化問題だが\(u = \theta \cdot x\) とおくと、方策勾配法の要領でコスト関数\(J\)のヤコビアンを使って最急降下法にて

\[ \theta \leftarrow \theta - \alpha \nabla_{\theta} J \]

を繰り返すことにより、最適な状態フィードバック制御のパラメータを求めることができる(できそう)。

問題設定

制御対象の離散システムとして

\[ x_{i+1} = f(x_i, u_i) \quad i \in [0, N] \]

を考える。このシステムにおいて\(u_i = \pi_{\theta}(x_i) \)の制御則を適用した上でコスト関数

\[ J^{\pi_{\theta}}(x_k, u_k; k, N) = \sum_{i=k}^{N-1} L(x_i, u_i) + L_f(x_N) \]

を考え、\(J^{\pi_{\theta}}(x_0, u_0)\)を最小化するパラメータ \(\theta^*\)を求めたい。ここで、\(L\)はステージコスト、\(L_f\)は終端コストである。確認のため、定義から

\[ J^{\pi_{\theta}}(x_k, u_k; k, N) = L(x_k, u_k) + J^{\pi_{\theta}}(x_{k+1}, u_{k+1}; k+1, N-1) \]

である。以降、簡単のため、単に各関数の添字を省略し\(\pi\)や\(J\)と書く。

コスト関数勾配の導出

状態変数勾配の導出

コスト関数の勾配導出に先立ち、\(\nabla_{\theta} x_k \)を計算する。\(\nabla_{\theta} x_0 = 0\)と\(\nabla_{\theta} u_0 = \nabla_{\theta} \pi (x_0)\)は自明。以下、\(k \geq 1\)について、

\[ \begin{aligned} \nabla_{\theta} x_k &= \nabla_{x} f \cdot \nabla_{\theta} x_{k-1} + \nabla_{u} f \cdot \nabla_{\theta} u_{k-1} \\ &= \nabla_{x} f \cdot \nabla_{\theta} x_{k-1} + \nabla_{u} f \cdot \nabla_{x} \pi \cdot \nabla_{\theta} x_{k-1} + \nabla_{u} f \cdot \nabla_{\theta} \pi\\ &= \nabla_{u} f \cdot \nabla_{\theta} \pi + (\nabla_{x} f + \nabla_{u} f \cdot \nabla_{x} \pi ) \cdot \nabla_{\theta} x_{k-1} \\ &= \sum_{i=0}^{k-1} \left( \prod_{j=i}^{k-2} \nabla_{x} f \big\vert_{x_{j+1}, \pi(x_{j+1})} + \nabla_{u} f \big\vert_{x_{j+1}, \pi(x_{j+1})} \cdot \nabla_{x} \pi \big\vert_{x_{j+1}} \right) \cdot \nabla_{u} f \big\vert_{x_{i}, \pi(x_{i})} \cdot \nabla_{\theta} \pi \big\vert_{x_{i}} \\ \end{aligned} \]

コスト関数勾配

先ほどと同様\(\nabla_{\theta} u_k = \nabla_{\theta} \pi \vert_{x_k} + \nabla_{x} \pi \vert_{x_k} \cdot \nabla_{\theta} x_{k} \)に注意して、

\[ \begin{aligned} \nabla_{\theta} J \vert_{x_0, \pi(x_0)} &= \nabla_{x} L \vert_{x_0, u_0} \cdot \nabla_{\theta} x_0 + \nabla_{u} L \vert_{x_0, u_0} \cdot \nabla_{\theta} u_0 + \nabla_{\theta} J \vert_{x_1, u_1} \\ &=\displaystyle \sum_{k=0}^{N-1} \Big( \nabla_{x} L \vert_{x_k, u_k} \cdot \nabla_{\theta} x_k + \nabla_{u} L \vert_{x_k, u_k} \cdot \nabla_{\theta} u_k \Big) + \nabla_{\theta} L_f \vert_{x_{N}} \\ &\begin{split} =& \displaystyle \sum_{k=0}^{N-1} \big( \nabla_{x} L \vert_{x_k, u_k} + \nabla_{u} L \vert_{x_k, u_k} \cdot \nabla_{x} \pi \vert_{x_k} \big) \cdot \nabla_{\theta} x_k \\ &+ \displaystyle \sum_{k=0}^{N-1} \nabla_{u} L \vert_{x_k, u_k} \cdot \nabla_{\theta} \pi \vert_{x_k} \\ &+ \nabla_{x} L_f \vert_{x_{N}} \cdot \nabla_{\theta} x_{N} \\ \end{split} \end{aligned} \]

線形時不変(LTI)システムでの検算

拘束条件なしの線形フィードバック問題を例にとって気休め程度に検算する。ただし、\(x \in \mathbb{R}^n, u \in \mathbb{R}^m \)とする。

\[ \begin{aligned} L(x, u) &= \frac{1}{2} x^\intercal Q x + \frac{1}{2} u^\intercal R u \\ L_f(x) &= \frac{1}{2} x^\intercal P x\\ x_{k+1} &= A x_k + B u_k \\ u_k &= - K x_k \end{aligned} \]

このとき、

\[ \begin{aligned} x_{k+1} &= (A-BK) x_{k} \\ &= (A - BK)^k x_0 \\ u_{k} &= -K(A - BK)^{k-1} x_0 \\ \end{aligned} \]

を用いて

\[ \begin{aligned} \nabla_{K} x_{i} &= (i-1)(A - BK)^{i-2} B \sum_{j=0,d=0}^{n, m} \frac{\partial K}{\partial k_{j,d}} x_0 \\ &= \sum_{l=0}^{k-1} (A-BK)^{k-l-2} B \sum_{j=0,d=0}^{n, m} \frac{\partial K}{\partial k_{j,d}} x_l \\ \nabla_{K} J &= \nabla_{K} \left( \sum_{i=0}^{N-1} \frac{1}{2} x_i^{\intercal} \left( Q + K^{\intercal} R K \right) x_i + \frac{1}{2} x_N^\intercal P x_N \right)\\ &= \sum_{i=0}^{N-1} \nabla_{K} x_i^{\intercal} \left( Q + K^{\intercal} R K \right) x_i + \sum_{i=0}^{N-1} \sum_{j=0,d=0}^{n, m} x_i^{\intercal} \frac{\partial K}{\partial k_{j,d}} R K x_i + \nabla_{K} x_N^\intercal P x_N \end{aligned} \]

であるので確かに一致している。

さらに、\(N \rightarrow \infty \)のとき、Global Convergence of Policy Gradient Methods for the Linear Quadratic Regulatorでも示されている通り、

\[ P_K = Q + K^\intercal R K + (A-BK)^\intercal P_K (A-BK) \]

と置くと、

\[ \begin{aligned} \nabla_{K} J &= \big((R+B^\intercal P_K B)K - B^\intercal P_K A \big) \sum_{i=0}^{\infty} x_i^\intercal x_i \end{aligned} \]

と求まる。上記論文ではモデルフリーの自然勾配法との収束速度比較などを行っているので参照すると面白い。