# Data-Driven Permissible Safe Control with Barrier Certificates

Rayan Mazouz∗1, John Skovbekk∗1, Frederik Baymler Mathiesen2, Eric Frew1,
Luca Laurenti2, and Morteza Lahijanian1 ∗Equal contributions.1Authors are with
the University of Colorado Boulder.2Authors are with Delft University of
Technology.

###### Abstract

This paper introduces a method of identifying a maximal set of safe strategies
from data for stochastic systems with unknown dynamics using barrier
certificates. The first step is learning the dynamics of the system via
Gaussian process (GP) regression and obtaining probabilistic errors for this
estimate. Then, we develop an algorithm for constructing piecewise stochastic
barrier functions to find a maximal permissible strategy set using the learned
GP model, which is based on sequentially pruning the worst controls until a
maximal set is identified. The permissible strategies are guaranteed to
maintain probabilistic safety for the true system. This is especially
important for learning-enabled systems, because a rich strategy space enables
additional data collection and complex behaviors while remaining safe. Case
studies on linear and nonlinear systems demonstrate that increasing the size
of the dataset for learning the system grows the permissible strategy set.

##  I Introduction

In an era defined by the increasing integration of Artificial Intelligence
(AI) into systems, ensuring the safety of stochastic systems with black-box
components has become a major concern. These systems, characterized by
uncertain and unpredictable dynamics, are ubiquitous across various domains,
from autonomous vehicles [1] to surgical robotics [2]. For such systems, the
notion of the _safety invariant set_ , which represents regions of the
system’s state space where non-myopic safety constraints are guaranteed to be
satisfied, emerges as a fundamental concept in providing safety guarantees.
Particularly, identifying the set of _permissible control strategies_ , within
these invariant sets holds critical importance, as it provides a comprehensive
understanding of the system’s operational boundaries and enables complex
behaviors while guaranteeing safety. However, computation of such a strategy
set poses a major challenge, especially if the system is unknown due to, e.g.,
AI or black-box components. In this work, we aim to provide a data-driven
method for identifying a maximal set of permissible strategies that guarantee
an unknown system remains inside a safe set.

Our approach is based on stochastic barrier functions (SBFs) [3, 4]. These
functions provide a systematic method to bounding the system’s behavior within
a safe set, even in the presence of uncertainties or disturbances. We also
utilize Gaussian process (GP) regression to learn the unknown dynamics, which
enable probabilistic bounding of the learning error [5, 6, 7]. Our key insight
in dealing with the computational challenge of identifying permissible
strategies is that local treatment of the uncertainty (namely, stochasticity
and learning error) is needed, and a particular form of functions that enables
such a treatment is piecewice (PW). Hence, we employ the recent results in PW-
SBFs [8] to reason about subsets of controls. Particularly, we show a
formulation of PW-SBFs that allow admissibility assessment of local control
subsets for each compact set of states via linear programming, achieving
computational efficiency. Based on these results, we propose an algorithm
that, given a set of input-output data on the dynamics, safe and initial sets,
and a lower-bound on the safety probability, returns a maximal permissible
strategies set by iteratively removing inadmissible control sets until a
fixed-point is reached. We show soundness of the algorithm. Our evaluations
validate the theoretical guarantees. They also illustrate that the approach is
general and works on both unknown systems with linear and nonlinear stochastic
dynamics.

In short, this work makes the following contributions:

  * •

a data-driven framework for computation of a maximal set of permissible
strategies,

  * •

extension of SBFs to provide safety control invariant sets for GP regressed
models on continuous control sets,

  * •

and a series of case studies on both linear and nonlinear systems that
demonstrate the efficacy of the method and validate the theoretical
guarantees.

###  I-A Related Work

Probabilistic safety invariance is an essential property to enforce the safety
of stochastic systems [9, 10, 11, 4]. For Markovian systems, these invariant
sets can be found using linear programs (LPs) and mixed integer LPs for finite
and infinite horizon probabilistic invariance, respectively [10]. The sound
application of these program requires knowledge of the transition kernel,
which is not available for unknown systems as in our setting. While a
sampling-based procedure for probabilistic invariant sets of deterministic,
stable black-box systems is available [9], it does not admit controlled
stochastic systems that may be inherently unstable.

Gaussian process (GP) regression is notable for its flexibility in learning
unknown systems and quantifying the uncertainty in safety-critical settings
[5, 12, 13, 11]. In addition to learning and estimating safe and stable
regions for unknown systems [5, 12], they have been used to learn state-
dependent modelling uncertainties for LTI control systems to subsequently find
probabilistic invariant sets under a fixed feedback controller [11]. However,
none of the previous can identify the full set of controls that leads to
safety invariance over a prescribed horizon.

Stochastic Barriers Functions (SBFs) can provide bounds on the probability of
remaining safe from an initial set [3, 4], even in the case a system is
learned from data [14, 15]. Synthesis based on GP regression over states and
controls can provide a single strategy with a lower-bound of safety [15]. In
addition, SBFs can be used to find invariant sets for uncontrolled systems
with bounded disturbances over an infinite horizon [4]. Commonly, SBF
synthesis relies on sum-of-squares optimization [3, 16] or neural barrier
training [17, 18]. While both methods are deemed useful, they each have their
own limitations. A novel and more efficient formulation relies on piecewise
(PW) SBFs, overcoming some of these major limitations [8]. In this work, we
adapt and extend PW-SBFs to enable computation of permissible controls for GP
regressed models.

##  II Problem Formulation

Consider a stochastic dynamical system with dynamics

$\begin{split}&\mathbf{x}_{k+1}=f(\mathbf{x}_{k},\mathbf{u}_{k})+\mathbf{w},\\\
k\in\mathbb{N}&,\quad\mathbf{x}_{k}\in\mathbb{R}^{n},\quad\mathbf{u}_{k}\in
U\subset\mathbb{R}^{m},\end{split}$ (1)

where $U$ is a compact set of controls, $f:\mathbb{R}^{n}\times
U\to\mathbb{R}^{n}$ is the vector field, and $\mathbf{w}$ represents noise
taking values in $W\subseteq\mathbb{R}^{n}$ and is assumed to be an
independent and identically distributed (i.i.d.) sub-Gaussian random variable
with known probability density function $p_{\mathbf{w}}$. The choice of the
control $\mathbf{u}_{k}\in U$ is performed by a stationary (feedback) control
strategy $\pi:\mathbb{R}^{n}\to U$. We denote the set of all control
strategies by $\Pi$.

In this paper, we assume that vector field $f$ is _unknown_. Instead, we
consider a given dataset $D=\\{(x_{i},u_{i},x^{\prime}_{i})\\}_{i=1}^{M}$,
consisting of $M$ input-output i.i.d. samples from System (LABEL:eq:system),
such that $x^{\prime}_{i}=f(x_{i},u_{i})+w_{i}$, where $w_{i}$ is a
realization of noise $\mathbf{w}$. We aim to learn $f$ using $D$. To that end,
we impose the following assumption, which ensures $f$ is a well-behaved
analytical function on a compact set.

######  Assumption 1 (RKHS Continuity).

Let a compact set $Y\subset\mathbb{R}^{n}\times U$ and
$\kappa:\mathbb{R}^{n+m}\times\mathbb{R}^{n+m}\to\mathbb{R}$ be a given
kernel. Further, define $H_{\kappa}(Y)$ as the reproducing kernel Hilbert
space (RKHS) of functions over $Y$ corresponding to $\kappa$ with norm
$\|\cdot\|\kappa$ [19]. Then, for each $i\in\\{1,\ldots,n\\}$, $f^{i}\in
H_{\kappa}(Y)$ and for a constant $C_{i}>0$, $\|f^{i}\|_{\kappa}\leq C_{i}$,
where $f^{i}$ is the $i$-th component of $f$.

Assumption 1 is standard in the literature and guarantees that $f$ can be
learned via Gaussian processes (GPs) with an appropriate kernel. For instance,
in the case where $\kappa$ is a universal kernel, such as the squared
exponential kernel, it includes a class of functions that is dense w.r.t. the
continuous functions over a compact set [20, 21].

Given these sources of uncertainty, we focus on probabilistic analysis of
System (LABEL:eq:system). We begin by defining a probability measure over the
trajectories of the system. For a measurable set $X\subset\mathbb{R}^{n}$, the
one-step transition kernel under a given strategy $\pi\in\Pi$ is

$T(X\mid
x,\pi):=\int_{\mathbb{R}^{\mathrm{w}}}\mathbf{1}_{X}(f(x,\pi(x))+w)p_{\mathbf{w}}(w)dw,$
(2)

where $\mathbf{1}_{X}$ is the indicator function for $X$ such that
$\mathbf{1}_{X}(x)=1$ if $x\in X$ and $\mathbf{1}_{X}(x)=0$ if $x\not\in X$
[22]. Then, for a time horizon $N\in\mathbb{N}$, a strategy $\pi,$ and an
initial condition $x_{0}\in\mathbb{R}^{n}$, we consider the probability
measure $\Pr$ over the trajectories of System (LABEL:eq:system) such that for
measurable sets $X_{0},X_{k}\subseteq X$, it follows that

$\displaystyle\Pr[\mathbf{x}_{0}\in X_{0}]=\mathbf{1}_{X_{0}}(x_{0}),$ (3)
$\displaystyle\Pr[\mathbf{x}_{k}\in
X_{k}\mid\mathbf{x}_{k-1}=x,\,\pi]=T(X_{k}\mid x,\pi).$

######  Definition 1 (Probabilistic Safety).

Given a measurable safe set $X_{\mathrm{s}}\subset\mathbb{R}^{n}$ and initial
set $X_{0}\subseteq X_{\mathrm{s}}$, the _probabilistic safety_ of System
(LABEL:eq:system) under a control strategy $\pi$ for $N$ time steps is defined
as

$P_{s}(X_{s},X_{0},N,\pi)=\\\ \inf_{x_{0}\in X_{0}}\Pr[\mathbf{x}_{k}\in
X_{s}\;\forall k\leq N\mid\mathbf{x}_{0}=x_{0},\pi].$

It is enough to find a single strategy that satisfies a threshold on $P_{s}$
to claim that a system is safe. However, identifying multiple strategies that
guarantee a minimum $P_{s}$ can unlock more complex behaviors for the system.
This is especially important for learning-enabled systems (e.g., the one
considered in this work) because under such strategies more data can be
collected in a safe manner to improve the learning model. To this end, we
define permissible strategies.

######  Definition 2 (Permissible strategy).

Given a safety threshold $p$, control strategy $\pi\in\Pi$ is called
_permissible_ if

$P_{s}(X_{s},X_{0},N,\pi)\geq p.$

Our goal is to find a maximal set of permissible strategies for System
(LABEL:eq:system).

######  Problem 1.

Given dataset $D$ obtained from System (LABEL:eq:system), safe set
$X_{s}\subset\mathbb{R}^{n}$, initial set $X_{0}\subseteq X_{s}$, time horizon
$N$, and safety probability threshold $p$, find a maximal set of permissible
strategies $\Pi_{s}\subseteq\Pi$ such that every $\pi\in\Pi_{s}$ is
permissible, i.e.,

$P_{s}(X_{s},X_{0},N,\pi)\geq p\qquad\forall\pi\in\Pi_{s}.$

Problem 1 poses several challenges. Firstly, it seeks a maximal set of
permissible strategies for an unknown System (LABEL:eq:system), which itself
is a stochastic process. Secondly, control set $U\subset\mathbb{R}^{m}$ is
continuous, making admissibility analysis of $\pi$ difficult due to its range
being an uncountable, infinite set. Furthermore, set $X_{s}$ can be non-
convex, possibly leading to an overall non-convex optimization problem.

######  Remark 1.

Note that a method that solves Problem 1 can also compute the probabilistic
control invariant set. This set consists of states from which, if the system
is initialized in it, it is possible to control the system to remain within
the set under probabilistic guarantees. The advantage of Problem 1 is that it
enables the control invariant set to be coupled with a maximal set of
permissible strategies.

Approach overview.  Our approach to Problem 1 is based on GP regression and
stochastic barrier functions (SBFs). GP regression allows learning of the
dynamics with formal learning-error quantification, and SFBs provide a formal
methodology to prove invariance of stochastic systems. To address issues
pertaining to continuity of $U$ and non-convexity of $X_{\mathrm{s}}$, we
utilize the recently-developed techniques for piecewise constant (PWC) SBFs
[8]. We show that based on these functions, an inner-approximation of the set
of permissible strategies can be obtained for System (LABEL:eq:system).

##  III Preliminaries

In this section, we provide a brief background on Gaussian process regression
and PWC-SBFs, which are core to our framework.

###  III-A Gaussian Process Regression

A Gaussian process (GP) is a collection of random variables $\mathrm{X}$, any
finite subset of which are jointly Gaussian [21]. Here, we describe GPs of
single-dimension output but note that multidimensional outputs are similar.
GPs are often interpreted as distributions over a function space, which are
defined completely by a mean function $\hat{f}:\mathrm{X}\to\mathbb{R}$ and
covariance (kernel) function $\kappa:\mathrm{X}\times\mathrm{X}\to\mathbb{R}$.
GP regression is a Bayesian approach to conditioning a prior GP model
$(\hat{f}_{0},\kappa_{0})$ on a dataset $D=(\mathrm{X},\mathrm{Y})$ to find a
posterior model $(\hat{f}_{D},\kappa_{D})$ defined as

$\displaystyle\hat{f}_{D}(x)$
$\displaystyle=\hat{f}_{0}(x)+\kappa_{0}(x,\mathrm{X})K^{-1}(\mathrm{Y}-\hat{f}_{0}(\mathrm{X})),$
$\displaystyle\kappa_{D}(x,x^{\prime})$
$\displaystyle=\kappa_{0}(x,x^{\prime})-\kappa_{0}(x,\mathrm{X})K^{-1}\kappa_{0}(\mathrm{X},x^{\prime}),$

where
$\kappa_{0}(\mathrm{X},x)=[\kappa_{0}(\mathrm{X}_{i},x),\dots,\kappa_{0}(\mathrm{X}_{M},x)]^{T}$
(similarly defined vectors of both arguments),
$K=\kappa_{0}(\mathrm{X},\mathrm{X})+\sigma^{2}I$, and $\sigma^{2}$ is the
noise variance on observations $\mathrm{Y}$. This posterior is a distribution
of functions that best describes the data $D$. When the input points to the
kernel are identical, $\kappa^{-1/2}_{D}(x,x)=\sigma_{D}(x)$.

###  III-B Piecewise Constant Stochastic Barrier Functions

Consider a stochastic process $\mathbf{z}_{k+1}=F(\mathbf{z}_{k},\mathbf{v}),$
where state $\mathbf{z}_{k}\in\mathbb{R}^{n_{z}}$, noise
$\mathbf{v}\in\mathbb{R}^{v}$ is a random variable with distribution
$p_{\mathbf{v}}$, and vector field
$F:\mathbb{R}^{n_{z}}\times\mathbb{R}^{v}\to\mathbb{R}^{n_{z}}$ is almost
everywhere continuous. For a measurable compact set
$Z\subset\mathbb{R}^{n_{z}}$, let $T_{z}(Z\mid z)$ and Pr be the stochastic
transition kernel and probability measure for process $\mathbf{z}_{k}$ defined
similarly to (2) and (3), respectively.

Further, consider partition $Z_{1},...,Z_{K_{z}}$ of $Z$ such that $F$ is
continuous over each $Z_{i}$, and for scalars $b_{i}\in\mathbb{R}_{\geq 0}$
with $i=\\{1,...,K_{z}\\}$, define PWC function
$B:\mathbb{R}^{n_{z}}\to\mathbb{R}_{\geq 0}$ as

$B(z)=\begin{cases}b_{i}&\text{if }z\in Z_{i}\\\
1&\text{otherwise.}\end{cases}$ (4)

The following theorem establishes the conditions under which $B$ is a _PWC-
SBF_ for process $\mathbf{z}$.

######  Theorem 1 ([8, Theorem 2]).

PWC function $B(z)$ in (4) is a stochastic barrier function for process
$\mathbf{z}$ with safe sets $Z\subset\mathbb{R}^{n_{z}}$ and initial set
$Z_{0}\subseteq Z$, if $\forall i\in\\{1,\ldots,K_{z}\\}$, there exist scalars
$\beta,\eta\geq 0$ such that

$\displaystyle b_{i}\leq\eta\qquad\qquad\qquad\qquad\qquad\qquad\forall
i:Z_{i}\cap Z_{0},$ (5a) $\displaystyle\sum_{j=1}^{K_{z}}b_{j}\cdot
T_{z}(Z_{j}\mid z)+T_{z}(\mathbb{R}^{n_{z}}\setminus Z\mid z)\leq$
$\displaystyle\qquad\qquad\qquad\qquad\qquad b_{i}+\beta\qquad\forall
i,\,\forall z\in Z_{i}.$ (5b)

Then, for a $N\in\mathbb{N}_{\geq 0}$, it follows that

$\displaystyle\Pr[\forall k\leq N,\;\mathbf{z}_{k}\in Z\mid\mathbf{z}_{0}\in
Z_{0}]\geq 1-(\eta+\beta N).$ (6)

The advantage of SBFs is that, by satisfying static conditions, they allow for
probabilistic reasoning about stochastic dynamical systems without having to
unfold their trajectories.

##  IV Piecewise Stochastic Barriers for Invariance via GP Regression

In this section, we introduce our main approach. First, we show that how GP
regression can be used to bound transition kernel $T$. Then, we propose a
method to use the kernel bounds with Theorem 1 to identify permissible
strategy sets with PWC-SBFs and provide soundness guarantees.

###  IV-A Learning Transition Kernel Bounds

#### Learning $f$

To begin, we learn function $f$ in System (LABEL:eq:system) using GP
regression on the state and control spaces using dataset $D$. For brevity, let

$z=(x,u)\text{ \quad with \quad}z\in Z=\mathbb{R}^{n}\times U.$

Then, we estimate each component of $f$ with a GP map
$\hat{f}_{D}^{i}:Z\to\mathbb{R}$ for $i=\\{1,\ldots,n\\}$. Assumption 1 allows
for uniform probabilistic error bounds of $\hat{f}_{D}^{i}$ [6]. Namely, for
scalar $\epsilon^{i}\geq 0$ for $1\leq i\leq n$ and compact set
$\tilde{Z}\subset Z$,

$\Pr\left[|f^{i}(z)-\hat{f}_{D}^{i}(z)|\leq\epsilon^{i}\right]\geq
1-\delta\qquad\forall z\in\tilde{Z},$ (7)

where $\delta$ satisfies $\epsilon^{i}=\alpha(\delta)\sigma^{i}_{D}(z)$ as
defined in [6, Theorem 2]. The term $\alpha(\delta)$ scales the posterior
according to the desired probability threshold and the RKHS properties of the
function $f$, including its RKHS norm that can be bounded using its Lipschitz
continuity [13]. This probabilistic error bound is the foundation for
computing bounds on the transition kernel $T$.

#### Bounding the Transition Kernel

Directly computing the transition kernel $T$ for every state in a region in
Theorem 1 is computationally intractable. Instead, we aim to compute bounds of
$T$ using $\hat{f}_{D}^{i}$ by employing techniques from [7, 13]. Consider a
partition $X_{1},\ldots,X_{K}$ of safe set $X_{\mathrm{s}}$ in $K$ compact
sets, and partition $U_{1},\ldots,U_{L}$ of control space $U$ in $L$ compact
sets, i.e.,

$\bigcup_{i=1}^{K}X_{i}=X_{\mathrm{s}}\quad\text{and}\quad\bigcup_{l=1}^{L}U_{l}=U$

all of which are non-overlapping. The product of these partitions is a
partition of the state-control space, i.e.,

$Z_{il}=X_{i}\times U_{l}.$

Then, the transition kernel from a state-control partition $Z_{il}$ to a
compact set $X_{j}\subset\mathbb{R}^{n}$ of states can be bounded by finding
the extrema of the transition kernel,

$\displaystyle\underline{p}^{l}_{ij}\leq\min_{x,u\in Z_{il}}T(X_{j}\mid x,u),$
(8) $\displaystyle\overline{p}^{l}_{ij}\geq\max_{x,u\in Z_{il}}T(X_{j}\mid
x,u).$

To compute bounds on these extrema, we partition the domains of the
uncertainty, i.e., those of the additive noise $W$ and the learning error
$\mathbb{R}_{\geq 0}^{n}$ in (7) [7]. This method seeks all the uncertainty
values to guarantee entering (avoiding) the target set to compute a non-
trivial lower (upper) bound of the extrema. Let $Q=W\times\mathbb{R}_{\geq
0}^{n}$ be the product space of the noise and learning error domains, and let
$\mathcal{C}=\\{C_{1},\dots,C_{p}\\}$ be non-overlapping partitions of $Q$.
Also, let the learned image, $Post(Z_{il},C)$, be all of the points in
$Z_{il}$ propagated through $\hat{f}_{D}$ with all uncertainties in
$C\in\mathcal{C}$ accounted for, i.e.,

$Post(Z_{il},C)=\\{\hat{f}_{D}(z)+c\mid z\in Z_{il},c\in C\\}.$

Finally, the extrema of $T$ are bounded using the following proposition.

######  Proposition 1 ([7, Theorem 1]).

Let $\mathcal{C}$ be a partition of the noise space and learning error domain
$Q=W\times\mathbb{R}^{n}_{\geq 0}$ of System (LABEL:eq:system), and $Z_{il}$
be a compact subset of $Z$. Then, the one-step transition kernel to a compact
set $X\subset\mathbb{R}^{n}$ from $Z_{il}$ is bounded by

$\displaystyle\min_{x,u\in Z_{il}}T(X\mid x,u)$
$\displaystyle\geq\sum_{C\in\mathcal{C}}\mathbf{1}(Post(Z_{il},C)\subseteq
X)\Pr(C),$ $\displaystyle\max_{x,u\in Z_{il}}T(X\mid x,u)$
$\displaystyle\leq\sum_{C\in\mathcal{C}}\mathbf{1}(Post(Z_{il},C)\cap
X=\emptyset)\Pr(C),$

where $\mathbf{1}(\cdot)$ returns one if the argument is true and zero
otherwise.

Proposition 1 bounds the transition kernel for any choice of uncertainty
partition $\mathcal{C}$.

######  Remark 2.

As the uncertainty structure of System (LABEL:eq:system) is relatively simple
(it is additive for both the noise and learning error), an optimal partition
$\mathcal{C}$ can be found that gives rise to the tightest bounds of $T$ in a
straightforward manner. In fact, computation of the tightest-possible bounds
requires only three partitions of the uncertainty space in each dimension. See
[7] for more details.

###  IV-B Permissible Strategy Set Synthesis

Using PWC-SBFs and the bounds on $T$, we extend Theorem 1 to find the maximal
permissible strategy set $\Pi_{s}$ for System (LABEL:eq:system). Let
$X_{\mathrm{u}}=\mathbb{R}^{n}\setminus X_{\mathrm{s}}$, and consider the
lower and upper bounds on the transition kernel $\underline{p}^{l}_{ij}$ and
$\overline{p}^{l}_{ij}$ in (8). Note that $\underline{p}^{l}_{ij}$ and
$\overline{p}^{l}_{ij}$ can be computed via Proposition 1. Then, the set of
all feasible values for $T(X_{j}\mid x,u)$ for all $x\in X_{i}$ and $u\in
U_{l}$ is

$\displaystyle\mathcal{P}^{l}_{i}=\Big{\\{}$ $\displaystyle
p^{l}_{i}=(p^{l}_{i1},\ldots,p^{l}_{iK},p^{l}_{i\mathrm{u}})\in[0,1]^{K+1}\quad
s.t.$ $\displaystyle\qquad\sum_{j=1}^{K}p^{l}_{ij}+p^{l}_{i\mathrm{u}}=1,$
$\displaystyle\qquad\underline{p}^{l}_{ij}\leq
p^{l}_{ij}\leq\overline{p}^{l}_{ij}\;\;\;\forall
j\in\\{1,\ldots,K,\mathrm{u}\\}\Big{\\}}.$ (9)

The following theorem sets up an optimization problem for synthesizing
permissible strategy set $\Pi_{s}$ using PWC-SBF.

######  Theorem 2 (PWC-SBF for Permissible Strategies).

Consider System (LABEL:eq:system) with safe set
$X_{\mathrm{s}}\subset\mathbb{R}^{n}$, initial set $X_{0}\subseteq
X_{\mathrm{s}}$, dataset $D$, and safety probability bound $p$. Given
$K$-partitions of $X_{\mathrm{s}}$, let the set of PWC functions
$\mathcal{B}_{K}$ be in the form (4), with $B_{i}(x)=b_{i}\in\mathbb{R}_{\geq
0}$, for every $i\in\\{1,\ldots,K\\}$. Further, given $L$-partitions of $U$,
let $\mathcal{P}^{l}=\prod_{i=1}^{K}\mathcal{P}_{i}^{l}$, where
$\mathcal{P}_{i}^{l}$ is the set of feasible transition kernel in (9), with
the bounds computed by Proposition 1, for every $l\in\\{1,\ldots,L\\}$.
Finally, let $B^{*}\in\mathcal{B}_{K}$ be a solution to the following
optimization problem

$\displaystyle
B^{*}=\arg\min_{B\in\mathcal{B}_{K}}\;\max_{\\{\mathcal{P}^{l}\\}_{l=1}^{L}}\;\eta+N\beta$
(10)

subject to

$\displaystyle b_{i}\leq\eta$ $\displaystyle\forall i:X_{i}\cap
X_{0}\neq\emptyset,$ (11a) $\displaystyle\sum_{j=1}^{K}b_{j}\cdot
p^{l}_{ij}+p^{l}_{i\mathrm{u}}\leq b_{i}+\beta_{i}^{l}$ $\displaystyle\forall
i,\leavevmode\nobreak\ \forall l,\leavevmode\nobreak\ \forall
p_{i}^{l}\in\mathcal{P}^{l}_{i},$ (11b) $\displaystyle
0\leq\beta_{i}^{l}\leq\beta,$ $\displaystyle\forall i,\leavevmode\nobreak\
\forall l.$ (11c)

Then,

  1. 1.

$B^{*}$ constitutes a stochastic barrier certificate for System
(LABEL:eq:system) that guarantees safety probability
$P_{s}(X_{s},X_{0},N,\pi)\geq 1-(\eta+\beta N)$ for every strategy
$\pi\in\Pi$;

  2. 2.

if $\beta_{i}^{l}\leq(1-\eta-p)/N$, then every choice of $u\in U_{l}$ is
permissible for every $x\in X_{i}$ of System (LABEL:eq:system).

Proof Sketch. The optimization is an exact solution to the problem in Theorem
1. Hence, it is straightforward that $B^{*}$ is a proper stochastic barrier
certificate, and that for every strategy $\pi\in\Pi$,
$P_{s}(X_{s},X_{0},N,\pi)\geq 1-(\eta+\beta N)$, where
$\beta=§\max\beta_{i}^{l}$. Further, it follows that if
$\beta_{i}^{l}\leq(1-\eta-p)/N$, state-control partition $Z_{il}$ is rendered
safe $\forall(x,u)\in X_{i}\times U_{l}$. Since the transition kernel bounds
of the learned system in (9) are contained in the distribution of the true
system [13, Theorem 1], the probabilistic safety guarantees hold for System
(LABEL:eq:system). $\square$

The above theorem provides a method of identifying permissible strategies,
which requires solving a minimax optimization problem. This problem includes
bilinear terms, i.e., $b_{j}\cdot p^{l}_{ij}$ in Condition (11b).
Nevertheless, by introducing dual variables, the minimax problem becomes a
simple linear program (LP) as shown in [8], which guarantees efficiency.
Furthermore, following the same approach developed in [8], the resulting LP
can be efficiently solved by a Counter Example Guided Synthesis (CEGS)
approach or a gradient descent method. Our evaluations in Section V use the
CEGS approach.

1

2

Input : Initial set $X_{0}$, state-control space $\mathcal{Z}$, time horizon
$N$, feasible transition kernel sets
$\tilde{\mathcal{P}}=\\{\mathcal{P}^{l}\\}_{l=1}^{L}$, and probability
threshold $p$.

Output : Permissible strategy set $\Pi_{s}$.

3

4$\eta,\beta\leftarrow\textsc{Barrier}(X_{0},\mathcal{Z},N,\tilde{\mathcal{P}}$)
$\triangleright$ Theorem 2

5

6while  _$1-(\eta+N\beta) <p$ _ do

7

8 $Z_{il}^{\prime},\mathcal{P}_{i}^{l}\leftarrow$ Identify region with
$\max\beta_{i}^{l}$

9

10 $\mathcal{Z}\leftarrow$ $\mathcal{Z}\setminus Z^{\prime}_{il}$
$\triangleright$ Remove worst region

11

12 $\tilde{\mathcal{P}}\leftarrow$
$\tilde{\mathcal{P}}\setminus\mathcal{P}_{i}^{l}$ $\triangleright$ Remove
kernel bounds

13

14
$\eta,\beta\leftarrow\textsc{Barrier}(X_{0},\mathcal{Z},N,\tilde{\mathcal{P}}$)
$\triangleright$ Theorem 2

15

16 if  _all control sets removed for $X_{i}$ _ then

17 return False

18 $\Pi_{s}\leftarrow\\{X_{i}\to
U_{l}\mid\forall(X_{i},U_{l})\in\mathcal{Z}\\}$

return $\Pi_{s}$

Algorithm 1 PWC-SBF based Control Invariant Sets

From the results of Theorem 2, we propose an algorithm that computes a maximal
set of permissible strategies $\Pi_{s}\subseteq\Pi$ w.r.t. to partitions of
the safe set and control space. An overview of the algorithm is shown in Alg.
1. For convenience, define set $\mathcal{Z}$ as

$\displaystyle\mathcal{Z}=\\{Z_{il}\mid Z_{il}=X_{i}\times U_{l},\;\forall
i,\;\forall l\\}$

The algorithm takes as input the initial set $X_{0}$, set $\mathcal{Z}$, time
horizon $N$, the feasible bounds on the transition kernel $T$, and the desired
probability of safety $p$. On Line 1, the algorithm computes an initial PWC-
SBF using Theorem 2. Next, the while loop is instantiated, with a termination
condition that the safety probability criteria $p$ must be met.

Using the PWC-SBF, the worst state-control pair ${Z}^{\prime}_{il}$
corresponding to $\max\beta_{i}^{l}$, along with the probability bounds
${P}_{i}^{l}$, are identified on Line 1. These are subsequently removed from
their respective sets on Lines 1 and 1. Intuitively, this means that the
algorithm sequentially prunes control region $U_{l}$ that is deemed not
permissible for a given state partition $X_{i}$, corresponding to
$\max\beta_{i}^{l}$. Then, a new PWC-SBF is synthesized on Line (1) using
Theorem 2. A check procedure on Line (1) is incorporated to make sure that for
all regions $X_{i}$, there exists at least one control region $U_{l}$ inside
updated space $\mathcal{Z}$. If the latter is true, the maximal set of
permissible strategies $\Pi_{s}$ is updated accordingly on Line (1).
Otherwise, the algorithm terminates as False. A successful algorithm
termination (Line 10) implies that the probability threshold $p$ is met.

######  Corollary 1.

Alg. 1 terminates in finite time. If it terminates with the return statement
on Line 1, then it returns $\Pi_{s}$, which is guaranteed to include only
permissible strategies.

The proof is a direct implication of Theorem 2.

######  Remark 3.

Inline with Remark 1, note that Alg. 1 can also compute the probabilistic
control invariant set. This is established by running the algorithm for every
$X_{i}\in X_{\mathrm{s}}$, $i\in\\{1,\ldots,K\\}$ as the initial set.

##  V Case Studies

We demonstrate the efficacy of Alg. 1 by finding permissible strategy sets for
linear and nonlinear systems learned from data. For each system, we perform GP
regression to estimate the system and then find the permissible sets. Each GP
uses zero mean and squared exponential kernel priors. To determine the quality
of our solution, we also compute the permissible strategy set of the known
linear system on the same partitions.

In all experiments, the initial state region of interest is
$X_{0}=[0.4,0.5]^{2}$, the safe set is
$X_{\mathrm{s}}=[0.0\leavevmode\nobreak\ 1.0]^{2}$, and the noise $\mathbf{w}$
is Gaussian zero-mean with variance $0.01^{2}$ and zero covariance in each
dimension. Permissible sets are constructed with Alg. 1 using a convergence
threshold of $p=1-1\times 10^{-4}$. Then, this set is used to generate bounds
on the probability of safety for $N=100$ time-step trajectories initialized
from $X_{0}$.

The safety probability bounds are validated by Monte Carlo simulations
randomly initialized in $X_{0}$ and propagated by random sampling of the
controls from the permissible set. To further validate the guarantees of the
approach, adversarial trajectories are generated by choosing control inputs
that drive the system towards exiting the safe set from both the full control
set and the permissible set.

Our implementation of Alg. 1 is written in Julia, and all the computations
were performed on an Apple M2 Pro MacBook with 16GB of RAM.

(a) Known System - $P_{\mathrm{s}}\geq 1-5.83\times 10^{-5}$

(b) GP 500 - $P_{\mathrm{s}}\geq 1-2.39\times 10^{-4}$

(c) GP 2000 - $P_{\mathrm{s}}\geq 1-1.93\times 10^{-5}$

Figure 1: Results for the linear system case studies using the known system
and two datasets. Regions with # of removed control regions and simulated
trajectories over 100 steps are shown in each plot. Adversarial trajectories
show permissible sets maintain safety.

###  V-A Linear System

First, consider the 2D system
$\mathbf{x}_{k+1}=A\mathbf{x}_{k}+B\mathbf{u}_{k}+\mathbf{w}_{k}$, where
$A=0.5I$, $B=(1,1)^{T}$, and $u\in[0.0,0.5]$. From the initial set, this
system primarily traverses around the diagonal $x_{1}=x_{2}$ line. The safe
set is discretized using a uniform grid with length 0.1, and the control range
is partitioned into five intervals. Computing each permissible set takes
approximately 70-100 seconds for this system.

Figure 1 shows the safe and the initial sets. It also shows the final results
for the known and learned systems. The number inside each cell is the number
of removed control sets. In each case, we simulated 1000 trajectories using
randomly sampled permissible strategies, which resulted in zero violations.
These trajectories are shown in gray color. Observe that, while the system can
reach towards the edges of the safe set, the permissible strategy set ensures
the system remains inside the safe set.

To further demonstrate the efficacy of the permissible strategy set, we
compare trajectories that choose the most adversarial actions according to the
full control set (shown in dashed black) and the permissible strategy set
(shown in solid black). Note that without restricting the strategies, the
system can always be driven outside of the safe set.

Figures 1(b) and 1(c) show the result using a GP regression with 500 and 2000
datapoints sampled from the system respectively. For the smaller dataset, more
control partitions are removed from the permissible set, resulting in
trajectories that cover a smaller range of the safe set. The learned model
with 2000 data points results in an permissible set and samples that more
resembles the known system, as fewer control partitions are removed. The
probabilities of exiting the safe set in 100 time steps is very small, and the
learned system with 2000 datapoints has a tighter upper bound than the known
system shown in Figure 1(a). This is likely due to convergence-related
hyperparameters in the optimization procedure used to synthesize the PWC-SBF
in Alg. 1.

The richness of the permissible strategy set is indicated by the proportion of
the control space it contains. The permissible strategy set for a known system
encompasses 93.6% of the system’s total possible actions, serving as a
benchmark for learning-based methods. With 500 data points, the procedure
recognizes 88.8% of the entire action space, while incorporating up to 2000
data points enhances this figure to 91.2%. This indicates that as more data
becomes available to estimate the system, the permissible strategy set
progressively approximates that of the known system.

(a) GP 500 - $P_{\mathrm{s}}\geq 1-8.96\times 10^{-3}$

(b) GP 1000 - $P_{\mathrm{s}}\geq 1-9.10\times 10^{-3}$

(c) GP 1500 - $P_{\mathrm{s}}\geq 1-9.95\times 10^{-3}$

Figure 2: Results for the nonlinear system case studies using various sizes of
datasets. Regions with # of removed control regions and simulated trajectories
over 100 steps are shown in each plot. Adversarial trajectories show
permissible control sets maintain safety.

###  V-B Nonlinear System

Consider the nonlinear dynamical system given by
$\mathbf{x}_{k+1}=\mathbf{x}_{k}+0.2[\cos(\mathbf{u}_{k}),\leavevmode\nobreak\
\sin(\mathbf{u}_{k})]^{T}+\mathbf{w}_{k}$, where $u\in[-\pi,\pi]$. This is a
simplified version of a Dubin’s car model, where the heading is set directly,
allowing movement in any direction. The state-space partition is the same as
for the system in Section V-A, and the control range is partitioned into 20
intervals. The permissible set computation takes approximately 1 to 1.5 hours
for each dataset.

Figure 2 illustrates the comparison of permissible strategy sets obtained for
the system using datasets of 500, 1000, and 1500 data points. In constructing
the permissible sets, numerous control intervals that would lead the system to
approach the safe set boundary are excluded. Despite the variations in dataset
size, each permissible set constrains the probability of exiting the safe set
within a similar bound.

The resulting permissible sets for the system learned with 500 and 1000
datapoints maintained 39.5% and 40.1% of all available controls, respectively.
With an additional 500 datapoints, the proportion of controls available rises
to 49.5%.

Out of the 1000 sampled trajectories using random permissible strategies, none
exited the safe set, and only the adversarial trajectory using the full
control set violated safety. The strict bound on $\beta$ during the synthesis
of the permissible set removed a majority of the controls to ensure safety. If
the bound on $\beta$ were relaxed, then additional strategies could be added
to the permissible set.

##  VI Conclusion

In this work, we introduce a data-driven framework for the computation of a
maximal set of permissible strategies for ensuring the safety of unknown
stochastic systems. The unknown model with a continuous control space is
regressed from data using Gaussian processes. Then, the regressed model is
used to synthesize a piecewise stochastic barrier function, which in turn
identifies theoretically-sound permissible control sets. The safety guarantees
are validated in case studies of linear and nonlinear systems. Future work
includes the incorporation of other regression approaches, extensions to
general temporal logic specifications, and applications to online safe
exploration and learning.

## References

  * [1] S. Shalev-Shwartz, S. Shammah, and A. Shashua, “On a Formal Model of Safe and Scalable Self-driving Cars,” _arXiv preprint arXiv:1708.06374_ , 2017.
  * [2] J. Guiochet, M. Machin, and H. Waeselynck, “Safety-critical Advanced Robots: A survey,” _Robotics and Autonomous Systems_ , vol. 94, pp. 43–52, 2017\.
  * [3] C. Santoyo, M. Dutreix, and S. Coogan, “A barrier function approach to finite-time stochastic system verification and control,” _Automatica_ , vol. 125, p. 109439, 2021.
  * [4] Y. Yu, T. Wu, B. Xia, J. Wang, and B. Xue, “Safe probabilistic invariance verification for stochastic discrete-time dynamical systems,” in _2023 62nd IEEE Conference on Decision and Control (CDC)_. IEEE, 2023, pp. 5804–5811.
  * [5] F. Berkenkamp, M. Turchetta, A. Schoellig, and A. Krause, “Safe model-based reinforcement learning with stability guarantees,” _Advances in neural information processing systems_ , vol. 30, 2017.
  * [6] S. R. Chowdhury and A. Gopalan, “On kernelized multi-armed bandits,” in _International Conference on Machine Learning_. PMLR, 2017, pp. 844–853.
  * [7] J. Skovbekk, L. Laurenti, E. Frew, and M. Lahijanian, “Formal abstraction of general stochastic systems via noise partitioning,” _IEEE Control Systems Letters_ , 2023.
  * [8] R. Mazouz, F. Baymler Mathiesen, L. Laurenti, and M. Lahijanian, “Piecewise Stochastic Barrier Functions,” _arXiv preprint arXiv:2404.16986_ , 2024\.
  * [9] Z. Wang and R. M. Jungers, “Data-driven computation of invariant sets of discrete time-invariant black-box systems,” _arXiv preprint arXiv:1907.12075_ , 2019.
  * [10] Y. Gao, K. H. Johansson, and L. Xie, “Computing probabilistic controlled invariant sets,” _IEEE Transactions on Automatic Control_ , vol. 66, no. 7, pp. 3138–3151, 2020.
  * [11] P. Griffioen, A. Devonport, and M. Arcak, “Probabilistic invariance for gaussian process state space models,” in _Learning for Dynamics and Control Conference_. PMLR, 2023, pp. 458–468.
  * [12] A. Lederer and S. Hirche, “Local asymptotic stability analysis and region of attraction estimation with gaussian processes,” in _2019 IEEE 58th Conference on Decision and Control (CDC)_. IEEE, 2019, pp. 1766–1771.
  * [13] J. Jackson, L. Laurenti, E. Frew, and M. Lahijanian, “Formal verification of unknown dynamical systems via gaussian process regression,” _arXiv preprint arXiv:2201.00655_ , 2021.
  * [14] P. Jagtap, S. Soudjani, and M. Zamani, “Formal synthesis of stochastic systems via control barrier certificates,” _IEEE Transactions on Automatic Control_ , vol. 66, no. 7, pp. 3097–3110, 2020.
  * [15] R. Wajid, A. U. Awan, and M. Zamani, “Formal synthesis of safety controllers for unknown stochastic control systems using gaussian process learning,” in _Learning for Dynamics and Control Conference_. PMLR, 2022, pp. 624–636.
  * [16] R. Mazouz, K. Muvvala, A. Ratheesh, L. Laurenti, and M. Lahijanian, “Safety Guarantees for Neural Network Dynamic Systems via Stochastic Barrier Functions,” _Advances in Neural Information Processing Systems_ , 2022.
  * [17] F. B. Mathiesen, S. C. Calvert, and L. Laurenti, “Safety Certification for Stochastic Systems via Neural Barrier Functions,” _IEEE Control Systems Letters_ , vol. 7, pp. 973–978, 2022.
  * [18] C. Dawson, S. Gao, and C. Fan, “Safe control with learned certificates: A survey of neural lyapunov, barrier, and contraction methods for robotics and control,” _IEEE Transactions on Robotics_ , 2023.
  * [19] N. Srinivas, A. Krause, S. M. Kakade, and M. W. Seeger, “Information-theoretic regret bounds for gaussian process optimization in the bandit setting,” _IEEE transactions on information theory_ , vol. 58, no. 5, pp. 3250–3265, 2012.
  * [20] I. Steinwart, “On the influence of the kernel on the consistency of support vector machines,” _Journal of machine learning research_ , vol. 2, no. Nov, pp. 67–93, 2001.
  * [21] C. E. Rasmussen and C. K. I. Williams, _Gaussian Processes for Machine Learning_. The MIT Press, 2006.
  * [22] D. P. Bertsekas and S. Shreve, _Stochastic optimal control: the discrete-time case_. Athena Scientific, 2004.