Revisiting LQR Control from the Perspective of Receding-Horizon Policy Gradient

TL;DR

This paper analyzes the discrete-time LQR problem using receding-horizon policy gradient, providing sample complexity bounds and demonstrating the method's effectiveness in control and estimation without prior stabilizing policies.

Contribution

It offers a detailed sample complexity analysis of RHPG for LQR, showing it can learn near-optimal policies without needing an initial stabilizing control.

Findings

RHPG can learn stabilizing policies close to optimal without initial stabilizing control.

The analysis applies to both control and estimation problems in linear systems.

The method demonstrates broad applicability and streamlined analysis in linear control tasks.

Abstract

We revisit in this paper the discrete-time linear quadratic regulator (LQR) problem from the perspective of receding-horizon policy gradient (RHPG), a newly developed model-free learning framework for control applications. We provide a fine-grained sample complexity analysis for RHPG to learn a control policy that is both stabilizing and $\epsilon$-close to the optimal LQR solution, and our algorithm does not require knowing a stabilizing control policy for initialization. Combined with the recent application of RHPG in learning the Kalman filter, we demonstrate the general applicability of RHPG in linear control and estimation with streamlined analyses.

Equations

$$\Vert X{\Vert }_W^2:=\underset{z\ne 0}{\sup }\frac{z^{\top }{X}^{\top }WXz}{z^{\top }Wz}.$$

$$x_{t+1}=A{x}_t+B{u}_t,$$

$$\displaystyle J_{\infty}:=\mathbb{E}_{x_{0}\sim\mathcal{D}}\left[\sum_{t=0}^{\infty}\big{(}x_{t}^{\top}Qx_{t}+u^{\top}_{t}Ru_{t}\big{)}\right],$$

$$K^{\ast }=(R+B^{\top }{P}^{\ast }B{)}^{-1}{B}^{\top }{P}^{\ast }A,$$

$$P=Q+A^{\top }PA-A^{\top }PB(R+B^{\top }PB{)}^{-1}{B}^{\top }PA.$$

$$\displaystyle\min_{K}~{}J_{\infty}(K)=\mathbb{E}_{x_{0}\sim\mathcal{D}}\left[\sum_{t=0}^{\infty}\big{(}x_{t}^{\top}(Q+K^{\top}RK)x_{t}\big{)}\right]$$

$$\text{s.t.}K\in \mathcal{K}:=\{K\mid \rho (A-BK)<1\}.$$

$$\displaystyle J\big{(}\{K_{t}\}\big{)}:=\mathbb{E}_{x_{0}\sim\mathcal{D}}\left[\sum_{t=0}^{N-1}x_{t}^{\top}(Q+K^{\top}_{t}RK_{t})x_{t}+x_{N}^{\top}Q_{N}x_{N}\right],$$

$$K_t^{\ast }=(R+B^{\top }{P}_{t+1}^{\ast }B{)}^{-1}{B}^{\top }{P}_{t+1}^{\ast }A,$$

$$P_t^{\ast }$$

$$-A^{\top }{P}_{t+1}^{\ast }B(R+B^{\top }{P}_{t+1}^{\ast }B{)}^{-1}{B}^{\top }{P}_{t+1}^{\ast }A.$$

$$\displaystyle N_{0}=\frac{1}{2}\cdot\frac{\log\big{(}\frac{\|Q_{N}-P^{*}\|_{*}\cdot\kappa_{P^{*}}\cdot\|A_{K}^{*}\|\cdot\|B\|}{\epsilon\cdot\lambda_{\min}(R)}\big{)}}{\log\big{(}\frac{1}{\|A_{K}^{*}\|_{*}}\big{)}}+1.$$

$$\underset{K_h}{\min }{J}_h(K_h)$$

$$\displaystyle+x_{h}^{\top}\big{(}Q+K_{h}^{\top}RK_{h}\big{)}x_{h}+x_{N}^{\top}Q_{N}x_{N}\bigg{]}.$$

$$\displaystyle\big{\|}\widetilde{K}_{h}-\widetilde{K}_{h}^{*}\big{\|}\sim\mathcal{O}(\epsilon)\mathcal{O}(\texttt{poly}(\text{system parameters})),$$

$$K_{h,i+1}=K_{h,i}-{\eta }_h\cdot \tilde{\nabla }{J}_h(K_{h,i})$$

$$x_{t+1}=A{x}_t+B{u}_t+w_t,w_t\sim \mathcal{N}(0,W),W>0.$$

$$\displaystyle J_{\infty}:=\limsup_{N\to\infty}\frac{1}{N}~{}\mathbb{E}_{x_{0},w_{t}}\left[\sum_{t=0}^{N-1}\big{(}x_{t}^{\top}Qx_{t}+u^{\top}_{t}Ru_{t}\big{)}\right].$$

$${\tilde{P}}_t:=P_t^{\ast }-P^{\ast },\tilde{R}:=R+B^{\top }{P}^{\ast }B,$$

$$\overline{A}:=A-B{\tilde{R}}^{-1}{B}^{\top }{P}^{\ast }A,$$

$${\tilde{P}}_t$$

$$\displaystyle=\overline{A}^{\top}\widetilde{P}_{t+1}^{1/2}\big{[}I+\widetilde{P}^{1/2}_{t+1}B\widetilde{R}^{-1}B^{\top}\widetilde{P}_{t+1}^{1/2}\big{]}^{-1}\widetilde{P}_{t+1}^{1/2}\overline{A}$$

$$\displaystyle\leq\big{[}1+\lambda_{\min}(\widetilde{P}^{1/2}_{t+1}B\widetilde{R}^{-1}B^{\top}\widetilde{P}_{t+1}^{1/2})\big{]}^{-1}\overline{A}^{\top}\widetilde{P}_{t+1}\overline{A}$$

$$=:{\mu }_t{\overline{A}}^{\top }{\tilde{P}}_{t+1}\overline{A},$$

$$\Vert {\tilde{P}}_t{\Vert }_{\ast }\le \Vert \overline{A}{\Vert }_{\ast }^2\cdot \Vert {\tilde{P}}_{t+1}{\Vert }_{\ast }.$$

$$\Vert {\tilde{P}}_t\Vert \le {\kappa }_{P^{\ast }}\cdot \Vert {\tilde{P}}_t{\Vert }_{\ast }\le {\kappa }_{P^{\ast }}\cdot \Vert \overline{A}{\Vert }_{\ast }^{2(N-t)}\cdot \Vert {\tilde{P}}_N{\Vert }_{\ast },$$

$$\displaystyle N\geq\frac{1}{2}\cdot\frac{\log\big{(}\frac{\|\widetilde{P}_{N}\|_{*}\cdot\kappa_{P^{*}}}{\epsilon}\big{)}}{\log\big{(}\frac{1}{\|\overline{A}\|_{*}}\big{)}}+1.$$

$$K_t^{\ast }-K^{\ast }$$

$$-(R+B^{\top }{P}^{\ast }B{)}^{-1}{B}^{\top }{P}^{\ast }A$$

$$\displaystyle=\big{[}(R+B^{\top}P^{*}_{t+1}B)^{-1}-(R+B^{\top}P^{*}B)^{-1}\big{]}B^{\top}P^{*}A$$

$$+(R+B^{\top }{P}_{t+1}^{\ast }B{)}^{-1}{B}^{\top }(P_{t+1}^{\ast }-P^{\ast })A$$

$$=(R+B^{\top }{P}_{t+1}^{\ast }B{)}^{-1}{B}^{\top }(P^{\ast }-P_{t+1}^{\ast })B{K}^{\ast }$$

$$-(R+B^{\top }{P}_{t+1}^{\ast }B{)}^{-1}{B}^{\top }(P^{\ast }-P_{t+1}^{\ast })A$$

$$=(R+B^{\top }{P}_{t+1}^{\ast }B{)}^{-1}{B}^{\top }(P^{\ast }-P_{t+1}^{\ast })(B{K}^{\ast }-A)$$

$$\Vert K_0^{\ast }-K^{\ast }\Vert \le \frac{\Vert \overline{A}\Vert \cdot \Vert B\Vert }{{\lambda }_{\min }(R)}\cdot \Vert {\tilde{P}}_1\Vert .$$

$$\displaystyle N=\frac{1}{2}\cdot\frac{\log\big{(}\frac{2\|Q_{N}-P^{*}\|_{*}\cdot\kappa_{P^{*}}\cdot\|A_{K}^{*}\|\cdot\|B\|}{\epsilon\cdot\lambda_{\min}(R)}\big{)}}{\log\big{(}\frac{1}{\|A_{K}^{*}\|_{*}}\big{)}}+1,$$

$$P_t^{\ast }$$

Full text

1 Introduction

Model-free policy gradient (PG) methods promise a universal end-to-end framework for controller designs. By utilizing input-output data of a black-box simulator, PG methods directly search the prescribed policy space until convergence, agnostic to system models, objective function, and design criteria/constraints. The general applicability of PG methods leads to countless empirical successes in continuous control, but the theoretical understanding of these PG methods is still in its early stage. Stemmed from the convergence theory of PG methods for general reinforcement learning tasks [1, 2], a recent thrust of research has specialized the analysis for the convergence and sample complexity of PG methods into several linear state-feedback control benchmarks [3, 4, 5, 6, 7, 8, 9]. However, incorporating imperfect-state measurements leads to a deficit of most, if not all, favorable landscape properties crucial for PG methods to converge (globally) in the state-feedback settings [10, 9]. Even worse, the control designer now faces several challenges unique to control applications: a) convergence might be toward a suboptimal stationary point without system-theoretic interpretations; b) provable stability and robustness guarantee could be lacking; c) convergence depends heavily on the initialization (e.g., the initial policy should be stabilizing), which would be challenging to hand-craft; and d) algorithm could be computationally inefficient. These bottlenecks blur the applicability of model-free PG methods in real-world control scenarios since the price for each of the above disadvantages could be unaffordable.

On the other hand, classic theories provide both elegant analytic solutions and efficient computational means (e.g., Riccati recursions) to a wide range of control problems [11, 12, 13]. They further reveal the intricate structure in various control settings and offer system-theoretical interpretations and guarantees to their characterized solutions. They suggest that, compared to viewing the dynamical system as a black box and studying the properties of PG methods from a (nonconvex) optimization perspective, it is better to incorporate those properties unique to decision and control into the design of learning algorithms.

In this work, we revisit the classical linear quadratic regulator (LQR) problem from the perspective of the newly-developed receding-horizon PG (RHPG) framework [14], which integrates Bellman’s principle of optimality into the development of a model-free PG framework. First, RHPG approximates infinite-horizon LQR using a finite-horizon problem formulation and further decomposes the finite-horizon problem into a sequence of one-step sub-problems. Second, RHPG solves each sub-problem recursively using model-free PG methods. To accommodate the inevitable computational errors in solving these sub-problems, we establish the generalized principle of optimality that bounds the accumulated bias by controlling the inaccuracies in solving each sub-problem. We characterize the convergence and sample complexity of RHPG in §3-D and emphasize that the RHPG algorithm does not require knowing a stabilizing initial control policy a priori.

1-A Literature Review

We mainly compare with [3, 4] and [15, 8, 16], where [3, 4] are the foundational work in applying policy optimization to LQR and [15, 8, 16] remove the assumption on an initial stabilizing point in LQR by adding a discount factor to the objective function as an extra parameter.

In contrast to [3, 4] that parametrizes LQR as a single nonconvex (constrained) optimization problem over the policy space, we provide a new parametrization and perspective to learning LQR control. The critical difference between the RHPG algorithm and [3, 4] is that RHPG incorporates existing theories into the design of model-free learning algorithms, which is more than just exploiting them for the convergence analysis. We view this work as an initial step toward the goal of “designing control-specific learning algorithms with performance guarantees” rather than “analyzing existing learning algorithms for control”. This line of research is motivated by the observation that viewing the dynamical system as a black box and directly searching in the policy space leads to a deficit of most, if not all, favorable landscape properties beyond LQR [10, 9]. This implies that the excellent properties in LQR following the parametrization in [3, 4], such as coercivity and gradient domination, are rare, problem-dependent, and hard to generalize. However, when the model information is known, existing theories have provided extremely efficient solutions (e.g., Riccati recursions) to these problems that seemed unsolvable in the PO paradigm. Hence, our rationale is that control settings inherently have more structures (that are problem-agnostic) than a black-box system. One should always exploit these structures in developing learning-based algorithms with performance guarantees. In our work and [14, 17], we have identified causality and the dynamic programming principle as the fundamental properties in all control settings and exploited them in the model-free learning paradigm. As demonstrated in our work and [14, 17], the RHPG framework efficiently solves LQR and the seemed-to-be-unsolvable Kalman filtering problem in the PO fashion, which serves as a fundamental benchmark in output feedback control. As a by-product, the RHPG framework also removes all assumptions inappropriate in model-free learning settings and is more consistent with existing theories in control and estimation. Our work and [14, 17] together lead to a promising path toward the theoretical foundation of model-free learning in partially observable settings and nonlinear control through the lens of RHPG.

In comparison with [15, 8, 16], we note that the $\gamma$-discounted LQR problems therein are equivalent to the standard non-discounted LQR with system matrices being $\sqrt{\gamma}A$ and $\sqrt{\gamma}B$. Then, for any $\gamma\in(0,1)$, the set of stabilizing policies with system matrices being $\sqrt{\gamma}A$ and $\sqrt{\gamma}B$ is strictly larger than that of the un-discounted case. Hence, when $\gamma$ is sufficiently small, one can initialize the PG algorithm with an arbitrary control policy. This removes the requirement of knowing a stabilizing policy in advance, but it comes with the price of solving multiple LQRs instead of only one. Moreover, the criterion for increasing $\gamma$ is more complex than our selection rule for $N$ in RHPG. Besides removing the assumption on the initial stabilizing point, it is more important to determine if the results/insights can be further generalized/extended to more complicated control problems, where the critical difference between our work and [15, 8, 16] appears. In [15, 8, 16], the landscape of discounted LQR is identical to those in [3, 4] for un-discounted LQR, with essentially scaled versions of system parameters. Thus, the same difficulties reported in [10, 9] will appear when considering the output-feedback setting with an additional discount factor. In contrast, RHPG can be directly extended to solve output-feedback problems [17].

1-B Notations

We use $\|X\|$, $\kappa_{X}$, and $\rho(X)$ to denote, respectively, the spectral norm, condition number, and the spectral radius of a square matrix $X$. If $X$ is symmetric, we use $X>0$ and $X\geq 0$ to denote that $X$ is positive definite (pd) and positive semi-definite (psd), respectively. For a pd matrix $W$ of appropriate dimensions, we define the $W$-induced norm of $X$ as

[TABLE]

2 Preliminaries

2-A Infinite-Horizon LQR

Consider the discrete-time linear dynamical system111For extensions to stochastic LQR with i.i.d. additive noises, as well as the setting with an arbitrary (deterministic) initial state, see §5.

[TABLE]

where $x_{t}\in\mathbb{R}^{n}$ is the state; $u_{t}\in\mathbb{R}^{m}$ is the control input; $A\in\mathbb{R}^{n\times n}$ and $B\in\mathbb{R}^{n\times m}$ are system matrices unknown to the control designer; and the initial state $x_{0}\in\mathbb{R}^{n}$ is sampled from a zero-mean distribution $\mathcal{D}$ that satisfies $\operatorname*{\rm Cov}(x_{0})=\Sigma_{0}>0$. The goal in the LQR problem is to obtain the optimal controller $u_{t}=\phi_{t}(x_{t})$ that minimizes the cost

[TABLE]

where $Q>0$ and $R>0$ are symmetric pd weightings chosen by the control designer. For the LQR problem as posed to admit a solution, we require $(A,B)$ to be stabilizable. Note that here $Q>0$ implies the observability of $(A,Q^{1/2})$. Then, the unique optimal LQR controller is linear state-feedback, i.e., $u^{*}_{t}=-K^{*}x_{t}$, and $K^{*}\in\mathbb{R}^{m\times n}$, which with a slight abuse of terminology we will call optimal control policy, can be computed by

[TABLE]

where $P^{*}$ is the unique positive definite (pd) solution to the algebraic Riccati equation (ARE)

[TABLE]

Moreover, the optimal control policy $K^{*}$ is guaranteed to be stabilizing, i.e., $\rho(A-BK^{*})<1$. Therefore, we can parametrize LQR as an optimization problem over the policy space $\mathbb{R}^{m\times n}$, subject to the stability condition [3]:

[TABLE]

Theoretical properties of model-free (zeroth-order) PG methods in solving (2.5) have been well understood [3, 18, 4]. In particular, the objective function (2.5), even though nonconvex, is coercive and (globally) gradient dominated [18]. Hence, if an initial control policy $K_{0}\in\mathcal{K}$ is known a priori, then any descent direction of the objective value (e.g., vanilla PG) suffices to ensure that all the iterates will remain in the interior of $\mathcal{K}$ while quickly converging toward the unique stationary point. Removing the assumption on $K_{0}$ (that an initial stabilizing policy can readily be found) has remained an active research topic [15, 8, 16].

2-B Finite-Horizon LQR

The finite-$N$-horizon version of the LQR problem is also described by the system dynamics (2.1), but with the objective function summing up only up to time $t=N$. Similar to (2.5), we can parametrize the finite-horizon LQR problem as \min_{\{K_{t}\}}J\big{(}\{K_{t}\}\big{)}, where

[TABLE]

and $Q_{N}$ is a symmetric psd terminal-state weighting to be chosen. The unique optimal control policy in the finite-horizon LQR is time-varying and can be computed by

[TABLE]

where $P^{*}_{t}$, for all $t\in\{0,\cdots,N-1\}$, are generated by the Riccati difference equation (RDE) starting with $P^{*}_{N}=Q_{N}$:

[TABLE]

Theoretical properties of zeroth-order PG methods in addressing (2.7) have been studied in [7, 6]. Compared to the infinite-horizon setting (2.5), the finite-horizon LQR problem (2.7) is also a nonconvex and gradient-dominated problem, but it does not naturally require the stability condition (2.6).

3 Receding-Horizon Policy Gradient

3-A LQR with Dynamic Programming

It is well known that the solution of the RDE (2.9) converges monotonically to the stabilizing solution of the ARE (2.4) exponentially [19]. It then readily follows that the sequence of time-varying LQR policies (2.8), denoted as $\{K_{t}\}_{t\in\{N-1,\cdots,0\}}$, converges monotonically to the time-invariant LQR policy $K^{*}$ as $N\to\infty$. Now, we formally present this non-asymptotic convergence result.

Theorem 3.1

Let $A_{K}^{*}:=A-BK^{*}$, use $\|\cdot\|_{*}$ to denote the $P^{*}$-induced norm, and define

[TABLE]

where it holds that $\|A_{K}^{*}\|_{*}<1$. Then, for all $N\geq N_{0}$, the control policy $K^{*}_{0}$ computed by (2.8) satisfies $\|K^{*}_{0}-K^{*}\|\leq\epsilon$ for any $\epsilon>0$.

We provide the proof of Theorem 3.1 in §-A. Theorem 3.1 demonstrates that if selecting $N\sim\mathcal{O}(\log(\epsilon^{-1}))$, then solving the finite-horizon LQR will result in a policy $K^{*}_{0}$ that is $\epsilon$-close to $K^{*}$, for any $\epsilon>0$. Furthermore, if one chooses a small enough $\epsilon$ such that an $\epsilon$-ball centered at $K^{*}$ lies entirely in $\mathcal{K}$, then this condition on $\epsilon$ constitutes a sufficient condition for $K^{*}_{0}$ to be stabilizing, i.e., $K^{*}_{0}\in\mathcal{K}$.

3-B Algorithm Design

We propose the RHPG algorithm (cf., Algorithm 1), which first selects $N$ by Theorem 3.1, and then sequentially decomposes the finite-$N$-horizon LQR problem backward in time. In particular, for every iteration indexed by $h\in\{N-1,\cdots,0\}$, the RHPG algorithm solves an LQR problem from $t=h$ to $t=N$, where we only optimize for the current policy $K_{h}$ and fix all the policies $\{K_{t}\}$ for $t\in\{h+1,\cdots,N-1\}$ to be the convergent solutions generated from earlier iterations. Concretely, for every $h$, the RHPG algorithm solves the following quadratic program in $K_{h}$:

[TABLE]

Due to the quadratic optimization landscape of (3.2) in $K_{h}$ for every $h$, applying any PG method with an arbitrary finite initial point (e.g., zero) would lead to convergence to the globally optimal solution of (3.2).

3-C Bias of Model-Free Receding-Horizon Control

The RHPG algorithm builds on Bellman’s principle of optimality, which requires solving each iteration to the exact optimal solution. However, PG methods can only return an $\epsilon$-accurate solution after a finite number of steps. To generalize Bellman’s principle of optimality, we analyze how computational errors accumulate and propagate in the (backward) dynamic programming process. In the theorem below, we show that if one solves every iteration of the RHPG algorithm to the $\mathcal{O}(\epsilon)$-neighborhood of the unique optimum, then the RHPG algorithm will output a policy that is $\epsilon$-close to the infinite-horizon LQR policy $K^{*}$.

Theorem 3.2

Choose $N$ according to Theorem 3.1 and assume that one can compute, for all $h\in\{N-1,\cdots,0\}$ and some $\epsilon>0$, a policy $\widetilde{K}_{h}$ that satisfies

[TABLE]

where $\widetilde{K}_{h}^{*}$ is the optimum of the LQR from $h$ to $N$, after absorbing errors in all previous iterations of Algorithm 1. Then, the RHPG algorithm outputs a control policy $\widetilde{K}_{0}$ that satisfies \big{\|}\widetilde{K}_{0}-K^{*}\big{\|}\leq\epsilon. Furthermore, if $\epsilon$ is sufficiently small such that $\epsilon<\frac{1-\|A-BK^{*}\|_{*}}{\|B\|}$, then $\widetilde{K}_{0}$ is stabilizing.

We illustrate Theorem 3.2 in Figure 1 and defer its proof to §-B. Additionally, we discuss the implications of Theorem 3.2 in the following remark.

Remark 3.3

Theorem 3.2 ensures that if each RHPG iteration satisfies a certain error tolerance level, then the RHPG output $\widetilde{K}_{0}$ will reach an $\epsilon$-neighborhood of $K^{*}$. Notably, upon selecting a sufficiently small $\epsilon$, such that the $\epsilon$-ball around $K^{*}$ lies entirely within $\mathcal{K}$, we can guarantee closed-loop stability of $\widetilde{K}_{0}$ solely based on its near-optimality. This approach differs significantly from existing PG for LQR literature (e.g., [3, 4]), which requires a stabilizing initial policy and focus on preserving the closed-loop stability during learning. Conversely, RHPG starts with an arbitrary initial point, potentially very far from the stabilizing region (see our numerical experiments in §4), yet still converges globally to $K^{*}$. The closed-loop stability certificate then comes for free, given near-optimality.

Now, it remains to establish the sample complexity for the convergence of (zeroth-order) PG methods in every iteration of the algorithm, which is done next.

3-D PG Update and Sample Complexity

We analyze here the sample complexity of the zeroth-order PG update in solving each iteration of the RHPG algorithm. Specifically, the zeroth-order PG update is defined as

[TABLE]

where $\eta_{h}>0$ is the stepsize to be determined and $\widetilde{\nabla}J_{h}(K_{h,i})$ is the estimated PG sampled from a (two-point) zeroth-order oracle. We formally present the sample complexity result in the following proposition.

Proposition 3.4

For all $h\in\{0,\cdots,N-1\}$, choose a constant smoothing radius $r_{h}\sim\mathcal{O}(\epsilon)$ and a constant stepsize $\eta_{h}\sim\mathcal{O}(\epsilon^{2})$. Then, the zeroth-order PG update (3.3) converges after $T_{h}\sim\mathcal{O}(\frac{1}{\epsilon^{2}}\log(\frac{1}{\delta\epsilon^{2}}))$ iterations in the sense that \big{\|}K_{h,T_{h}}-\widetilde{K}_{h}^{*}\big{\|}\leq\epsilon with a probability of at least $1-\delta$.

For completeness, we provide a supplementary proof of Proposition 3.4 in §-C, which mostly follows existing results in the literature [4]. Combining Theorem 3.2 with Proposition 3.4, we conclude that if we spend $\widetilde{\mathcal{O}}(\epsilon^{-2}\log(\delta^{-1}))$ iterations in solving every subproblem to $\mathcal{O}(\epsilon)$-accuracy with a probability of $1-\delta$, for all $h\in\{0,\cdots,N-1\}$, then the RHPG algorithm will output a $\widetilde{K}_{0}$ that satisfies $\|\widetilde{K}_{0}-K^{*}\|\leq\epsilon$ with a probability of at least $1-N\delta$. By (3.1), this implies that the total iteration complexity of RHPG is also $\widetilde{\mathcal{O}}(\epsilon^{-2}\log(\delta^{-1}))$ with the dependence on various system parameters being polynomial.

We discuss the tradeoffs in selecting $N$ to balance minimizing finite-to-infinite error and minimizing errors through the backward propagation in §5-A. To compare our sample complexity bound with the sharpest result in the literature [4], our dependence on $\epsilon$ matches that of [4]222Note that the $\widetilde{\mathcal{O}}(\epsilon^{-1})$ sample complexity presented in [4] is for the convergence in objective value (e.g., $f(K)-f(K^{*})\leq\epsilon$), and is equivalent to an $\widetilde{\mathcal{O}}(\epsilon^{-2})\cdot\mathcal{O}(\texttt{poly}(\text{system parameters}))$ sample complexity for the convergence in policy (i.e., $\|K-K^{*}\|\leq\epsilon$).. Our sample complexity and that of [4] have polynomial dependence on system parameters. However, it is not clear how to compare the polynomial dependencies between our bounds and those of [4], and these polynomial factors might affect the overall computational efficiency of both algorithms in a substantial way. We leave this comparison as an important future research topic.

4 Numerical Experiments

We verify our theories on a scalar linear system studied in [4], where $A=5$, $B=0.33$, $Q=R=1$, and the unique optimal LQR policy is $K^{*}=14.5482$. For the PG method in [3, 4] to converge in this simple setting, one has to initialize with a policy $K_{0}$ that satisfies $K_{0}\in\underline{\mathcal{K}}:=\{K\mid 12.12<K\pm r<18.18\}$, which is necessary to prevent the zeroth-order oracle with a smoothing radius of $r$ from perturbing $K_{0}$ outside of the stabilizing region $\mathcal{K}$ during the first iteration of the PG update. In contrast, we initialize the PG updates in Algorithm 1 with a zero policy $K_{h}=0$, set $Q_{N}=3$, and choose $N=\texttt{ceil}(\log(\epsilon^{-1}))$ according to (3.1). Furthermore, we choose $r_{h}=\sqrt{\epsilon}$, select a constant stepsize in each iteration of the RHPG algorithm, and run the algorithm to solve the LQR problem under six different $\epsilon$, namely $\epsilon\in\{10^{-3},3.16\times 10^{-3},10^{-2},3.16\times 10^{-2},10^{-1},3.16\times 10^{-1}\}$. We apply the zeroth-order PG update in solving every subproblem to \big{\|}\widetilde{K}_{h}-\widetilde{K}^{*}_{h}\big{\|}\leq\epsilon. As shown in Figure 2, the empirical observation of the iteration complexity of RHPG (right) for the convergence in policy (left) is around $\mathcal{O}(\epsilon^{-2})$ under varying $\epsilon$, which corroborates our theoretical findings.

5 Discussions

In this section, we provide two discussion paragraphs. The first paragraph discusses the tradeoffs in selecting the problem horizon $N$. The second paragraph covers extensions of our results to the stochastic LQR setting that incorporates a zero-mean, independent stochastic disturbance in the dynamics and, additionally, another setting where the initial state $x_{0}$ could be arbitrary.

5-A Tradeoffs in Selecting the Problem Horizon $N$

There exist two potential usages of the RHPG algorithm in practice. In most real-world scenarios, the user wants to address the finite-horizon LQR problem with time-varying system parameters and does not necessarily care whether the resulting control policy is close to the infinite-horizon solution or not. Model-predictive or receding-horizon control is the most prevalent choice when detailed model information is available. The RHPG solution, on the other hand, extends receding-horizon control to the setting with an unknown model, and the user could specify a desired horizon $N$.

In more theory-oriented cases, the user wants to address the infinite-horizon LQR problem with time-invariant system parameters. Then the user should choose the problem horizon $N$ carefully, balancing between i) reducing the finite-to-infinite horizon error by increasing $N$, and ii) potentially creating more computational error in the approximate dynamic programming when $N$ is large. The RHPG framework performs these two tasks sequentially: first dealing with the finite-to-infinite horizon error and fixing an $N$, and then addressing the finite-$N$-horizon problem step-by-step in time. The main rationale comes from the different rates of $N$ contributing to the two errors. In the first step (cf., Theorem 3.1), the convergence of the finite-horizon solution toward the infinite-horizon solution is exponential, meaning that it suffices to choose a problem horizon $N=\mathcal{O}(\log(\epsilon^{-1}))$ for an accurate approximation (the dependence to other system parameters are also logarithmic). After deciding on an $N$ based on the desired $\epsilon$, the user solves the finite-$N$-horizon problem iteratively, where the subproblem in each iteration is strongly convex and smooth. In the sample complexity analysis, we have considered the choice of $N$ since one would require a higher solution accuracy in each iteration if $N$ is large (cf., §-B). Moreover, $N$ also appears in the probability term of $1-N\delta$ due to the utilization of Boole’s inequality. In solving each sub-problem, the complexity in the failure probability $\delta$ is $\mathcal{O}(\log(\delta^{-1}))$. Replacing $\delta$ with $\delta/N$ yields the rate of $\mathcal{O}(\log(N/\delta))=\mathcal{O}(\log(\delta^{-1}))+\mathcal{O}(\log\log(\epsilon^{-1}))$. Hence, with a fixed $N$ chosen according to Theorem 3.1, the impact of $N$ on the total complexity is an additional $\mathcal{O}(\log\log(\epsilon^{-1}))$ factor and is therefore minor. In summary, the user in the infinite-horizon setting should prioritize the choice of $N$ toward reducing the finite-to-infinite horizon error according to Theorem 3.1 and then solve the finite-$N$-horizon problem.

5-B Extensions to Arbitrary $x_{0}$ and Stochastic LQR

We discuss here extensions of the RHPG framework to i) the setting where $x_{0}\in\mathbb{R}^{n}$ is an arbitrary (deterministic) vector that is unknown to the designer, and ii) stochastic LQR with $x_{0}\sim\mathcal{D}$ and the system dynamics being

[TABLE]

Note that in the stochastic LQR setting, the objective function (2.2) should be replaced with the time-average cost

[TABLE]

In both settings and under the stabilizability and detectability assumptions, a unique stationary state-feedback control policy exists and is identical to (2.3), which also stabilizes the closed-loop system. The difference is that in setting i), the LQR problem is deterministic, which allows implementing RHPG with a two-point zeroth-order oracle (cf., Algorithm 1). In contrast, setting ii) involves additive noises in the system dynamics, which necessitates using a one-point zeroth-order oracle, and thus, the gradient sampling will be noisier.

We note that the procedure of RHPG is identical whether solving deterministic or stochastic problems and with either finite or infinite problem horizons. The only difference is in the choice of the inner-loop oracles (two-point v.s. one-point). In both settings, the receding-horizon parametrization, the convergence of the Riccati equations in Theorem 3.1, the analysis of dynamic programming in Theorem 3.2, and the quadratic optimization landscape in each subproblem are the same as in the presentation in the main body of this paper. Simply combining Theorems 3.1-3.2 with the corresponding inner-loop convergence result (for quadratic minimization and as a replacement of Proposition 3.3) yields the overall complexity in these two extended settings.

6 Conclusion

We have revisited discrete-time LQR from the perspective of RHPG and provided a fine-grained sample complexity analysis for RHPG to learn a control policy that is stabilizing and $\epsilon$-close to the optimal LQR policy. Our result demonstrates the potential of RHPG in addressing various tasks in linear control and estimation with streamlined analyses.

-A Proof of Theorem 3.1

This proof is dual to the proof of Theorem 3.1 in [14]. We first identify one-to-one correspondences between system parameters in LQR and those in Kalman filtering [14]:

[TABLE]

We also identify direct correspondences between our $P_{t}$, $P^{*}$, $K_{t}$, and $K^{*}$ and [14]’s $\Sigma_{N-t}$, $\Sigma^{*}$, $L_{N-1-t}$, and $L^{*}$, respectively. Then, letting

[TABLE]

and following equations (A.2)-(A.3) of [14], we have

[TABLE]

where $\widetilde{P}_{t+1}^{1/2}$ denotes the unique positive semi-definite (psd) square root of the psd matrix $\widetilde{P}_{t+1}$, $0<\mu_{t}\leq 1$ for all $t$, and $\overline{A}$ satisfies $\rho(\overline{A})<1$. We now use $\|\cdot\|_{*}$ to represent the $P^{*}$-induced matrix norm and invoke Theorem 14.4.1 of [19], where our $\widetilde{P}_{t}$, $\overline{A}^{\top}$ and $P^{*}$ correspond to $P_{i}-P^{*}$, $F_{p}$ and $W$ in [19], respectively. By Theorem 14.4.1 of [19] and (-A.1), we obtain $\|\overline{A}\|_{*}<1$ and given that $\mu_{t}\leq 1$,

[TABLE]

Therefore, the convergence is exponential such that $\|\widetilde{P}_{t}\|_{*}\leq\|\overline{A}\|_{*}^{2(N-t)}\cdot\|\widetilde{P}_{N}\|_{*}$. As a result, the convergence of $\widetilde{P}_{t}$ to [math] in spectral norm can be characterized as

[TABLE]

where we have used $\kappa_{X}$ to denote the condition number of $X$. That is, to ensure $\|\widetilde{P}_{1}\|\leq\epsilon$, it suffices to require

[TABLE]

Lastly, we show that the (monotonic) convergence of $K^{*}_{t}$ to $K^{*}$ follows from the convergence of $P^{*}_{t}$ to $P^{*}$. Similar to (A.5) of [14], this can be verified through:

[TABLE]

Hence, we have $\|K^{*}_{t}-K^{*}\|\leq\frac{\|\overline{A}\|\cdot\|B\|}{\lambda_{\min}(R)}\cdot\|P^{*}_{t+1}-P^{*}\|$ and

[TABLE]

Substituting $\epsilon$ in (-A.2) with $\frac{\epsilon\cdot\lambda_{\min}(R)}{\|\overline{A}\|\cdot\|B\|}$ completes the proof.

-B Proof of Theorem 3.2

This proof is dual to the proof of Theorem 3.3 in [14]. First, according to Theorem 3.1, we select

[TABLE]

where $A_{K}^{*}:=A-BK^{*}$. This ensures that $K^{*}_{0}$ is stabilizing and $\|K^{*}_{0}-K^{*}\|\leq\epsilon/2$. Then, it remains to show that the output $\widetilde{K}_{0}$ satisfies $\|\widetilde{K}_{0}-K^{*}_{0}\|\leq\epsilon/2$.

Recall that the RDE (2.9) is a backward iteration starting with $P^{*}_{N}=Q_{N}\geq 0$, and can also be represented as:

[TABLE]

Moreover, for any $K_{t}$, we introduce the Lyapunov equation:

[TABLE]

Furthermore, for clarity of the proof, we define/recall:

[TABLE]

We argue that $\|\widetilde{K}_{0}-K^{*}_{0}\|\leq\epsilon/2$ can be achieved by carefully controlling $\delta_{t}$ for all $t$. At $t=0$, it holds that

[TABLE]

where substituting $K^{*}_{t}$ and $K^{*}$ in (-A.3), respectively, with $\widetilde{K}^{*}_{0}$ and $K^{*}_{0}$ leads to

[TABLE]

Hence, the error size $\|\widetilde{K}^{*}_{0}-K^{*}_{0}\|$ could be bounded by

[TABLE]

Define the helper constants

[TABLE]

Next, we require $\|\delta_{0}\|\leq\epsilon/4$ and $\|\widetilde{K}^{*}_{0}-K^{*}_{0}\|\leq\epsilon/4$ to fulfill $\|\widetilde{K}_{0}-K^{*}_{0}\|\leq\epsilon/2$. We select a fixed scalar $a>0$ that is independent of system parameters and $\epsilon$, and additionally require $\|P^{*}_{1}-\widetilde{P}_{1}\|\leq a$ to upper-bound the pd solutions of (-B.7). Then, by (-B.8), in order to fulfill $\|\widetilde{K}^{*}_{0}-K^{*}_{0}\|\leq\epsilon/4$, it suffices to require

[TABLE]

Subsequently, by (-B.7), we have

[TABLE]

The first difference term on the RHS of (-B.10) is

[TABLE]

where applying (-B.8) in deriving (-B.12) and (-B.13) is due to the matrix inversion lemma. Moreover, the second term on the RHS of (-B.10) is

[TABLE]

where (-B.14) follows from completion of squares. Thus, combining (-B.10), (-B.11), and (-B.15) yields

[TABLE]

where (-B.16) is due to that $\widetilde{P}_{2}BR^{-1}B^{\top}\geq 0$ and thus $\|(I+\widetilde{P}_{2}BR^{-1}B^{\top})^{-1}\|\leq 1$. Now, we require

[TABLE]

where $C_{2}$ and $C_{3}$ are positive constants defined as333As the scalar $a>0$ increases, the constant $C_{3}$ grows correspondingly.

[TABLE]

Then, conditions (-B.17) and (-B.18) are sufficient for (-B.9) (and thus for $\|\widetilde{K}_{0}-K^{*}_{0}\|\leq\epsilon/2$) to hold. Subsequently, we can propagate the required accuracies in (-B.17) and (-B.18) forward in time. Specifically, we iteratively apply the arguments in (-B.16) (i.e., by plugging quantities with subscript $t$ into the LHS of (-B.16) and plugging quantities with subscript $t+1$ into the RHS of (-B.16)) to obtain the result that if at all $t\in\{2,\cdots,N-1\}$, we require

[TABLE]

then (-B.17) holds true and therefore (-B.9) is satisfied.

We now compute the required accuracy for $\delta_{N-1}$. Note that $P^{*}_{N-1}=\widetilde{P}^{*}_{N-1}$ since no prior computational errors happened at $t=N$. By (-B.16), the distance between $P^{*}_{N-1}$ and $\widetilde{P}_{N-1}$ can be bounded as

[TABLE]

To fulfill the requirement (-B.19) for $t=N-1$, which is

[TABLE]

it suffices to let

[TABLE]

Finally, we analyze the worst-case complexity of RHPG by computing, at the most stringent case, the required size of $\|\delta_{t}\|$. When $C_{2}\leq 1$, the most stringent dependence of $\|\delta_{t}\|$ on $\epsilon$ happens at $t=0$, which is of the order $\mathcal{O}(\epsilon)$, and the dependences on system parameters are $\mathcal{O}(1)$. We then analyze the case where $C_{2}>1$, where the requirement on $\|\delta_{0}\|$ is still $\mathcal{O}(\epsilon)$. Note that in this case, $\|\delta_{N-1}\|\leq\|\delta_{t}\|$ for all $t\in\{1,\cdots,N-1\}$ and by (-B.20):

[TABLE]

Since we require $N$ to satisfy (-B.4), the dependence of $\|\delta_{N-1}\|$ on $\epsilon$ in (-B.21) becomes $\|\delta_{N-1}\|\sim\mathcal{O}(\epsilon^{\frac{3}{4}})$ with additional polynomial dependences on system parameters, but one can observe that the dependence on $\epsilon$ is still milder than the requirement for $\|\delta_{0}\|$. Therefore, it suffices to require error bound for all $t$ to be $\|\delta_{t}\|\sim\mathcal{O}(\epsilon)\mathcal{O}(\texttt{poly}(\text{system parameters)})$ to reach the $\epsilon$-neighborhood of the infinite-horizon LQR policy. Lastly, for $\widetilde{K}_{0}$ to be stabilizing, it suffices to select a sufficiently small $\epsilon$ such that the $\epsilon$-ball centered at the infinite-horizon LQR policy $K^{*}$ lies entirely in the set of stabilizing policies. A crude bound that satisfies this requirement is

[TABLE]

This completes the proof.

-C Proof of Proposition 3.4

Recall that for all $h$, the objective function $J_{h}$ is $L_{h}$-smooth and $\alpha_{h}$-strongly-convex. Define $\varsigma_{h}:=\frac{\epsilon^{2}\alpha_{h}}{2}$ and $\varsigma:=\min_{h}\varsigma_{h}>0$. We argue that if with a probability of at least $1-\delta$, it holds that

[TABLE]

then $\|K_{h,T_{h}}-\widetilde{K}^{*}_{h}\|\leq\epsilon$ also holds with a probability of at least $1-\delta$ and the proof of Proposition 3.4 is complete. This is due to the $\alpha_{h}$-strong convexity and $\nabla J_{h}(\widetilde{K}^{*}_{h})=0$. Thus,

[TABLE]

As a result, we will focus on proving (-C.22) with a high probability of at least $1-\delta$. First, define the cost difference $\Delta_{t}:=J_{h}(K_{h,t})-J_{h}(\widetilde{K}_{h}^{*})$ and the stopping time $\tau:=\min\{t\mid\Delta_{t}>10\delta^{-1}\Delta_{0}\}$. Let $\mathbb{E}^{t}[\cdot]$ denote the expectation conditioned on all the randomness up to $t$. Then, we state the following helper lemma and defer its proof to §-D.

Lemma 0.1

For all $h$, choose the parameters of Algorithm 1 according to

[TABLE]

Then, for all $t$, it holds that

[TABLE]

where $\alpha_{h}$ and $L_{h}$ are the strong convexity and smoothness constants of $J_{h}$, respectively, and $G_{2}$ is a uniform constant to be introduced shortly.

Following the proof of Theorem 8 in [4, 20], we first consider the case of $\tau>T_{h}$. In this case, we can bound $\mathbb{E}^{t}[\Delta_{t+1}]$ using Lemma 0.1 directly. When $\tau\leq T_{h}$, it implies that $\mathbb{E}^{t}[\Delta_{t+1}]1_{\tau>t}=0$. We require $\eta_{h}\leq\frac{\varsigma\delta\alpha_{h}}{40L_{h}G_{2}}$ and show

[TABLE]

Setting $t+1=T_{h}$, it suffices to let $T_{h}=\frac{4}{\eta_{h}\alpha_{h}}\log(\frac{4\Delta_{0}}{\delta\varsigma})$ to ensure that

[TABLE]

Next, we prove that the event $\tau\leq T_{h}$ has a probability smaller than $\frac{\delta}{2}$. For all $t$, we define the stopping process as

[TABLE]

By Eq. (20)-(21) of [4], $Y_{t}$ is a super-martingale. Applying Doob’s maximal inequality yields

[TABLE]

Imposing the condition that $\varsigma\log(\frac{4\Delta_{0}}{\varsigma\delta})\leq 16\delta^{-1}\Delta_{0}$, we can prove that P\big{(}\max_{t=1,\cdots,T_{h}}Y_{t}\geq\frac{10\Delta_{0}}{\delta}\big{)}\leq\frac{\delta\mathbb{E}[Y_{0}]}{10\Delta_{0}}\leq\frac{\delta}{2}. We can now conclude that $\mathbb{E}[\Delta_{T_{h}}1_{\tau>T_{h}}]\leq\frac{\delta\varsigma}{2}$ and the event $\tau$ occurs after $T_{h}$ with probability at least $1-\frac{\delta}{2}$. As a result,

[TABLE]

where we have used Markov’s inequality. This verifies (-C.22), and thus $\|K_{h,T_{h}}-\widetilde{K}^{*}_{h}\|\leq\epsilon$ is satisfied with a probability of at least $1-\delta$. Lastly, we analyze the constant $G_{2}$ following Corollary 10 of [4] and [21], where

[TABLE]

By Corollary 10 of [4], it holds almost surely that $G_{2}\leq(mn)\lambda^{2}$, where $\lambda:=\max_{h}\lambda_{h}$ and $\lambda_{h}$ is the Lipschitz continuity constant of $J_{h}$ taken over the compact domain $\Phi_{h}$. In summary, for $\|K_{h,T_{h}}-\widetilde{K}^{*}_{h}\|\leq\epsilon$ to hold with a probability of at least $1-\delta$, we need to choose the (constant) algorithmic parameters according to $\eta_{h}\sim\mathcal{O}(\epsilon^{2})$, and $r_{h}\sim\mathcal{O}(\epsilon)$. Then, the iteration complexity for the convergence of the zeroth-order PG method is $T_{h}\sim\mathcal{O}(\frac{1}{\epsilon^{2}}\log(\frac{1}{\delta\epsilon^{2}}))$.

-D Proof of Lemma 0.1

The proof of this lemma mostly follows the steps in Section 4.1.1 of [4]. First, define the smoothed version of $J_{h}$ as $J_{h}^{r_{h}}(K_{h}):=\mathbb{E}[J_{h}(K_{h}+r_{h}U)]$, where the expectation is taken over $U$ that is uniformly drawn from the surface of a unit sphere. Then, we use $J(K_{h};x_{0})$ to denote an instantiation of the objective value $J_{h}(K_{h})$ given $x_{0}$, and define the two-point zeroth-order estimate of $\nabla J_{h}^{r_{h}}$ as

[TABLE]

where $mn$ is the dimension of the policy space. Subsequently, we invoke the $L_{h}$-smoothness property to derive

[TABLE]

where the inequalities are due to Lemma 14 of [4]. Moreover,

[TABLE]

Again by the $L_{h}$-smoothness property, we have

[TABLE]

Then, letting $\eta_{h}\in(0,\frac{1}{2L_{h}}]$, we can derive

[TABLE]

where the second inequality is due to that the $\alpha_{h}$ strong-convexity implies the $\alpha_{h}$ gradient domination property, the choice of stepsize $\eta_{h}\leq\frac{1}{2L_{h}}$, and $2ab\leq a^{2}+b^{2}$ for any $a,b$. Recall the choices of algorithmic parameters as follows:

[TABLE]

Then, using the bounds on algorithmic parameters and rearranging terms lead to

[TABLE]

which completes the proof.