Sufficiently Accurate Model Learning for Planning

Dynamics models play an essential role in many modern controllers and planners. This is for instance the case in Model Predictive Control (MPC) [1]. They can also be used to compute feedforward terms for feedback controllers [2], or be used with a planner to find a trajectory from which stabilizing controllers can be generated [3, 4]. While model-free Reinforcement Learning techniques can solve similar problems without the use of a dynamics model [5, 6, 7], they generally do not scale to new problems or changes in problem parameters. Model learning methods perform admirably when the models can approximate the dynamics of the system accurately. However, the performance of these controllers can be degraded with uncertainty in the model. While robust controllers can be designed to attempt to alleviate these issues [8], these methods typically perform conservatively due to having to cater to worst case approximations of the model. In addition, there may be effects that a robust controller designer may not be aware. For example, consider the problem of landing a quadrotor precisely at a target as shown in Figure 9. There are complex aerodynamic effects associated when nearby surfaces cause disturbances to the airflow. This may result in large torques when the quadrotor is hovering close to the ground and hamper precise landings. This is known as the “ground effect.” These aerodynamic effects can be hard to model from just prior knowledge and may show up as a highly correlated, state dependent, non-zero mean noise. A common method that has been suggested to model similar effects is to learn or adjust a dynamics model with real data taken from running the system.

For linear systems, this method of System Identification (SI) or Model Learning has many results [9, 10, 11, 12, 13, 14] —such as recoverability of linear system dynamics when the data contains sufficient excitation. System identification methods have been proposed for non-linear systems [15, 16, 17]. However, they suffer from issues such as the need for a large amount of data or the requirement of a special problem structure such as block-based systems.

While generic methods for non-linear systems do not provide the same theoretical guarantees, there has been some success in practice. In robotics, Gaussian Processes, Gaussian Mixture Models, or Neural Networks have been used to learn models of dynamics [18, 19, 20, 21]. A typical process for learning these models involves selecting a parameterized model, such as a neural network with a fixed number of layers and neurons, and choosing a loss function that penalizes the output of the model for not matching the data gathered from running the real system. Then, one optimizes the parameters by minimizing the empirical risk using, for instance, stochastic gradient descent like algorithms. This formulation assumes that all transitions are equally important, since it penalizes the mismatch between model and data uniformly on all portions of the state-action space. While this formulation has shown success in some applications, prior knowledge about the task and system can inform better learning objectives. A control designer may know that a certain part of the state space requires a certain accuracy for a robust controller to work well, or that some part of the state space is more important and should have hard constraints on the model accuracy. For example, to precisely land a quadrotor, a designer may note that the accuracy of modeling the complex ground effect forces is more important near the landing site. Incorporating this prior knowledge can lead to better performing models.

To address the problem of incorporating prior knowledge into model learning, we introduced the idea of Sufficiently Accurate Model Learning in [22]. This formulation is based on the inclusion of constraints in the optimization problem whose role is to introduce prior-knowledge about the importance of different state-control subsets. In the example of the quadrotor, notice that when the quadrotor is away from the surfaces, the ground effect is minor and thus, it is important to focus the learned model’s expressiveness in the region of the state-space that is most heavily affected. This can be easily captured by a constraint that the average error in the important state-space regions is smaller than a desired value. These constraints will allow models with finite expressiveness concentrate on modeling important aspects of a system. One point to note is that this constrained objective can be used orthogonally to many existing methods. For example, the constrained objective can replace the unconstrained objective in [20, 21], and all other aspects of the methods can remain the same. The data can be collected the same way. The idea of using extra sensors during training for [20] need not change. While not trivial, the idea of constraints can also be applied to Gaussian process models such as [18, 19].

In its most generic form, the problem of model learning is an infinite dimensional non-convex optimization problem that involves the computation of expectations with respect to an unknown distribution over the state-action space. In addition, the formulation proposed here introduces constraints which seems to make the learning process even more challenging. However, in this work we show that solving this problem accurately is not more challenging than solving an unconstrained parametric learning problem. To reach this conclusion we solve a relaxation of the problem of interest with three common modifications: (i) function parameterization, (ii) empirical approximation, and (iii) dual problem solving. Function parameterization turns the infinite dimensional problem into one over finite function parameters. Empirical approximation allows for efficient computation of approximate expectations, and solving the dual problem leads to a convex unconstrained optimization problem. The three approximations introduced however may not yield solutions that are good approximations of the original problem. To that end, we establish a bound on the difference of the value of these solutions. This gap between the original and approximate problem depends on the number of samples of data as well as the expressiveness of the function approximation (Theorem 1). In particular, the bound can be made arbitrarily small with sufficient number of samples and with the selection of a rich enough function approximator. This implies that solving the functional constrained problem is nearly equivalent to solving a sequence of unconstrained approximate problems using primal-dual methods.

This paper extends [22] to the case of empirical samples and presents Theorem 1 that relates number of samples, function approximation expressiveness, and loss function smoothness to approximation error. In addition, there is experimental validation of Theorem 1 as well as different examples to showcase the framework of Sufficiently Accurate model learning. This paper is organized as follows: Section II introduces the Sufficiently Accurate model learning framework, Section III presents in detail the three approximations introduced to solve the problem as well as a result that bounds the error on the approximate problem (Theorem 1). Section IV provides the proof for the main theorem. Section V presents a simple primal-dual method to solve the constrained problem, while experimental results are presented in Section VI on a double integrator system with friction, a ball bouncing task, and a quadrotor landing with ground effect in simulation. Theoretical results are experimentally tested for the double integrator system. Section VII presents the paper’s conclusions.

In this paper we consider the problem of learning a discrete time dynamical system. Let us denote $t\in\mathbb{Z}$ as the time index and let $x_{t}\in\mathbb{R}^{n}$, $u_{t}\in\mathbb{R}^{p}$ be the state and input of the system at time $t$. The dynamical system of interest is represented by a function $f\mathrel{\mathop{\ordinarycolon}}\mathbb{R}^{n}\times\mathbb{R}^{p}\to\mathbb{R}^{n}$ that relates the state and input at time $t$ to the state at time $t+1$

One approach to System Identification or Model Learning consists of fitting an estimated dynamical model, $\phi\mathrel{\mathop{\ordinarycolon}}\mathbb{R}^{n}\times\mathbb{R}^{p}\rightarrow\mathbb{R}^{n}$, to state transition data [18]. This state transition data consists of tuples of $s=(x_{t},u_{t},x_{t+1})$ drawn from a distribution $\mathcal{S}_{\mathcal{D}}$ with the sample space $\mathcal{S}\subseteq\mathbb{R}^{n}\times\mathbb{R}^{p}\times\mathbb{R}^{n}$. The estimated model $\phi$ belongs to a class of functions, $\Phi$, of which $f$ is an element. $\Phi$ could be, for example, the space of all continuous functions. Then, the problem of model learning reduces to finding the function $\phi$ in the class $\Phi$ that best fits the transition data. The figure of merit is a loss function $\ell\mathrel{\mathop{\ordinarycolon}}\mathcal{S}\times\Phi\to\mathbb{R}$. With these definitions, the problem of interest can be written as the following stochastic optimization problem

The loss function needs to be selected to encourage the model $\phi$ to match its output with the state transition data. An example of a loss function is the p-norm, $l(s,\phi)=\mathinner{\!\left\lVert\phi(x_{t},u_{t})-x_{t+1}\right\rVert}_{p}$. For $p=1$, this reduces to a sum of absolute differences between each output dimension of $\phi$ and the true next state, $x_{t+1}$. When $p=2$, this is simply a euclidean loss. Other common losses can include the Huber loss and, in discrete state settings, a 0-1 loss.

Often times, one can derive, from first principles, a model $\hat{f}$ of the dynamical system of interest. Depending on the complexity of the system of, these models may be inaccurate since they may ignore hard to model dynamics or higher order effects. For instance, one can derive a model $\hat{f}$ for a quadrotor from rigid body dynamics where forces and torques are functions of each motor’s propeller speed. However, the accuracy of this model will depend on other effects that are harder to model such as aerodynamic forces near the ground. If these aerodynamic effects are ignored, it can result in a failure to control the system or in poor performance [23]. In these cases, the target model, denoted by $\tilde{\phi}$, can decomposed as the sum of an analytic model and an error

The learning of the error term—or residual model—fits the framework described in (2). For instance, for the $p$ norm loss can be modified to take the form $l(s,\phi)=\mathinner{\!\left\lVert\tilde{\phi}(x_{t},u_{t})-x_{t+1}\right\rVert}_{p}=\mathinner{\!\left\lVert\phi(x_{t},u_{t})-(x_{t+1}-\hat{f}(x_{t},u_{t}))\right\rVert}_{p}$.

A characteristic of the classic model learning problem defined in (2) is that errors are uniformly weighted across the state-input space. In principle, one can craft the loss so to represent different properties on different subsets of the state-input space. However, this design is challenging, system dependent, and dependant on the transition data available. In contrast, our approach aims to exploit prior knowledge about the errors and how they impact the control performance. For instance, based on the analysis of robust controllers one can have bounds on the error required for successful control. This information can be used to formulate the Sufficiently Accurate model learning problem, where we introduce the prior information in the form of constraints. Formally, we encode the prior information by introducing $K\in\mathbb{N}$ functions $g_{k}\mathrel{\mathop{\ordinarycolon}}\mathcal{S}\rightarrow\mathbb{R}$. Define in addition, a collection of subsets of transition tuples where this prior information is relevant, $\mathcal{S}_{k}\subset\mathcal{S}$ and corresponding indicator functions $\mathbb{I}_{k}(s)\mathrel{\mathop{\ordinarycolon}}\mathcal{S}\rightarrow\{0,1\}$ taking the value one if $s\in\mathcal{S}_{k}$ and zero otherwise. With these definitions the sufficiently accurate model learning problem is defined as

Note that the sets $\mathcal{S}_{k}$ that define the indicator variables $\mathbb{I}_{k}(s)$ are not necessarily disjoint. In fact, in practice, they are often not. The sets can be arbitrary and have no relation to each other. Examples of how these sets are used are given in examples at the end of this section. Notice that the (4) is an infinite dimensional problem since the optimization variable is a function and it involves the computation of expectations with respect to a possibly unknown distribution. An approximation to this problem is presented in Section III. For technical reasons, the functions $g_{k}$ and $l$ should be expectation-wise Lipschitz continuous.

The expectation-wise Lipschitz assumption is a weaker assumption than Lipschitz-continuity, as any Lipschitz-continuous function with a Lipschitz constant $L$ is also expectation-wise Lipschitz-continuous with a constant of $L$. In particular, the loss functions in Example 1 and 2 are expectation-wise Lipschitz-continuous with some constant (cf., Appendix A). There is no assumption that the functions should be convex or continuous.

Before we proceed, we present two examples of Sufficiently Accurate model learning. For notational brevity, when an expectation does not have a subscript, it is always taken over $s\sim\mathcal{S}_{\mathcal{D}}$.

The unconstrained problem (2) and the constrained problem (4) are functional optimization problems. In general, these are infinite dimensional and usually intractable. Instead of optimizing over the entire function space $\Phi$, one may look at function spaces, $\Phi_{\theta}\subset\Phi$, parameterized by a $d$-dimensional vector $\theta\in\Theta=\mathbb{R}^{d}$. Examples of these classes of functions are linear functions of the form $\phi_{\theta}(x,u)=\theta_{x}^{\top}x+\theta_{u}^{\top}u$ where $\theta=[\theta_{x},\theta_{u}]$ is a vector of weights for the state and control input. More complex function approximators, such as neural networks, may be used to express a richer class of functions [21, 24]. Restricting the function space poses a problem in that the optimal solution to (4) may no longer exist in the set $\Phi_{\theta}$. The goal under these circumstances should be to find the closest solution in $\Phi_{\theta}$ to the true optimal solution $\phi^{\star}$. Additionally, the expectations of the loss and constraint functions are in general intractable. The distributions can be unknown or hard to compute in closed form. In practice, the expectation is approximated with a finite number of data samples $s_{i}\sim\mathcal{S}_{\mathcal{D}}$ with $i=1,\ldots,N$. This yields the following empirical parameterized risk minimization problem

While both function and empirical approximations are common ways to simplify the problem, the approximate problem introduced in (8) is still a constrained optimization problem and can be difficult to solve in general as it can be nonconvex in the parameters $\theta$. This is the case for instance when the function approximator is a neural network. One approach to solve this problem is to solve the dual problem associated with (8). To aid in the definition of the dual problem, we first define the dual variables (also known as Lagrange multipliers), $\lambda\in\mathbb{R}_{+}^{K}$, along with the Lagrangian associated with (8)

Here, the symbol, $\ell_{0}(s_{i},\phi_{\theta})$, is defined as $l(s_{i},\phi_{\theta})\mathbb{I}_{0}(s)$ to condense the notation. Similarly, the bolded vector, $\bm{g}(s_{i},\phi_{\theta})$ is a vector where the $k^{\text{th}}$ entry is defined as $g_{k}(s_{i},\phi_{\theta})\mathbb{I}_{k}(s)$. The dual problem is now defined as

Notice that (10) is similar to a regularized form of (8) where each constraint is penalized by a coefficient $\omega_{k}$

Adding this type of regularization can weight certain state-action spaces more. In fact, if $\omega_{k}$ is chosen to be $\lambda_{N}^{\star}$, solving (11) would be equivalent to solving (10). However, arbitrary choices of $\omega_{k}$ provide no guarantees on the constraint function values. By defining the constraint functions directly, constraint function values are determined independent of any tuning factor. For problems where strong guarantees are required or easier to define, the Sufficiently Accurate framework will satisfy them by design. An alternative interpretation is that (10) provides a principled way of selecting the regularization coefficients. In Section V, we discuss an implementation of a primal dual algorithm to do so.

The dual problem has two important properties that hold regardless of the structure of the optimization problem (8). For any $\theta$ that minimizes the Lagrangian $\mathcal{L}_{N}$, the resulting function—termed the dual function—is concave on the multipliers, since it is the point-wise minimum of linear functions (see e.g. [25]). Therefore, its maximization is tractable and it can be done efficiently for instance using stochastic gradient descent. In addition, the dual function is always a lower bound on the value $P_{N}^{\star}$ and in that sense solving the dual problem (10) provides the tightest lower bound. In the case of convex problems (that fulfill Slater’s Condition), it is well known that the problems have zero duality gap, and therefore $P_{N}^{\star}=D_{N}^{\star}$ [25, Section 5.2.3]. However, the problem (8) is non-convex and a priori we do not have guarantees on how far the values of the primal and the dual are. Moreover, recall that the primal problem in (8) is an empirical approximation of the problem that we are actually interested in solving (4).

The previous discussion leads to the question about the quality of the solution (10) as an approximation to (4). The duality gap is defined as the difference between the primal and dual solutions of the same problem. Here, the gap is the difference between the primal and the dual of different but closely related problems. Hence, the quantity we are interested in bounding is the surrogate duality gap defined as

To provide specific bounds for the difference in the previous expression we consider the family of function classes $\Phi_{\theta}$ termed $\epsilon$-universal function approximators. We define this notion next.

To provide some intuition on the definition consider the case where $\Phi$ is the space of all continuous function, the above property is satisfied by some neural network architecture. That is, for any $\epsilon$, there exists a class of neural network architectures, $\Phi_{\theta}$ such that $\Phi_{\theta}$ is an $\epsilon$-universal approximator for the set of continuous functions [24, Corollary 2.1]. Thus, for any dynamical system with continuous dynamics, this assumption is mild. Other parameterizations, such as Reproducing Kernel Hilbert Spaces, are $\epsilon$-universal as well [26]. Notice that the previous definition is an approximation on the total norm variation and hence it is a milder assumption than the universal approximation property that fully connected neural networks exhibit [24].

Next, we define an intermediate problem on which the surrogate duality gap depends: a perturbed version of problem (4) where the constraints are relaxed by $L\epsilon\geq 0$ where $L$ is the constant defined in Assumption 1 and $\epsilon$ the universal approximation constant in Definition 1

$\bm{1}$ is a vector of ones. The perturbation results in a problem whose constraints are tighter as compared to (4). The set of feasible solutions for the perturbed problem (13) is a subset of the feasible solutions for the unperturbed problem (4). The perturbed problem accounts for the approximation introduced by the parameterization. In the worst case scenario, if the problem (13) is infeasible, the parameterized approximation of (8) may turn infeasible as the number of samples increases.

Let $\lambda_{L\epsilon}^{\star}$ be the solution to the dual of (13)

With these definitions, we can present the main theorem that bounds the surrogate duality gap.

The intuition behind the theorem is that given some acceptable surrogate duality gap, $\delta$, there exists a neural network architecture, $\Phi_{\theta}$, and a number of samples, $N$ such that the probability that the solution to (10) is within $\delta$ to the solution to (4) is very high. The choice of neural network will influence the value of $\epsilon$ and $\lambda_{L\epsilon}^{\star}$. These in turn will decide the duality gap, $\delta$, as the quantity $\delta-\epsilon L(\mathinner{\!\left\lVert\lambda_{L\epsilon}^{\star}\right\rVert}_{1}+1)$ must be positive. A larger neural network will correspond to a smaller $\epsilon$ which will also has an impact on the perturbed problem (13). A smoother function and smaller $\epsilon$ will lead to smaller perturbations. Smaller perturbations can lead to a smaller dual variable, $\lambda_{L\epsilon}^{\star}$. Thus, larger neural networks and smoother dynamic systems will have smaller duality gaps. If $L\epsilon$ is large, then the perturbed problem may be infeasible. In theses cases, $\lambda_{L\epsilon}^{\star}$ will be infinite. This corresponds to problems where the function approximation simply can not satisfy the original constraint functions. For example, using constant functions to approximate a complicated system may violate the constraint functions $g_{k}$ for all possible constant functions. Thus, no $\delta$ exists to bound the solution as the parameterized empirical problem (8) has no feasible solution. This theorem suggests that with a good enough function approximation and large enough $N$, solving (10) is a good approximation to solving (4) with large probability.

There are some details to point out in Theorem 1. First, the function $H^{\Phi_{\theta}}$ a complicated function that will usually scale with the size of the neural network. A larger neural network will lead to a smaller $\epsilon$, but may require a larger number of samples $N$ to adequately converge to an acceptable solution. The assumption on the limiting behavior of $H^{\Phi_{\theta}}$ is fufilled by some neural network architectures [29], but the general behavior of this function for all neural network architectures is still a topic of research. Additionally, we assume the space of Lagrange multipliers is a compact set. This will imply, along with compact state-action space, that $\mathcal{L}$ is bounded. A finite Lagrange multiplier is a reasonable assumption as the problem (4) is required to be feasible [27].

The bound established in Theorem 1 depends on quantities that are in general difficult to estimate, These include $H^{\Phi_{\theta}}$, $\mathinner{\!\left\lVert\lambda_{L\epsilon}^{\star}\right\rVert}_{1}$, $L$, $\epsilon$. Thus, while this theorem provides some insights on how these quantities influence the gap between solutions, it is mainly a statement of the existence of such values that can provide a desired result. In practice, this result can be achieved by choosing increasing the sizes of neural networks as well as data samples until the desired performance is reached. Note that the theorem follows our intuition that larger neural networks and more data will give us more accurate result. However, this theorem formalizes not only that it is more accurate, but that the error will tend to 0 as number of samples and number of parameters increase.

To begin, we define an intermediate problem

Note that this is the unperturbed version of (13). As a reminder, this problem uses a class of parameterized functions, but does not use data samples to approximate the expectation. Thus, it can be seen as a step in between (4) and (8). As with the dual problem to (8), we can define the Lagrangian associated with (17)

and the dual problem

Using this intermediate problem, we can break the bound $\mathinner{\!\left\lvert P^{\star}-D_{N}^{\star}\right\rvert}$ into two components.

As a reminder, $P^{\star}$ is the solution to the problem we want to solve in (4), $D_{N}^{\star}$ is the solution to the problem (10) we can feasibly solve, and $D_{\theta}^{\star}$ is the solution to an intermediate problem (19). The first half of this bound, $\mathinner{\!\left\lvert P^{\star}-D_{\theta}^{\star}\right\rvert}$, is the error that arises from using a parameterized function and dual approximation. The second half of this bound, $D_{\theta}^{\star}-D_{N}^{\star}$, is the error that arises from using empirical data samples. It can be seen as a kind of generalization error. The proof will now be split into two parts that will find a bound for each of these errors.

Section III has shown that problem (10) can approximate (4) given a large enough neural network and enough samples. This section will discuss how to compute a solution (10). There are many primal-dual methods [32, 33, 34] in the literature to solve this exact problem, and Algorithm 1 is an example of a simple primal-dual algorithm. One way to approach this problem is to consider the optimal dual variable, $\lambda_{N}^{\star}$. Given knowledge of $\lambda_{N}^{\star}$, the problem reduces to the following unconstrained minimization

A possible solution method is to start with an estimate of $\lambda_{N}^{\star}$, and solve the minimization problem. Then holding the primal variables fixed, update the dual variables by solving the outer maximization. This method can be seen as solving a sequence of unconstrained minimization problems. This method can be further approximated; instead of fully minimizing with respect to the primal variables, a gradient descent step can be taken. And instead of fully maximizing with respect to the dual variables, a gradient ascent step can be taken. This leads to Algorithm 1 where we iterate between the inner minimization step and the outer maximization step. At each iteration, dual variables are projected onto the positive orthant of $\mathbb{R}^{K}$, denoted by the projection operator, $[\lambda]_{+}$. This is to ensure non-negativity of the dual variables.

In many cases, the full gradient of $\nabla_{\theta}\mathcal{L}_{N}(\theta,\lambda)$ and $\nabla_{\lambda}\mathcal{L}_{N}(\theta,\lambda)$ can be too expensive to compute. This is due to the possibly large number of samples $N$. An alternative is to take stochastic gradient descent and ascent steps. The gradients can be approximated by taking $M$ random samples of the whole dataset $S={s_{1},\ldots,s_{N}}$. The samples will be denoted as $\hat{S}={s_{i_{1}},\ldots,s_{i_{M}}}$ where $i_{1}$ is an integer index into whole dataset $S$. Using $\hat{S}$, we obtain

The gradients $\nabla_{\theta}\ell_{0}(s_{i_{j}},\phi_{\theta})$ and $\nabla_{\theta}\bm{g}(s_{i_{j}},\phi_{\theta})$ can be computed easily using backpropogation. Similarly, for $\nabla_{\lambda}\mathcal{L}_{N}(\theta,\lambda)$,

The dual gradient can be estimated as simply the average of the constraint functions over the sampled dataset.

In the simplest form of the primal-dual algorithm, the variables are updated with simple gradient ascent/descent steps. These updates can be replaced with more complicated update schemes, such as using momentum [35] or adaptive learning rates [36]. Higher order optimization methods such as Newton’s method can be used to replace the gradient ascent and descent steps. For large neural networks, this can be unfeasible as it requires the computation of Hessians with respect to neural network weights. The memory complexity for the Hessian is quadratic with the number of neural network weights.

The primal-dual algorithm presented here is not guaranteed to converge to the global optimum. With proper choice of learning rate, it can converge to a local optimum or saddle point. This issue is present in unconstrained formulations like (2) as well. An example of the evolution of the loss and constraint functions is shown in Figure 1.

This section shows examples of the Sufficiently Accurate model learning problem. First, experiments are performed using a simple double integrator experiencing unknown dynamic friction. The simplicity of this problem along with the small state space allows us to explore and visualize some of the properties of the approximated solution. Next, two more interesting examples are shown. One example learns how a ball bounces on a paddle with unknown paddle orientation and coefficient of restitution. The other example mitigates ground effects which can disturb the landing sequence of a quadrotor. The experiments will compare the sufficiently accurate problem (4) with the unconstrained problem (2) which will be denoted as the Uniformly Accurate problem. Each experimental subsection will be broken down into three parts, 1) System and Task introduction, 2) Experimental details, and 3) Results.