D2C 2.0: Decoupled Data-Based Approach for Learning to Control Stochastic Nonlinear Systems via Model-Free ILQR

It is well-known that a general solution to stochastic optimal control problem is computationally intractable. Hence, most methods rely on approximated dynamic programming methods (See, for e.g., [1]). In case of an unknown system model, the problem is further complicated. This problem of learning to optimally control an unknown general nonlinear system has been well-studied in the stochastic adaptive control literature [2] [3]. A very similar problem is also addressed in the field of reinforcement learning (RL) except that it specifically aims to solve a discrete-time stochastic control process (Markov Decision Process), more often with an infinite horizon formulation [4]. Also, RL has attracted a great attention with its recent success, especially in its sub-domain of deep RL [5] [6] [7].
Most methods in the RL can be divided into model-based and model-free approaches, based on whether a model is constructed before solving for an optimal policy or not, respectively. Among them, model-based RL typically requires lesser number of samples and can be global or local model-based. Global models are typically approximated with Gaussian Processes (GPs) or Neural Networks, and hence can be inaccurate in some portions of its state space. Local models (such as in [8][9]) in contrast offer more sample efficient and accurate solutions. On the other hand, model-free RL, since has to bypass learning the system model requires substantially more number of training samples to train the policy [10]. While an overwhelming number of successful works in the recent past have relied on a complex parameterizations such as by using neural networks [7] [11][12] [13], they are challenged with non-determinism in reproducing the results and in benchmarking against existing methods [14][15]. Recent works such as by [16] showed that a simple linear parameterization can be very effective and can result in performances competitive to that of the former methods. The current work, aimed at deriving a simple but structured, sample-efficient and thereby a fast policy, falls under the local model-based regime, though we won’t explicitly build the system model but rather deal with its gradients in the neighborhood of a nominal trajectory.

In the prequel of this work, called D2C [17], instead of solving for the optimal closed-loop policy via dynamic programming as in typical RL techniques, our approach searches for a computationally tractable policy that is nonetheless near-optimal to second order in a noise parameter $\epsilon$. This is done via decoupling of the open and the closed loops, wherein we first solve for an optimal open-loop trajectory, and then design an optimal linear feedback law to stabilize the linear perturbations about the nominal trajectory. The success of this approach was previously demonstrated in [18][19] for model-based problems. This way of decoupling results in a compact parameterization of the feedback law in terms of an open loop control sequence and a linear feedback policy, and thereby, constitutes a novel approach to the RL problem via the decoupled search. The contributions of this sequel paper are as follows:
i) We incorporate a data-efficient, and hence, significantly faster way to compute the gradients in the ILQR framework in a model free fashion. This, in turn, is incorporated in our D2C framework to solve for the optimal open-loop trajectory. This is in contrast to the ‘direct’ gradient-descent based approach introduced in the earlier version of our algorithm (D2C 1.0).
ii) We make a comparison with a state-of-the-art reinforcement learning algorithm - DDPG [7] as well as D2C 1.0, and show the significant merit of the current method, termed D2C 2.0, in terms of its computational efficiency.
iii) D2C 2.0 also constitutes a novel approach to solving complex robotic motion planning problems under uncertainty in a data efficient model-free fashion.

Differential Dynamic Programming (DDP) [20] is a class of iterative algorithms for trajectory optimization. One salient aspect of DDP is that it exhibits quadratic convergence for a general non-linear optimal control problem and can be solved in a single step for linear systems. In spite of its origins in the 1960s, it gained popularity only in the last decade due to the success of the modified algorithm called ILQR [21]. Though DDP (theoretically) guarantees a quadratic convergence, some of the terms in it involve computing second order derivatives of the system dynamical models. Since the dynamical models of most systems are multivariable vector-valued functions and their second order derivatives being third order tensors, DDP in its original form was not effective for practical implementation. ILQR [22] [23] dropped these terms and introduced regularization schemes [24][25] to handle the absence of these computationally expensive terms. This resulted in a faster and stable convergence to the optimal solution.

One advantage that ILQR has in terms of solving it in a model-free manner is that this algorithm is explicit in terms of system model and their first order derivatives. Earlier works such as by [25] and [26] employed finite differencing in computing the Jacobians in ILQR. Typically, a forward Euler scheme is chosen to independently determine each element (i.e, gradient) in a Jacobian matrix. [9] presented ILQR-LD by learning the dynamics through Locally Weighted Projection Regression (LWPR) approximation of the model first and then obtaining a corresponding analytical Jacobian. Their work also demonstrated a better computational efficiency over finite differences.

Spall introduced the ‘Simultaneous Perturbation Stochastic Approximation’ (SPSA) [27] method that evaluates the function only twice to calculate the Jacobian of the cost function. In this paper, a similar formulation is derived to compute the Jacobians of the system model online, in a model-free fashion, through the least squares method in a central-difference scheme.

The rest of the paper is organized as follows. Section 3 introduces the terminology and provides the problem statement. In section 4, we present a near-optimal decoupling principle that forms the basis of our algorithm. Section 5 describes our algorithm - D2C 2.0 in detail. In Sections 6 and 7, we provide results and a discussion of the application of the D2C 2.0 to several benchmark examples, along with comparisons to DDPG and the D2C technique.

Let ${\bf x_{t}}$ $\in$ $\mathcal{X}\in$ $\mathbb{R}^{n_{x}}$ represent the state of a system and ${\bf u_{t}}\in\mathcal{X}\in\mathbb{R}^{n_{u}}$ be the control signal at time $t$ respectively. Let $f:\mathbb{R}^{n_{x}}\times\mathbb{R}^{n_{u}}\rightarrow\mathbb{R}^{n_{x}}$ denote the deterministic state transition model of the system. Let $\pi_{t}:\mathbb{R}^{n_{x}}\rightarrow\mathbb{R}^{n_{u}}$ be the policy to be applied on the system. Let $c_{t}:\mathbb{R}^{n_{x}}\times\mathbb{R}^{n_{u}}\rightarrow\mathbb{R}$ be the time-indexed incremental cost function for all $t\in\{0,1,\cdots,N-1\}$ and $c_{N}:\mathbb{R}^{n_{x}}\rightarrow\mathbb{R}$ be the terminal cost function. Let $J_{t}^{\pi_{t}}:\mathbb{R}^{n_{x}}\rightarrow\mathbb{R}$ represent the cost-to-go function under the policy $\pi_{t}$ and $Q_{t}^{\pi_{t}}:\mathbb{R}^{n_{x}}\times\mathbb{R}^{n_{u}}\rightarrow\mathbb{R}$ be the corresponding action-value function, both at time $t$. Finally, let us denote the derivatives of a variable at time $t$ w.r.t. the state ${\bf x_{t}}$ or the control ${\bf u_{t}}$ by placing respective variables in the subscript. For example, $J_{\bf x_{t}}=\frac{\partial J_{t}}{\partial{\bf x_{t}}}|_{\bf\bar{x}_{t}}$. Unless stated otherwise, they are evaluated at the nominal state or control at time $t$.

Now, let us define the notion of process noise in the system through an important parameter called ‘noise scaling factor’, denoted by $\epsilon$. We assume that the noise is additive and white, zero-mean Gaussian. Hence, the state evolution equation of our stochastic system is represented as ${\bf x_{t+1}}=f({\bf x_{t}},{\bf u_{t}})+\epsilon\omega_{t}$, where $\epsilon<1$ and $\omega_{t}$ is zero-mean Gaussian distributed.
Let us define the notion of a ‘nominal’ trajectory, given a policy $\pi_{t}(.)$, as follows - the nominal control actions of the policy are the control actions output by the policy when all the noise variables $\omega_{t}=0$. Let us represent the nominal trajectory variables with bars over them $i.e,$ ${\bf\bar{x}_{t}}$ and ${\bf\bar{u}_{t}}$, both at time $t$. So, the nominal trajectory is given by, $\mathbb{T}_{nom}=\{\bf\bar{x}_{0:N},\bar{u}_{0:N-1}\}$, where $N$ is the time horizon, and $\bar{\mathbf{x}}_{t+1}=f(\bar{\mathbf{x}_{t}},\bar{\mathbf{u}}_{t})$, where $\bar{\mathbf{u}}_{t}=\pi_{t}(\bar{\mathbf{x}}_{t})$.

Using the above definitions, the problem of stochastic optimal control that we consider is to find the set of time-varying optimal control policies $\pi=\{\pi_{1},\pi_{2},\cdots,\pi_{N-1}\}$ (subscripts denote the time-indices) that minimize the finite horizon total expected cumulative cost, $i.e,$

such that the system model i.e, ${\bf x_{t+1}}=f({\bf x_{t}},\pi_{t}({\bf x_{t}}))+\epsilon\omega_{t}$, relevant constraints and boundary conditions are satisfied. Since the current approach deals with model-free problem, the system model and its constraints are assumed to be implicitly satisfied in the data that the system/agent gathers from its environment (through a simulator or a real world experiment).

Note : The results provided in this particular section are restated from our prequel paper on D2C primarily for the reader’s convenience and for the sake of completeness.
In this section, we describe a near-optimal decoupling principle which provides the theoretical foundation for our D2C technique.
Assuming that $f(.)$ and $\pi_{t}(.)$ are sufficiently smooth, let us describe the linearized dynamics about the nominal trajectory $\mathbb{T}_{nom}$ corresponding to the policy $\pi_{t}(.)$ as follows. The state and control perturbations, given by $\delta{\bf x_{t}}={\bf x_{t}}-{\bf\bar{x}_{t}}$ and $\delta{\bf u_{t}}={\bf u_{t}}-{\bf\bar{u}_{t}}$ respectively, evolve as follows:

where $A_{t}=\frac{\partial f}{\partial{\bf x}}|_{{(\bf\bar{x}_{t}},{\bf\bar{u}_{t}})}$, $B_{t}=\frac{\partial f}{\partial{\bf u}}|_{({\bf\bar{x}_{t}},{\bf\bar{u}_{t}})}$, $K_{t}=\frac{\partial\pi({\bf x_{t}})}{\partial x}|_{\bf\bar{x}_{t}}$ and $S_{t}(.)$, $\tilde{S}_{t}(.)$, $R_{t}(.)$ are higher order terms in respective Taylor’s expansions. By substituting (3) in (2), we have

where $\bar{S}_{t}(\delta{\bf x_{t}})=B_{t}\tilde{S}_{t}(\delta{\bf x_{t}})+S_{t}(\delta{\bf x_{t}})+R_{t}(K_{t}\delta{\bf x_{t}}+\tilde{S}_{t}(\delta{\bf x_{t}}))$. Now, similarly expanding the cost function $c_{t}({\bf x_{t}},{\bf u_{t}})$ under the policy $\pi_{t}({\bf x_{t}})$ $\forall t\in\{0,1,\cdots,N-1\}$ about the nominal trajectory $\mathbb{T}_{nom}$ as $c_{t}({\bf x_{t}},\pi_{t}({\bf x_{t}}))=\bar{c}_{t}+\bar{C}_{t}\delta{\bf x_{t}}+\bar{H}_{t}(\delta{\bf x_{t}})$, where $\bar{c}_{t}=c_{t}({\bf\bar{x}_{t}},\pi_{t}({\bf\bar{x}_{t}}))$ and $\bar{C}_{t}=\frac{\partial c_{t}({\bf x_{t}},\pi_{t}({\bf x_{t}}))}{\partial{\bf x_{t}}}|_{\bf\bar{x}_{t}}$ results in

Lemma 1: The state perturbation equation (4) can be equivalently separated into $\delta{\bf x_{t}}=\delta{\bf x_{t}^{l}}+\bar{S}_{t}$ and $\delta{\bf x_{t+1}^{l}}=\bar{A}_{t}\delta{\bf x_{t}^{l}}+\epsilon\omega_{t}$, and $\bar{S}_{t}$ is $O(\epsilon^{2})$.
Proof: Proof is provided in the supplementary file.

From (5) and lemma-1, we obtain

Proposition 1: Given a feedback policy $\pi$ satisfying the smoothness requirements,

Proof: Proof is provided in the supplementary file. Noting that ${\bf u_{t}}=\pi_{t}({\bf x_{t}})={\bf\bar{u}_{t}}+K_{t}\delta{\bf x_{t}}+\tilde{S}_{t}(\delta{\bf x_{t}})$, the above proposition entails the following observations:
i) Observation 1: The cost-to-go along the nominal trajectory $\mathbb{T}_{nom}=\{{\bf\bar{x}_{0:N}},{\bf\bar{u}_{0:N-1}}\}$ given by $\bar{J}^{\pi}$ is within $O(\epsilon^{2})$ of the expected cost-to-go $\tilde{J}^{\pi}$ of the policy $\pi$.
ii) Observation 2: Given the nominal control sequence ${\bf\bar{u}_{0:N-1}}$, the variance of the cost-to-go is overwhelmingly determined by the linear feedback law $K_{t}\delta{\bf x_{t}}$ $i.e,$ by the variance of the linear perturbation of the cost-to-go, $\delta J_{1}^{\pi}$, under the linear closed loop system dynamics $\delta{\bf x_{t+1}}=\bar{A}_{t}\delta{\bf x_{t}^{l}}+\epsilon\omega_{t}$.

Proposition 1 and the above observations suggest that an optimal open-loop control sequence with a suitable linear feedback law (for the perturbed linear system) wrapped around it is approximately optimal. The following subsection summarizes the problems to be solved for each of them.

This section presents our decoupled data-based control (D2C-2.0) algorithm. There are two main steps involved in solving for the requisite policy:

1. 1.

   Design an open-loop nominal trajectory by solving the optimization problem via a model-free ILQR based approach $i.e,$ without explicitly making use of an analytic model of the system. This is described in subsection 5.1. This is the primary difference from the prequel D2C where we used a first order gradient based technique to solve the optimal open loop control problem.

2. 2.

   Determine the parameters of the time-varying linear perturbation system about its optimal nominal trajectory and design an LQR controller based on the system identification. This step is detailed in subsection 5.2.

This section reports the implementation details of D2C-2.0 in simulation on various systems and their training results. We train the systems with varying complexities in terms of the dimension of the state and the action spaces. The physical models of the system are deployed in the simulation platform ‘MuJoCo-2.0’ [28] as a surrogate to their analytical models. The models are imported from the OpenAI gym [29] and Deepmind’s control suite [30].

Figure 1 (left) shows the open-loop training plots of various systems under consideration. The data is averaged over 5 experimental runs to reflect the training variance and the repeatability of D2C-2.0 over multiple training sessions. One might notice that the episodic cost is not monotonically decreasing during training. This is because of the flexibility provided in training for the episodic cost to stay within a ‘band’ around the cost at the previous iteration. Though this might momentarily increase the cost in subsequent training iterations, it is observed to have improved the overall speed of convergence.

In this section, we compare D2C-2.0 with D2C and DDPG, and discuss the merits and demerits of each of the methods.

The manuscript on D2C provided an in-depth analysis on the comparison of D2C with DDPG. We have shown that training in D2C is much more sample-efficient and fares well w.r.t. the training time and the reproducibility of the resulting policy.

Table 1 shows the training times of D2C, DDPG and D2C-2.0. It is evident from the table that for simple and lower dimensional systems such as pendulum and cart pole, the policy can be almost calculated online (parallelization could make it much faster). It is due to the tendency of ILQR/DDP to quickly converge if the linearization of system models has bigger basin of attraction (such as pendulum, which is, in fact, often approximated to be linear about its resting state). Moreover, for a higher dimensional system such as in the case of a robotic fish, D2C-2.0 took 54.8 seconds while the original D2C took around 6096 seconds, which is approximately 90 times faster. The table makes it evident that D2C-2.0 clearly exhibits much faster convergence when compared to both D2C and DDPG, while D2C is much faster than DDPG. The primary reason for the much larger training times in DDPG is due to the complex parametrization of the feedback policy in terms of deep neural networks which runs into hundreds of thousands of parameters while D2C has a parametrization, primarily the open loop control sequence which is of dimension $pT$, where $p$ is the number of inputs and $T$ is the time horizon of the problem, that runs from tens to a few thousands (see Table 2).

A major reason behind the differences in training efficiencies stems from their formulation. D2C is based on gradient descent, which is linearly convergent. It is a direct and generic method of solving an optimization problem. On the other hand, the equations involved in D2C-2.0 arise from the optimality conditions that inherently exploit the recursive optimal structure and the causality in the problem. Moreover, the equations are explicit in the Jacobians of the system model. As a result, it is possible to solve for it independently and in a sample efficient manner, at every time-step of each iteration. However, one should expect that as we consider very high dimensional systems, the algorithmic complexity of D2C-2.0 could outgrow that of D2C. This is because it is $O(n_{x}^{2})$ in terms of its computational complexity. Also, D2C, by the nature of its formulation does not need the state of the system, and hence, is generalizable to problems with partial state observability.

The closed-loop performance of D2C-2.0 w.r.t. that of DDPG is shown in the figure 1 (right). Noise is injected into the system via the control channel by varying the noise scaling factor-$\epsilon$ and then, measuring the performance. We note here that DDPG fails to converge for the robotic fish problem, and thus, the nominal performance ($\epsilon=0$) is poor. We see that the performance of D2C is better than that of DDPG when the noise levels are lower, however, the performance of DDPG remains relatively similar for higher levels of noise while that of D2C degrades. In our opinion, this closed loop performance aspect of the methods needs to be studied in more detail in order to make a definitive judgement, and is left for future work. For instance, by the time DDPG overtakes D2C in terms of the average performance, the actual performance of both systems is quite poor, and thus, it is not clear if this can be construed as an advantage.

Now, as mentioned in section 2, earlier works have attempted to tackle the idea of model-free ILQR and have mostly confined to finite-differences for Jacobian computation. Table 3 shows the comparison of per-iteration computational times between ‘finite-differences’ and ‘linear least squares with central difference formulation’. It is clearly evident that, as the dimension of the state space increases, the method of finite-differences requires many more function evaluations, and hence, our LLS-CD method is much more efficient.

In order to conclude the discussion that started the global vs local policy approximation dilemma in section 1, how well deep RL methods extrapolate their policy to the entire continuous state and action spaces remains an open question. Now, consider that the results shown in the Table 1 are based on serial implementations on an off-the-shelf computer and D2C-2.0 can be highly parallelizable. By looking at the order of magnitude of numbers in the same table, we expect with reasonable augmentation of computational power by parallelization, D2C-2.0 could offer a near-real time solution in high dimensional problems. In such cases, one could rely on the near-optimal policy described in this paper for low noises, and re-solve for the open-loop trajectory online by the approach presented in this paper along with the attendant feedback, whenever the noise takes the method out of the region of attraction of the linear feedback policy. This will be our motivation for the future work along with addressing the following limitation.
Limitations: D2C-2.0, in its current form, is not reliable for hybrid models such as legged robots that involve impacts. Primarily, the issue is in the feedback design which requires smoothness of the underlying dynamics while legged systems are only piece-wise smooth. Future work shall model impacts to tractably incorporate them in order to estimate the required gradients. This seems plausible as iLQR has been used to design controllers for such systems in the literature [25].

Proof of Lemma 1:

We proceed by induction. The first general instance of the recursion occurs at $t=3$. It can be shown that:

Noting that $\bar{S}_{1}(.)$ and $\bar{S}_{2}(.)$ are second and higher order terms, it follows that $\bar{\bar{S}}_{3}$ is $O(\epsilon^{2})$.w
Suppose now that $\delta{\bf x_{t}}=\delta{\bf x_{t}^{l}}+\bar{\bar{S}}_{t}$ where $\bar{\bar{S}}_{t}$ is $O(\epsilon^{2})$. Then:

Noting that $\bar{S}_{t+1}$ is $O(\epsilon^{2})$ and that $\bar{\bar{S}}_{t+1}$ is $O(\epsilon^{2})$ by assumption, the result follows.

Proof of Proposition 1:

From (6), we get,

The first equality in the last line of the equations before follows from the fact that $\mathbb{E}[\delta{\bf x_{t}^{l}}]=0$, since its the linear part of the state perturbation driven by white noise and by definition $\delta{\bf x_{1}^{l}}=0$. The second equality follows from the fact that $\delta J_{2}^{\pi}$ is an $O(\epsilon^{2})$ function. Now,