Reach-SDP: Reachability Analysis of Closed-Loop Systems with Neural Network Controllers via Semidefinite Programming

Haimin Hu, Mahyar Fazlyab, Manfred Morari, George J. Pappas

Introduction

Deep neural networks (DNN) have seen renewed interest in recent years due to the proliferation of data and access to more computational power. In autonomous systems, DNNs are either used as feedback controllers , motion planners , perception modules, or end-to-end controllers . Despite their high performance, DNN-driven autonomous systems lack formal safety and stability guarantees. Indeed, recent studies show that DNNs can be vulnerable to small perturbations or adversarial attacks . This issue is more pronounced in closed-loop systems, as a small perturbation in the loop can dramatically change the behavior of the closed-loop system over time. Therefore, it is of utmost importance to develop tools for verification of DNN-driven control systems. The goal of this paper is to develop a methodology, based on semidefinite programming, for safety verification and reachability analysis of linear dynamical systems in feedback interconnection with DNNs.

Safety verification or reachability analysis aims to show that starting from a set of initial conditions, a dynamical system cannot evolve to an unsafe region in the state space. Methods for reachability analysis can be categorized into exact (complete) or approximate (incomplete), which compute the reachable sets exactly and approximately, respectively. Verification of dynamical systems has been extensively studied in the past . More recently, the problem of output range analysis of neural networks has been addressed in , mainly motivated by robustness analysis of DNNs against adversarial attacks . Compared to these bodies of work, verification of closed-loop systems with neural network controllers has been less explored . In , a method for verification of sigmoid-based neural networks in feedback with a hybrid system is proposed, in which the neural network is transformed into a hybrid system and then a standard verification tool for hybrid systems is invoked. In , a new reachability analysis approach based on Bernstein polynomials is proposed that can verify DNN-controlled systems with Lipschitz continuous activation functions. Dutta et al. use a flow pipe construction scheme to over approximate the reachable sets. A piecewise polynomial model is used to provide an approximation of the input-output mapping of the controller and an error bound on the approximation. This approach, however, is only applicable to Rectified Linear Unit (ReLU) activation functions.

Contributions. In this paper, we propose a semidefinite program (SDP) for reachability analysis of linear time-varying dynamical systems in feedback interconnection with neural network controllers equipped with a projection operator, which projects the output of the neural network (the control action) onto a specified set of control inputs. Our technical approach relies on abstracting the nonlinear activation functions as well as the projection operator by quadratic constraints , which leads to an outer-approximation of forward reachable sets of the closed-loop system. We show that we can compute these approximate reachable sets using semidefinite programming. Our approach can be used to analyze control policies learned by neural networks in, for example, model predictive control (MPC) and constrained reinforcement learning . To the best of our knowledge, our result is the first convex-optimization-based method for reachability analysis of closed-loop systems with neural networks in the loop. We illustrate the utility of our approach in two numerical examples, in which we certify finite-time reachability and constraint satisfaction of a double integrator and a quadrotor.

which can be equivalently expressed as the quadratic inequality

Problem Formulation

We consider a discrete-time linear time-varying system

We denote the closed-loop system with dynamics (5) and the projected neural network control policy (9) by

which is a non-smooth nonlinear system because of the presence of the nonlinear activation functions in the neural network and the projection operator.

In this paper, we only consider box constraints for the input and leave more sophisticated constraints such as polytopes to future work. These constraints are useful in, for example, MPC design problems. In a robust MPC controller is approximated by a neural network equipped with a projection operator to ensure satisfaction of polytopic input constraints.

2 Finite-Time Reach-Avoid Verification Problem

holds true for (10). There exist efficient methods and software implementations for testing set inclusion (12a) and set intersection (12b). However, computing exact reachable sets for the nonlinear closed-loop system (10) is, in general, computationally intractable. Therefore, we resort to finding outer-approximations of the closed-loop reachable sets, Rˉt(X0)⊇Rt(X0)\bar{\mathcal{R}}_{t}(\mathcal{X}_{0})\supseteq\mathcal{R}_{t}(\mathcal{X}_{0}), and use them to test if,

is true. Note that (13) are sufficient conditions for (12). To obtain meaningful certificates, we want the approximations Rˉt(X0)\bar{\mathcal{R}}_{t}(\mathcal{X}_{0}) to be as tight as possible. Thus our goal is to compute the tightest outer-approximations of the tt-step reachable sets. This can be addressed, for example, by solving the following optimization problem,

The solution to the above problem is the minimum-volume outer-approximation of the tt-step reachable set of the closed-loop system (10). In the following sections, we will derive a convex relaxation to the optimization problem (14).

Problem Abstraction via Quadratic Constraints

We begin with a formal definition of QCs for sets .

holds for all x∈Xx\in\mathcal{X}. Then we say X\mathcal{X} satisfies the QC defined by Q\mathcal{Q}. The vector [x⊤ 1]⊤\begin{bmatrix}x^{\top}\ 1\end{bmatrix}^{\top} is called the basis of this QC.

Here, for each fixed QQ, the set of xx’s satisfying (15) is a superset of X\mathcal{X}. Indeed, we have that

In this paper, we mainly use polytopes and ellipsoids as the initial set X0\mathcal{X}_{0}. Nonetheless, as addressed by , other types of sets such as hyper-rectangles and zonotopes are also applicable in this setting.

where mm is the dimension of bb. The basis is [x⊤ 1]⊤\begin{bmatrix}x^{\top}\ 1\end{bmatrix}^{\top}.

with basis [x⊤ 1]⊤\begin{bmatrix}x^{\top}\ 1\end{bmatrix}^{\top}.

To this end, we have over-approximated the initial set X0\mathcal{X}_{0} and represented it as QCs. We will see in Section 4 that the matrix P∈PP\in\mathcal{P} appears as a decision variable in the SDP. It provides an extra degree of freedom towards a less conservative estimation of the reachable set Rˉ(X0)\bar{\mathcal{R}}(\mathcal{X}_{0}).

2 Reachable Set

In order to facilitate the relaxation of the volume-minimization problem in (14), we assume the candidate set Rˉ(X0)\bar{\mathcal{R}}(\mathcal{X}_{0}) that over-approximates R(X0)\mathcal{R}(\mathcal{X}_{0}) is represented by the intersection of finitely many quadratic inequalities:

2.2 Ellipsoidal reachable set

For polytopic reachable sets, if the facets aia_{i} are properly chosen, the resulting outer-approximations can be very tight, as we will see in Section 5.1. However, finding facets of a higher dimensional polytope can be prohibitively challenging. Ellipsoidal reachable sets scale better and therefore are more suitable in this case.

3 ReLU Neural Networks with Projection

In this subsection, we show how to abstract the nonlinear activation functions of the neural network by QCs. We will focus on the ReLU function throughout the paper. Other types of activation functions such as sigmoid and tanh can also be represented by QCs. See for details. Recall that the input constraint sets are Ut={ut∣u‾t≤ut≤uˉt}\mathcal{U}_{t}=\{u_{t}\mid\underline{u}_{t}\leq u_{t}\leq\bar{u}_{t}\}. Consider the following recursion,

for all i=1,⋯ ,di=1,\cdots,d. In addition, the ReLU ϕ(x)\phi(x) is also slope-restricted on $withrepeatednonlinearity.Thisallowsustowrite(4)withwith repeated nonlinearity. This allows us to write (4) with\alpha=0andand\beta=1$,

By taking a weighted sum of the constraints (23) and (24), we obtain the single quadratic constraint,

As we will see in Section 4, the QQ matrix in (26) will appear as a decision variable in the SDP problem. Note that there are many ways to refine (26) to yield a tighter relaxation for a specific region in the state-space, such as using interval arithmetic or linear programming .

Reach-SDP: Computing Forward Reachable Sets via Semidefinite Programming

In this section, we propose Reach-SDP, an optimization-based approach that uses the QC abstraction developed in the previous section to estimate the reachable set of the closed-loop system (10). Specifically, the NN-step reachable set is estimated using the following recursive computations,

for t=0,⋯ ,N−1t=0,\cdots,N-1. In the sequel, we discuss how to implement Reach_SDP⁡\operatorname{Reach\_{SDP}} in detail.

In the previous section, we have abstracted the initial set, the reachable set and the projected neural network with QCs or quadratic inequalities, each with a different basis vector. For Reach-SDP, we unify those quadratic terms with the same basis vector:

Substitute Eξb=xbE\xi_{b}=x_{b} into xb⊤Qxb≥0x_{b}^{\top}Qx_{b}\geq 0 and we have xb⊤E⊤QExb=ξb⊤Mξb≥0x_{b}^{\top}E^{\top}QEx_{b}=\xi_{b}^{\top}M\xi_{b}\geq 0. ∎

The quadratic inequalities (19) defined by matrices SiS_{i} in (20) or (21), describing the candidate set Rˉ(X0)\bar{\mathcal{R}}(\mathcal{X}_{0}), are equivalent to the quadratic inequalities defined by matrices Mout(Si)=Eout⊤SiEoutM_{\rm{out}}(S_{i})=E_{\rm{out}}^{\top}S_{i}E_{\rm{out}} with basis [x⊤ 1]⊤[\mathbf{x}^{\top}\ 1]^{\top}. The change-of-basis matrix is

2 Over-Approximating the One-Step Reachable Set

In the next theorem, we state our main result for over-approximating the one-step reachable set R(X0)\mathcal{R}(\mathcal{X}_{0}) for the closed-loop system (10).

Since the initial set X0\mathcal{X}_{0} satisfies the QC defined by Min\mathcal{M}_{\rm{in}}, we have,

for all x0∈X0x_{0}\in\mathcal{X}_{0}. Similarly, the projected neural network controller Proj⁡Ut(π(⋅))\operatorname{Proj}_{\mathcal{U}_{t}}(\pi(\cdot)) satisfying the QC defined by Mmid\mathcal{M}_{\rm{mid}} implies that,

By left- and right-multiplying both sides of (34) by [x⊤ 1][\mathbf{x}^{\top}\ 1] and [x⊤ 1]⊤[\mathbf{x}^{\top}\ 1]^{\top}, the basis vector in (28), we have

The first two quadratic terms in (37) are nonnegative for any x0∈X0x_{0}\in\mathcal{X}_{0} by (35) and (36), respectively. Consequently, the last quadratic term must be nonpositive for all x0∈X0x_{0}\in\mathcal{X}_{0},

By Proposition 4, the above condition is equivalent to

for all yˉ∈{y∣y=fπ(x0), x0∈X0}=R(X0)\bar{y}\in\{y\mid y=f_{\pi}(x_{0}),\ x_{0}\in\mathcal{X}_{0}\}=\mathcal{R}(\mathcal{X}_{0}). Recall from (19) that the set of all points yy that satisfies (39) is the candidate set Rˉ(X0)\bar{\mathcal{R}}(\mathcal{X}_{0}). Therefore, we conclude that Rˉ(X0)\bar{\mathcal{R}}(\mathcal{X}_{0}) must be a superset of the exact one-step reachable set R(X0)\mathcal{R}(\mathcal{X}_{0}), i.e. R(X0)⊆Rˉ(X0)\mathcal{R}(\mathcal{X}_{0})\subseteq\bar{\mathcal{R}}(\mathcal{X}_{0}). ∎

3 Minimum-Volume Approximate Reachable Set

Theorem 1 provides a sufficient condition for over-approximating the one-step reachable set R(X0)\mathcal{R}(\mathcal{X}_{0}). Now we can use this result to reformulate problem (14), which finds a minimum-volume approximate reachable set Rˉ(X0)\bar{\mathcal{R}}(\mathcal{X}_{0}).

If the approximate reachable set is parametrized by an ellipsoid as in (21), we can easily obtain a minimum-volume ellipsoid that encloses R(X0)\mathcal{R}(\mathcal{X}_{0}) by solving,

Note that (34) is not convex in AA and bb. Nonetheless we can find a convex constraint equivalent to (34) using Schur complement. See for detail.

In , an MILP-based approach is proposed for estimating forward reachable sets of dynamical systems in closed-loop with neural network controllers. If the facets are given with the same directions as the ones of the true reachable sets, then the estimated reachable sets are exact. However, this method only works for polytopic reachable sets and does not scale well with the size of the neural network, the volume, and the dimension of the initial set.

Numerical Experiments

In this section, we demonstrate our approach with two application examples. The controllers used to generate training data were implemented in YALMIP . All neural network controllers were trained with ReLU activation functions and the Adam algorithm in PyTorch. We used MATLAB, CVX and Mosek to solve the Reach-SDP problems.

We first consider a double integrator system

discretized with sampling time ts=1t_{s}=1s and subject to state and input constraints, A∁=×\mathcal{A}^{\complement}=\times and U=\mathcal{U}=, respectively. We implemented a standard linear MPC following with a prediction horizon NMPC=10N_{\text{MPC}}=10, weighting matrices Q=I2Q=I_{2}, R=1R=1, the terminal region O∞LQR\mathcal{O}^{LQR}_{\infty} and the terminal weighting matrix P∞P_{\infty} synthesized from the discrete-time algebraic Riccati equation. The MPC is designed as a stabilizing controller which steers the system to the origin while satisfying the constraints. We then use the MPC to generate 2420 samples of state and input pairs (x,πMPC(x))(x,\pi_{\text{MPC}}(x)) for learning. The neural network has 2 hidden layers with 10 and 5 neurons, respectively. Our goal is to verify if all initial states in X0=[2.5,3]×[−0.25,0.25]\mathcal{X}_{0}=[2.5,3]\times[-0.25,0.25] can reach the set G=[−0.25,0.25]×[−0.25,0.25]\mathcal{G}=[-0.25,0.25]\times[-0.25,0.25], a region near the origin, in N=6N=6 steps while avoiding A\mathcal{A} at all times. We computed the minimum-volume polytopic approximate reachable sets Rˉ1(X0),⋯ ,Rˉ6(X0)\bar{\mathcal{R}}^{1}(\mathcal{X}_{0}),\cdots,\bar{\mathcal{R}}^{6}(\mathcal{X}_{0}) using the Reach-SDP introduced in Section 4. The matrix that determines the direction of the facets is chosen as

As shown in Figure 2, our approach yielded a tight outer-approximation of the reachable sets and successfully verified the safety properties sought.

2 6D Quadrotor

In the second example, we apply Reach-SDP to the 6D quadrotor model in . In order to have a linear model as in (5), we rewrite the nonlinear quadrotor dynamics in as follows,

where gg is the gravitational acceleration, the state vector x=[px,py,pz,vx,vy,vz]⊤x=[p_{x},p_{y},p_{z},v_{x},v_{y},v_{z}]^{\top} include positions and velocities of the quadrotor in the 3D space and the control vector uu is a function of θ\theta (pitch), ϕ\phi (roll) and τ\tau (thrust), which are the control inputs of the original model. The control task is to steer the quadrotor to the origin while respecting the state constraints A∁=×××××\mathcal{A}^{\complement}=\times\times\times\times\times and the actuator constraints [θ,ϕ,τ]⊤∈[−π/9,π/9]×[−π/9,π/9]×[0,2g][\theta,\phi,\tau]^{\top}\in[-\pi/9,\pi/9]\times[-\pi/9,\pi/9]\times[0,2g]. We implemented a nonlinear MPC with a prediction horizon NMPC=30N_{\text{MPC}}=30, a least-squares objective function that penalizes both states and inputs with weighting matrices Q=I6Q=I_{6}, R=I4R=I_{4} and the terminal constraint xNMPC=0x_{N_{\text{MPC}}}=0. We used the original nonlinear dynamics in as the prediction model in MPC, which is discretized with a sampling time ts=0.1t_{s}=0.1s using the Runge–Kutta 4th order method. The nonlinear MPC problems were solved using SNOPT . A total number of 4531 feasible samples of state and input pairs (x,πMPC(x))(x,\pi_{\text{MPC}}(x)) were generated and used to train a neural network with 2 hidden layers and 32 neurons in each layer. The initial set is given as an ellipsoid X0=E(q0,Q0)\mathcal{X}_{0}=\mathcal{E}(q_{0},Q_{0}), where q0=[4.7 4.7 3 0.95 0 0]⊤q_{0}=[4.7\ 4.7\ 3\ 0.95\ 0\ 0]^{\top} is the center and Q0=diag⁡(0.052,0.052,0.052,0.012,0.012,0.012)Q_{0}=\operatorname{diag}(0.05^{2},0.05^{2},0.05^{2},0.01^{2},0.01^{2},0.01^{2}) is the shape matrix. Here, we want to verify if all initial states in X0\mathcal{X}_{0} can reach the set G=[3.7,4.1]×[2.5,3.5]×[1.2,2.6]\mathcal{G}=[3.7,4.1]\times[2.5,3.5]\times[1.2,2.6], which is defined in the (px,py,pz)(p_{x},p_{y},p_{z})-space, in t=1t=1 second subject to the state and input constraints. We approximated the ellipsoidal forward reachable sets Rˉ1(X0),⋯ ,Rˉ10(X0)\bar{\mathcal{R}}^{1}(\mathcal{X}_{0}),\cdots,\bar{\mathcal{R}}^{10}(\mathcal{X}_{0}) using the Reach-SDP. The resulting approximate reachable sets are plotted in Figure 3 and 4 which shows that our method is able to verify the given reach-avoid specifications.

Conclusions

In this paper, we propose the first convex-optimization-based reachability analysis method for linear systems in feedback interconnection with neural network controllers. Our approach relies on abstracting the nonlinear components of the closed-loop system by quadratic constraints. Then we show that we can compute the approximate reachable sets via semidefinite programming. Future work includes extending the current approach to incorporate nonlinear dynamics and to approximate backward reachable sets, which is useful for certifying invariance properties.

Appendix A Proof of Proposition 5

Assuming the same activation function for all neurons throughout the entire network, we can write (22) compactly as

By Lemma 2 the neural network (22) satisfies the QC defined by Mmid\mathcal{M}_{\rm{mid}} with basis [x⊤ 1]⊤[\mathbf{x}^{\top}\ 1]^{\top}, which concludes the proof.

References