# Chance constrained sets approximation:
A probabilistic scaling approach - EXTENDED VERSION

M. Mammarella<EMAIL_ADDRESS>V. Mirasierra<EMAIL_ADDRESS>M. Lorenzen<EMAIL_ADDRESS>T. Alamo<EMAIL_ADDRESS>F. Dabbene
<EMAIL_ADDRESS>CNR-IEIIT; c/o Politecnico di Torino; C.so Duca
degli Abruzzi 24, Torino; Italy.  Universidad de Sevilla, Escuela Superior de
Ingenieros, Camino de los Descubrimientos s/n, Sevilla; Spain.
Systemwissenschaften, TTI GmbH, Nobelstr. 15, 70569 Stuttgart, Germany

###### Abstract

In this paper, a sample-based procedure for obtaining simple and computable
approximations of chance-contrained sets is proposed. The procedure allows to
control the complexity of the approximating set, by defining families of
simple-approximating sets of given complexity. A probabilistic scaling
procedure then allows to rescale these sets to obtain the desired
probabilistic guarantees. The proposed approach is shown to be applicable in
several problem in systems and control, such as the design of Stochastic Model
Predictive Control schemes or the solution of probabilistic set membership
estimation problems.

##  1 Introduction

In real-world applications, the complexity of the phenomena encountered and
the random nature of data makes dealing with uncertainty essential. In many
cases, uncertainty arises in the modeling phase, in some others it is
intrinsic to both the system and the operative environment, as for instance
wind speed and turbulence in aircraft or wind turbine control [1]. Hence, it
is crucial to include underlying stochastic characteristic of the framework
and eventually accept a violation of constraints with a certain probability
level, in order to improve the coherence of the model and reality. Deriving
results in the presence of uncertainty is of major relevance in different
areas, including, but not limited to, optimization [2] and robustness analysis
[3]. However, with respect to robust approaches, where the goal is to
determine a feasible solution which is optimal in some sense for all possible
uncertainty instances , the goal in the stochastic framework is to find a
solution that is feasible for almost all possible uncertainty realizations,
[4, 5]. In several applications, including engineering and finance, where
uncertainties in price, demand, supply, currency exchange rate, recycle and
feed rate, and demographic condition are common, it is acceptable, up to a
certain safe level, to relax the inherent conservativeness of robust
constraints enforcing probabilistic constraints. More recently, the method has
been used also in unmanned autonomous vehicle navigation [6, 7] as well as
optimal power flow [8, 9].

In the optimization framework, constraints involving stochastic parameters
that are required to be satisfied with a pre-specified probability threshold
are called chance constraints (CC). In general, dealing with CC implies facing
two serious challenges, that of stochasticity and of nonconvexity [10].
Consequently, while being attractive from a modeling viewpoint, problems
involving CC are often computationally intractable, generally shown to be NP-
hard, which seriously limits their applicability. However, being able to
efficiently solve CC problems remains an important challenge, especially in
systems and control, where CC often arise, as e.g. in stochastic model
predictive control (SMPC) [11, 12]. The scientific community has devoted large
research in devising computationally efficient approaches to deal with chance-
constraints. We review such techniques in Section 3, where we highlight three
mainstream approaches: i) exact techniques; ii) robust approximations and iii)
sample-based approximations . In this paper, we present what we consider an
important step forward in the sample-based approach. We propose a simple and
efficient strategy to obtain a probabilistically guaranteed inner
approximation of a chance constrained set, with given confidence.

In particular, we describe a two step procedure the involves: i) the
preliminary approximation of the chance constraint set by means of a so-called
Simple Approximating Set (SAS), ii) a sample-used scaling procedure that
allows to properly scale the SAS so to guarantee the desired probabilistic
properties. The proper selection of a low-complexity SAS allows the designer
to easily tune the complexity of the approximating set, significantly reducing
the sample complexity. We propose several candidate SAS shapes, grouped in two
classes: i) sampled-polytopes; and ii) norm-based SAS.

The probabilistic scaling approach was presented in the conference papers [13,
14]. The present work extends these in several directions: first, we performe
here a thorough mathematical analysis the results, providing of all results.
Second, the use of norm-based SAS is extended to comprise more general sets
(as e.g. , and More importantly, we consider here joint chance constraints.
This choice is motivated by the fact that enforcing joint chance constraints,
which have to be satisfied simultaneously, adheres better to some
applications, despite the inherent complexity. Finally, we present here a
second application, besides SMPC, related to probabilistic set-membership
identification.

The paper is structured as follows. Section 2 provides a general preamble of
the problem formulation and of chance constrained optimization, including two
motivating examples. An extensive overview on methods for approximating chance
constrained sets is reported in Section 3 whereas the probabilistic scaling
approach has been detailed in Section 4. Section 5 and Section 6 are dedicated
to the definition of selected candidate SAS, i.e. sampled-polytope and norm-
based SAS, respectively. Last, in Section 7, we validate the proposed approach
with a numerical example applying our method to a probabilistic set membership
estimation problem. Main conclusions and future research directions are
addressed in Section 8.

###  1.1 Notation

Given an integer $N$, $[N]$ denotes the integers from 1 to $N$. Given
$z\in\mathbb{R}^{s}$ and $p\in[1,\infty)$, we denote by $\|z\|_{p}$ the
$\ell_{p}$-norm of $z$, and by
$\mathbb{B}^{s}_{p}\doteq\\{\;z\in\mathbb{R}^{s}\;:\;\|z\|_{p}\leq 1\;\\}$
$\ell_{p}$-norm ball of radius one. Given integers $k,N$, and parameter
$p\in(0,1)$, the Binomial cumulative distribution function is denoted as

$\mathbf{B}(k;N,p)\doteq\sum\limits_{i=0}^{k}\left(\begin{array}[]{c}N\\\ i\\\
\end{array}\right)p^{i}(1-p)^{N-i}.$ (1)

The following notation is borrowed from the field of order statistics [15].
Given a set of $N$ scalars $\gamma_{i}\in\mathbb{R}^{N}$, $i\in[N]$, we denote
$\gamma_{1:N}$ the smallest one, $\gamma_{2:N}$ the second smallest one, and
so on and so forth until $\gamma_{N:N}$, which is equal to the largest one. In
this way, given $r\geq 0$ we have that $\gamma_{r+1:N}$ satisfies that no more
than $r$ elements of $\\{\gamma_{1},\gamma_{2},\ldots,\gamma_{N}\\}$ are
strictly smaller than $\gamma_{r+1:N}$.
The Chebyshev center of a given set $\mathbb{X}$, denoted as
$\mathsf{Cheb}(\mathbb{X})$, is defined as the center of the largest ball
inscribed in $\mathbb{X}$, i.e.

$\mathsf{Cheb}(\mathbb{X})\doteq\arg\min_{\theta_{c}}\max_{\theta\in\mathbb{X}}\left\\{\|\theta-\theta_{c}\|^{2}\right\\}.$

Given an $\ell_{p}$-norm $\|\cdot\|_{p}$, its dual norm $\|\cdot\|_{p^{*}}$ is
defined as

$\|c\|_{p^{*}}\doteq\sup\limits_{z\in\mathbb{B}^{s}_{p}}c^{\top}z,\;\forall
c\in\mathbb{R}^{s}.$

In particular, the couples $(p,p^{*})$: $(2,2)$, $(1,\infty)$, $(\infty,1)$
give raise to dual norms.

##  2 Problem formulation

Consider a robustness problem, in which the controller parameters and
auxiliary variables are parametrized by means of a decision variable vector
$\theta$, which is usually referred to as design parameter and is restricted
to a set $\Theta\subseteq\mathbb{R}^{n_{\theta}}$. Furthermore, the
uncertainty vector $w\in\mathbb{R}^{n_{w}}$ represents one of the admissible
uncertainty realizations of a random vector with given probability
distribution $\mathsf{Pr}_{\mathbb{W}}$ and (possibly unbounded) support
$\mathbb{W}$.

This paper deals with the special case where the design specifications can be
decoded as a set of $n_{\ell}$ uncertain linear inequalities

$F(w)\theta\leq g(w),$ (2)

where

$F(w)=\begin{bmatrix}f_{1}^{\top}(w)\\\ \vdots\\\
f_{n_{\ell}}^{\top}(w)\end{bmatrix}\in\mathbb{R}^{n_{\ell}\times{n_{\theta}}},\quad
g(w)=\begin{bmatrix}g_{1}(w)\\\ \vdots\\\
g_{n_{\ell}}(w)\end{bmatrix}\in\mathbb{R}^{n_{\ell}},$

are measurable functions of the uncertainty vector $w\in\mathbb{R}^{n_{w}}$.
The inequality in (2) is to be interpreted component-wise, i.e.

$f_{\ell}(w)\theta\leq g_{\ell}(w),\forall\ell\in[n_{\ell}].$

Furthermore, we notice that each value of $w$ gives raise to a corresponding
set

$\mathbb{X}(w)=\\{\;\theta\in\Theta\;:\;F(w)\theta\leq g(w)\;\\}.$ (3)

Due to the random nature of the uncertainty vector $w$, each realization of
$w$ corresponds to a different set of linear inequalities. Consequently, each
value of $w$ gives raise to a corresponding set

$\mathbb{X}(w)=\\{\;\theta\in\Theta\;:\;F(w)\theta\leq g(w)\;\\}.$ (4)

In every application, one usually accepts a risk of violating the constraints.
While this is often done by choosing the set $\mathbb{W}$ appropriately, we
can find a less conservative solution by choosing the set $\mathbb{W}$ to
encompass all possible values and characterizing the region of the design
space $\Theta$ in which the fraction of elements of $\mathbb{W}$, that violate
the constraints, is below a specified level. This concept is rigorously
formalized by means of the notion of _probability of violation_.

######  Definition 1 (Probability of violation).

Consider a probability measure ${\rm Pr}_{\mathbb{W}}$ over $\mathbb{W}$ and
let $\theta\in\Theta$ be given. The probability of violation of $\theta$
relative to inequality (2) is defined as

$\mathsf{Viol}(\theta)\doteq\mathsf{Pr}_{\mathbb{W}}\,\\{\,F(w)\theta\not\leq
g(w)\,\\}.$

Given a constraint on the probability of violation, i.e.
$\mathsf{Viol}(\theta)\leq\varepsilon$, we denote as (joint) _chance
constrained set_ of probability $\varepsilon$ (shortly, $\varepsilon$-CCS) the
region of the design space for which this probabilistic constraint is
satisfied. This is formally stated in the next definition.

######  Definition 2 ($\varepsilon$-CCS).

Given $\varepsilon\in(0,1)$, we define the chance constrained set of
probability $\varepsilon$ as follows

$\mathbb{X}_{\varepsilon}=\\{\;\theta\in\Theta\;:\;\mathsf{Viol}(\theta)\leq\varepsilon\;\\}.$
(5)

Note that the $\varepsilon$-CCS represents the region of the design space
$\Theta$ for which this probabilistic constraint is satisfied and it is
equivalently defined as

$\mathbb{X}_{\varepsilon}\doteq\Bigl{\\{}\theta\in\Theta\;:\;\mathsf{Pr}_{\mathbb{W}}\left\\{F(w)\theta\leq
g(w)\right\\}\geq 1-\varepsilon\Bigr{\\}}.$ (6)

######  Remark 1 (Joint vs. individual CCs).

The constraint $\theta\in\mathbb{X}_{\varepsilon}$, with
$\mathbb{X}_{\varepsilon}$ defined in (6), describes a joint chance
constraint. That is, it requires that the joint probability of satisfying the
inequality constraint

$F(w)\theta\leq g(w)$

is guaranteed to be greater than the probabilistic level $1-\varepsilon$. We
remark that this constraint is notably harder to impose than individual CCs,
i.e. constraints of the form

$\displaystyle\theta\in\mathbb{X}_{\varepsilon_{\ell}}^{\ell}\,\,$
$\displaystyle\\!\\!\\!\\!\doteq\\!\\!\\!$
$\displaystyle\Bigl{\\{}\theta\in\Theta\,:\,\mathsf{Pr}_{\mathbb{W}}\left\\{f_{\ell}(w)^{\top}\theta\leq
g_{\ell}(w)\right\\}\geq 1-\varepsilon_{\ell}\Bigr{\\}},$
$\displaystyle\qquad\ell\in[n_{\ell}],$

with $\varepsilon_{\ell}\in(0,1)$. A discussion on the differences and
implications of joint and individual chance constraints may be found in
several papers, see for instance [10, 16] and references therein.

######  Example 1.

A simple illustrating example of the set $\varepsilon$-CCS is shown in Figure
1. The dotted circle is the region of the design space that satisfies all the
constraints (the so called robust region), which are tangent to the dotted
circle at points uniformly generated. The outer red circle represents the
chance constrained set $\mathbb{X}_{\varepsilon}$ for the specific value
$\varepsilon=0.15$. That is, the red circle is obtained in such a way that
every point in it has a probability of violating a random constraint no larger
than $0.15$. Note that in this very simple case, the set
$\mathbb{X}_{\varepsilon}$ can be computed analytically, and turns out to be a
scaled version of the robust set. We observe that the $\varepsilon$-CCS is
significantly larger than the robust set.

Figure 1: Red circle = $\mathbb{X}_{\varepsilon}$, dotted circle = unit
circle, blue lines = constraint samples.

Hence, while there exist simple examples for which a closed-form computation
of $\mathbb{X}_{\varepsilon}$ is possible, as the one re-proposed here and
first used in [13], we remark that this is not the case in general. Indeed, as
pointed out in [10], typically the computation of the $\varepsilon$-CCS is
extremely difficult, since the evaluation of the probability
$\mathsf{Viol}(\theta)$ amounts to the computation of a multivariate integral,
which is NP-Hard [17].

Moreover, the set $\varepsilon$-CCS is often nonconvex, except for very
special cases. For example, [1, 18] show that the solution set of separable
chance constraints can be written as the union of cones, which is nonconvex in
general.

######  Example 2 (Example of nonconvex $\varepsilon$-CCS).

To illustrate these inherent difficulties, we consider the following three-
dimensional example ($n_{\theta}=3$) with $w=\left\\{w_{1},w_{2}\right\\}$,
where the first uncertainty $w_{1}\in\mathbb{R}^{3}$ is a three-dimensional
normal-distributed random vector with zero mean and covariance matrix

$\Sigma=\left[\begin{array}[]{ccc}4.5&2.26&1.4\\\ 2.26&3.58&1.94\\\
1.4&1.94&2.19\end{array}\right],$

and the second uncertainty $w_{2}\in\mathbb{R}^{3}$ is a three-dimensional
random vector whose elements are uniformly distributed in the interval
$[0,1]$. The set of viable design parameters is given by $n_{\ell}=4$
uncertain linear inequalities of the form

$F(w)\theta\leq\mathbf{1}_{4},\quad
F(w)=\left[\begin{array}[]{cccc}w_{1}&w_{2}&(2w_{1}-w_{2})&w_{1}^{2}\end{array}\right]^{\top}.$
(7)

The square power $w_{1}^{2}$ is to be interpreted element-wise.

In this case, to obtain a graphical representation of the set
$\mathbb{X}_{\varepsilon}$, we resorted to gridding the set $\Theta$ and, for
each point $\theta$ in the grid, to approximate the probability through a
Monte Carlo computation. This procedure is clearly unaffordable for higher
dimensions frameworks. In Figure 2 we report the plot of the computed
$\varepsilon$-CCS set for different values of $\varepsilon$. We observe that
the set is indeed nonconvex.

Figure 2: The $\varepsilon$-CCS set for $\varepsilon=0.15$ (smaller set),
$\varepsilon=0.30$ (intermediate set), and $\varepsilon=0.45$ (larger set). We
observe that all sets are nonconvex, but the nonconvexity is more evident for
larger values of $\varepsilon$, corresponding to larger levels of accepted
violation, while the set $\mathbb{X}_{\varepsilon}$ appears “almost convex”
for small values of $\varepsilon$. This kind of behaviour is in accordance
with a recent result that prove convexity of the $\varepsilon$-CCS for values
of $\varepsilon$ going to zero, and it is usually referred to as eventual
convexity [19].

###  2.1 Chance constrained optimization

Finding an optimal $\theta\in\mathbb{X}_{\varepsilon}$ for a given cost
function $J:~{}\mathbb{R}^{n_{\theta}}\rightarrow\mathbb{R}$, leads to the
chance constrained optimization (CCO) problem

$\min_{\theta\in\mathbb{X}_{\varepsilon}}J(\theta),$ (8)

where the cost-function $J(\theta)$ is usually assumed to be a convex, often
even a quadratic or linear function.

We remark that the solution of the CCO problem (8) is in general NP-hard, for
the same reasons reported before. We also note that several stochastic
optimization problems arising in different application contexts can be
formulated as a CCO. Typical examples are for instance the reservoir system
design problem proposed in [20], where the problem is to minimize the total
building and penalty costs while satisfying demands for all sites and all
periods with a given probability, or the cash matching problem [21], where one
aims at maximizing the portfolio value at the end of the planning horizon
while covering all scheduled payments with a prescribed probability. CCO
problems also frequently arise in short-term planning problems in power
systems. These optimal power flow (OPF) problems are routinely solved as part
of the real-time operation of the power grid. The aim is determining minimum-
cost production levels of controllable generators subject to reliably
delivering electricity to customers across a large geographical area, see e.g.
[8] and references therein.

In the next subsections, we report two control-related problems which served
as motivation of our study.

###  2.2 First motivating example: Stochastic MPC

To motivate the proposed approach, we consider the Stochastic MPC framework
proposed in [12, 11]. We are given a discrete-time system

$x_{k+1}=A(\sigma_{k})x_{k}+B(\sigma_{k})u_{k}+a_{\sigma}(\sigma_{k}),$ (9)

subject to generic uncertainty $\sigma_{k}\in\mathbb{R}^{n_{\sigma}}$, with
state $x_{k}\in\mathbb{R}^{n_{x}}$, control input
$u_{k}\in\mathbb{R}^{n_{u}}$, and the vector valued function
$a_{\sigma}(\sigma_{k})$ representing additive disturbance affecting the
system state. The system matrices $A(\sigma_{k})$ and $B(\sigma_{k})$, of
appropriate dimensions, are (possibly nonlinear) functions of the uncertainty
$\sigma_{k}$ at step $k$. For $k=1,2,\ldots$, the disturbances $\sigma_{k}$
are modeled as realizations of a stochastic process. In particular,
$\sigma_{k}$ are assumed to be independent and identically distributed (iid)
realizations of zero-mean random variables with support
$\mathcal{S}\subseteq\mathbb{R}^{n_{\sigma}}$. Note that the presence of both
additive and multiplicative uncertainty, combined with the nonlinear
dependence on the uncertainty, renders the problem particularly arduous.
Furthermore, we remark that the system representation in (9) is very general,
and encompasses, among others, those in [11, 12, 22].

Given the model (9) and a realization of the state $x_{k}$ at time $k$, state
predictions $t$ steps ahead are random variables as well and are denoted
$x_{t|k}$, to differentiate it from the realization $x_{t+k}$. Similarly
$u_{t|k}$ denotes predicted inputs that are computed based on the realization
of the state $x_{k}$.

Contrary to [11, 12, 22], where the system dynamics were subject to individual
state and input chance constraints, here we take a more challenging route, and
we consider joint state and input chance constraints of the form 111The case
where one wants to impose hard input constraints can be also be formulated in
a similar framework, see e.g. [11].

$\mathsf{Pr}_{\boldsymbol{\sigma}}\left\\{H_{x}x_{t|k}+H_{u}u_{t|k}\leq\mathbf{1}_{n_{t}}|x_{k}\right\\}\geq
1-\varepsilon,$ (10)

with $t\in\\{0,\ldots,T-1\\}$, $\varepsilon\in(0,1)$, and
$H_{x}\in\mathbb{R}^{n_{\ell}\times n_{x}}$,
$H_{u}\in\mathbb{R}^{n_{\ell}\times n_{u}}$.

The probability $\mathsf{Pr}_{\boldsymbol{\sigma}}$ is measured with respect
to the sequence ${\boldsymbol{\sigma}}=\\{\sigma_{t}\\}_{t>k}$. Hence,
equation (10) states that the probability of violating the linear constraint
$H_{x}x+H_{u}u\leq 1$ for any future realization of the disturbance should not
be larger than $\varepsilon$.

The objective is to derive an asymptotically stabilizing control law for the
system (9) such that, in closed loop, the constraint (10) is satisfied.
Following the approach in [12], a stochastic MPC algorithm is considered to
solve the constrained control problem. The approach is based on repeatedly
solving a stochastic optimal control problem over a finite, moving horizon,
but implementing only the first control action. The design parameter $\theta$
is then given by the control sequence
$\mathbf{u}_{k}=(u_{0|k},u_{1|k},...,u_{T-1|k})$ and the prototype optimal
control problem to be solved at each sampling time $k$ is defined by the cost
function

$\displaystyle J_{T}(x_{k},\mathbf{u}_{k})=$
$\displaystyle\mathbb{E}\left\\{\sum_{t=0}^{T-1}\left(x_{t|k}^{\top}Qx_{t|k}+u_{t|k}^{\top}Ru_{t|k}\right)+x_{T|k}^{\top}Px_{T|k}~{}|~{}x_{k}\right\\},$

with $Q\in\mathbb{R}^{n_{x}\times n_{x}}$, $Q\succeq 0$,
$R\in\mathbb{R}^{n_{u}\times n_{u}}$, $R\succ 0$, and appropriately chosen
$P\succ 0$, subject to the system dynamics (9) and constraints (10).

The online solution of the stochastic MPC problem remains a challenging task
but several special cases, which can be evaluated exactly, as well as methods
to approximate the general solution have been proposed in the literature. The
approach followed in this work was first proposed in [11, 12], where an
offline sampling scheme was introduced. Therein, with a prestabilizing input
parameterization

$u_{t|k}=Kx_{t|k}+v_{t|k},$ (12)

with suitably chosen control gain $K\in\mathbb{R}^{n_{u}\times n_{x}}$ and new
design parameters $v_{t|k}\in\mathbb{R}^{n_{u}}$, equation (9) is solved
explicitly for the predicted states $x_{1|k},\ldots,x_{T|k}$ and predicted
inputs $u_{0|k},\ldots,u_{T-1|k}$. In this case, the expected value of the
finite-horizon cost (2.2) can be evaluated offline, leading to a quadratic
cost function of the form

$J_{T}(x_{k},\mathbf{v}_{k})=\begin{bmatrix}x_{k}^{\top}&\textbf{v}_{k}^{\top}&\textbf{1}_{n_{x}}^{\top}\end{bmatrix}\tilde{S}\begin{bmatrix}x_{k}\\\
\textbf{v}_{k}\\\ \textbf{1}_{n_{x}}\\\ \end{bmatrix}$ (13)

in the deterministic variables
$\mathbf{v}_{k}=(v_{0|k},v_{1|k},...,v_{T-1|k})$ and $x_{k}$.

Focusing now on the constraint definition, we notice that by introducing the
uncertainty sequence
$\boldsymbol{\sigma}_{k}=\\{\sigma_{t}\\}_{t=k,...,k+T-1}$, we can rewrite the
joint chance constraint defined by equation (10) as

$\displaystyle\mathbb{X}_{\varepsilon}^{\textsc{smpc}}=\left\\{\
\begin{bmatrix}x_{k}\\\
\mathbf{v}_{k}\end{bmatrix}\in\mathbb{R}^{n_{x}+n_{u}T}~{}:~{}\right.$
$\displaystyle\Bigl{.}\mathsf{Pr}_{\boldsymbol{\sigma}_{k}}\left\\{\begin{bmatrix}f_{\ell}^{x}(\boldsymbol{\sigma}_{k})\\\
f_{\ell}^{v}(\boldsymbol{\sigma}_{k})\end{bmatrix}^{\top}\begin{bmatrix}x_{k}\\\
\mathbf{v}_{k}\end{bmatrix}\leq 1,\ell\in[n_{\ell}]\right\\}\geq
1-\varepsilon\Bigr{\\}},$ (14)

with
$f_{\ell}^{x}:\mathbb{R}^{n_{\sigma}}\to\mathbb{R}^{n_{x}},f_{\ell}^{v}:\mathbb{R}^{n_{\sigma}}\to\mathbb{R}^{n_{u}T}$
being known functions of the sequence of random variables
$\boldsymbol{\sigma}_{k}$. We remark that, in the context of this paper,
neither the detailed derivation of the cost matrix $\tilde{S}$ in (13) nor
that of $f_{\ell}^{v},f_{\ell}^{x}$ are relevant for the reader, who can refer
to [12, Appendix A] for details. Note that, by defining
$\theta=[x_{k}^{\top},\mathbf{v}_{k}^{\top}]^{\top}$, (14) is given in the
form of (5) .

As discussed in [11], obtaining a good and simple enough approximation of the
set $\mathbb{X}_{\varepsilon}^{\textsc{smpc}}$ is extremely important for
online implementation of SMPC schemes. In particular, if we are able to
replace the set $\mathbb{X}_{\varepsilon}^{\textsc{smpc}}$ by a suitable inner
approximation, we would be able to guarantee probabilistic constraint
satisfaction of the ensuing SMPC scheme. On the other hand, we would like this
inner approximation to be simple enough, so to render the online computations
fast enough.

###  2.3 Second motivating example: probabilistic set membership estimation

Suppose that there exists $\bar{\theta}\in\Theta$ such that

$|y-\bar{\theta}^{T}\varphi(x)|\leq\rho,\;\forall(x,y)\in\mathbb{W}\subseteq\mathbb{R}^{n_{x}}\times\mathbb{R},$

where $\varphi:\mathbb{R}^{n_{x}}\to\mathbb{R}^{n_{\theta}}$ is a (possibly
non-linear) regressor function, and $\rho>0$ accounts for modelling errors.
The (deterministic) set membership estimation problem, see [23], [24],
consists of computing the set of parameters $\theta$ that satisfy the
constraint

$|y-\theta^{T}\varphi(x)|\leq\rho$

for all possible values of $(x,y)\in\mathbb{W}$. In the literature, this set
is usually referred to as the feasible parameter set, that is

${\mathsf{FPS}}\doteq\\{\;\theta\in\Theta\;:\;|y-\theta^{T}\varphi(x)|\leq\rho,\;\forall(x,y)\in\mathbb{W}\;\\}.$
(15)

If, for given $w=(x,y)$, we define the set

$\mathbb{X}(w)=\\{\;\theta\in\Theta\;:\;|y-\theta^{T}\varphi(x)|\leq\rho\;\\},$

then the feasible parameter set ${\mathsf{FPS}}$ can be rewritten as

${\mathsf{FPS}}=\\{\;\theta\in\Theta\;:\;\theta\in\mathbb{X}(w),\;\forall
w\in\mathbb{W}\;\\}.$

The deterministic set membership problem suffers from the following
limitations in real applications: i) due to the possible non-linearity of
$\varphi(\cdot)$, checking if a given $\theta\in\Theta$ satisfies the
constraint $\theta\in\mathbb{X}(w)$, for every $w\in\mathbb{W}$, is often a
difficult problem; ii) in many situations, only samples of $\mathbb{W}$ are
available: thus, the robust constraint cannot be checked and only outer bounds
of ${\mathsf{FPS}}$ can be computed; and iii) because of outliers and possible
non finite support of $\mathbb{W}$, set ${\mathsf{FPS}}$ is often empty
(especially for small values of $\rho$).

If a probability distribution is defined on $\mathbb{W}$, the probabilistic
set membership estimation problem is that of characterizing the set of
parameters $\theta$ that satisfy

$\mathsf{Pr}_{\mathbb{W}}\\{|y-\theta^{T}\varphi(x)|\leq\rho\\}\geq
1-\epsilon,$

for a given probability parameter $\epsilon\in(0,1)$. Hence, we can define
${\mathsf{FPS}}_{\epsilon}$ the set of parameters that satisfy the previous
probabilistic constraint, that is,

${\mathsf{FPS}}_{\epsilon}=\\{\;\theta\in\Theta\;:\;\mathsf{Pr}_{\mathbb{W}}\\{\theta\in\mathbb{X}(w)\\}\geq
1-\epsilon\;\\}.$

It is immediate to notice that this problem fits in the formulation proposed
in this section: It suffices to define

$F(w)=\left[\begin{array}[]{c}\varphi^{T}(x)\\\
-\varphi^{T}(x)\end{array}\right],\;g(w)=\left[\begin{array}[]{c}\rho+y\\\
\rho-y\end{array}\right].$

###  2.4 Chance constrained approximations

Motivated by the discussion above, we are ready to formulate the main problem
studied in this paper.

######  Problem 1 ($\varepsilon$-CCS approximation).

Given the set of linear inequalities (2), and a violation parameter
$\varepsilon$, find an inner approximation of the set
$\mathbb{X}_{\varepsilon}$. The approximation should be: i) simple enough, ii)
easily computable.

A solution to this problem is provided in the paper. In particular, regarding
i), we present a solution in which the approximating set is represented by few
linear inequalities. Regarding ii), we propose a computationally efficient
procedure for its construction (see Algorithm 1).

Before presenting our approach, in the next section we provide a brief
literature overview of different methods presented in the literature to
construct approximations of the $\varepsilon$-CCS set.

##  3 Overview on different approaches to $\varepsilon$-CCS approximations

The construction of computational efficient approximations to
$\varepsilon$-CCS is a long-standing problem. In particular, the reader is
referred to the recent work [10], which provides a rather complete discussion
on the topic, and covers the most recent results. The authors distinguish
three different approaches, which we very briefly revisit here.

###  3.1 Exact techniques

In some very special cases, the $\varepsilon$-CCS is convex and hence the CCO
problem admits a unique solution. This is the case, for instance, of
individual chance constraints with $w$ being Gaussian [25]. Other important
examples of convexity of the set $\mathbb{X}_{\varepsilon}$ involve log-
concave distribution [1, 26]. General sufficient conditions on the convexity
of chance constraints may be found in [27, 28, 29, 19]. However, all these
cases are very specific and hardly extend to joint chance constraints
considered on this work.

###  3.2 Robust techniques

A second class of approaches consist in finding deterministic conditions that
allow to construct a set $\underline{\mathbb{X}}$, which is a guaranteed inner
convex approximation of the probabilistic set $\mathbb{X}_{\varepsilon}$. The
classical solution consists in the applications of Chebyshev-like
inequalities, see e.g. [30, 31]. More recent techniques, which are proved
particularly promising, involve robust optimization [3], as the convex
approximations introduced in [32]. A particular interesting convex relaxation
involves the so-called Conditional Value at Risk (CVaR), see [33] and
references therein. Finally, we point out some recent techniques based on
polynomial moments relaxations [34, 35]. Nonetheless, it should be remarked
that these techniques usually suffer from conservatism and computational
complexity issues, especially in the case of joint chance constraints.

###  3.3 Sample-based techniques

In recent years, a novel approach to approximate chance constraints, based on
random sampling of the uncertain parameters, has gained popularity, see e.g.
[4, 5] and references therein. Sampling-based techniques are characterized by
the use of a finite number $N$ of iid samples of the uncertainty
$\left\\{w^{(1)},w^{(2)},\ldots,w^{(N)}\right\\}$ drawn according to a
probability distribution $\mathsf{Pr}_{\mathbb{W}}$. To each sample
$w^{(i)},i\in[N]$, we can associate the following sampled set

$\mathbb{X}(w^{(i)})=\\{\;\theta\in\Theta\;:\;F(w^{(i)})\theta\leq
g(w^{(i)})\;\\},$ (16)

sometimes referred to as scenario, since it represents an observed instance of
our probabilistic constraint.

Then, the scenario approach considers the CCO problem (8) and approximates its
solution through the following scenario problem

$\displaystyle\theta^{*}_{sc}=\arg\min J(\theta)$ (17)
$\displaystyle\text{subject to }\theta\in\mathbb{X}(w^{(i)}),i\in[N].$

We note that, if the function $J(\theta)$ is convex, problem (17) becomes a
linearly constrained convex program, for which very efficient solution
approaches exist. A fundamental result [36, 37, 38, 39] provides a
probabilistic certification of the constraint satisfaction for the solution to
the scenario problem. In particular, it is shown that, under some mild
assumptions (non-degenerate problem), we have

$\mathsf{Pr}_{\mathbb{W}^{N}}\left\\{\mathsf{Viol}(\theta^{*}_{sc})>\varepsilon\right\\}\leq\mathbf{B}(n_{\theta}-1;N,\varepsilon),$
(18)

where the probability in (18) is measured with respect to the samples
$\\{w^{(1)},w^{(2)},\ldots,w^{(N)}$}. Moreover, the bound in (18) is shown to
be tight. Indeed, for the class of so-called fully-supported problems, the
bound holds with equality, i.e. the Binomial distribution
$\mathbf{B}(n_{\theta}-1;N,\varepsilon)$ represents the exact probability
distribution of the violation probability [37].

A few observations are at hand regarding the scenario approach and its
relationship with Problem 1. First, if we define the sampled constraints set
as

$\mathbb{X}_{N}\doteq\bigcap_{i=1}^{N}\mathbb{X}(w^{(i)}),$ (19)

we see that the scenario approach consists in approximating the constraint
$\theta\in\mathbb{X}_{\varepsilon}$ in (8) with its sampled version
$\theta\in\mathbb{X}_{N}$. On the other hand, it should be remarked that the
scenario approach cannot be used to derive any guarantee on the relationship
existing between $\mathbb{X}_{N}$ and $\mathbb{X}_{\varepsilon}$. Indeed, the
nice probabilistic property in (18) holds only for the optimum of the scenario
program $\theta^{*}_{sc}$. This is a fundamental point, since the scenario
results build on the so-called support constraints, which are defined for the
optimum point $\theta^{*}_{sc}$ only.

On the contrary, in our case we are interested in establishing a direct
relation (in probabilistic terms) between the set $\mathbb{X}_{N}$ and the
$\varepsilon$-CCS $\mathbb{X}_{\varepsilon}$. This is indeed possible, but
needs to resort to results based on Statistical Learning Theory [40],
summarized in the following lemma.

######  Lemma 1 (Learning Theory bound).

Given probabilistic levels $\delta\in(0,1)$ and $\varepsilon\in(0,0.14)$, if
the number of samples $N$ is chosen so that $N\geq N_{LT}$, with

$N_{LT}\doteq\frac{4.1}{\varepsilon}\Big{(}\ln\frac{21.64}{\delta}+4.39n_{\theta}\,\log_{2}\Big{(}\frac{8en_{\ell}}{\varepsilon}\Big{)}\Big{)},$
(20)

then
$\mathsf{Pr}_{\mathbb{W}^{N}}\left\\{\mathbb{X}_{N}\subseteq\mathbb{X}_{\varepsilon}\right\\}\geq
1-\delta$.

The lemma, whose proof is reported in Appendix A.1, is a direct consequence of
the results on VC-dimension of the so-called $(\alpha,k)$-Boolean Function,
given in [41].

######  Remark 2 (Sample-based SMPC).

The learning theory-based approach discussed in this section has been applied
in [11] to derive an _offline_ probabilistic inner approximation of the chance
constrained set $\mathbb{X}_{\varepsilon}^{\textsc{smpc}}$ defined in (14),
considering individual chance constraints. In particular, the bound (2) is a
direct extension to the case of joint chance constraints of the result proved
in [11]. Note that since we are considering multiple constraints at the same
time (like in (2)), the number of constraints $n_{\ell}$ enters into the
sample size bound. To explain how the SMPC design in [11] extends to the joint
chance constraints framework, we briefly recall it.

First, we extract offline (i.e. when designing the SMPC control) $N$ iid
samples of the uncertainty, $\boldsymbol{\sigma}_{k}^{(i)}$ of
$\boldsymbol{\sigma}_{k}$, and we consider the sampled set

$\displaystyle\mathbb{X}^{\textsc{smpc}}(\boldsymbol{\sigma}_{k}^{(i)})=\Biggl{\\{}\
\begin{bmatrix}x_{k}\\\
\mathbf{v}_{k}\end{bmatrix}:\begin{bmatrix}f_{\ell}^{x}(\boldsymbol{\sigma}_{k}^{(i)})\\\
f_{\ell}^{v}(\boldsymbol{\sigma}_{k}^{(i)})\end{bmatrix}^{\top}\begin{bmatrix}x_{k}\\\
\mathbf{v}_{k}\end{bmatrix}\leq 1,\Biggl{.}\ell\in[n_{\ell}]\Biggr{\\}},$

and
$\mathbb{X}_{N}^{\textsc{smpc}}\doteq\bigcap_{i=1}^{N}\mathbb{X}^{\textsc{smpc}}(\boldsymbol{\sigma}_{k}^{(i)})$.
Then, applying Lemma 1 with $n_{\theta}=n_{x}+n_{u}T$, we conclude that if we
extract $N\geq N_{LT}^{\textsc{smpc}}$ samples, it is guaranteed that, with
probability at least $1-\delta$, the sample approximation
$\mathbb{X}_{N}^{\textsc{smpc}}$ is a subset of the original chance constraint
$\mathbb{X}_{\varepsilon}^{\textsc{smpc}}$. Exploiting these results, the SMPC
problem can be approximated conservatively by the linearly constrained
quadratic program

$\displaystyle\min_{\mathbf{v}_{k}}~{}J_{T}(x_{k},\mathbf{v}_{k})\textrm{
subject to }(x_{k},\mathbf{v}_{k})\in\mathbb{X}_{N}^{\textsc{smpc}}.$ (21)

Hence the result reduces the original stochastic optimization program to an
efficiently solvable quadratic program. This represents an undiscussed
advantage, which has been demonstrated for instance in [12]. On the other
hand, it turns out that the ensuing number of linear constraints, equal to
$n_{\ell}\cdot N_{LT}^{\textsc{smpc}}$ may still be too large. For instance,
even for a moderately sized MPC problem with $n_{x}=5$ states, $n_{u}=2$
inputs, prediction horizon of $T=10$, simple interval constraints on states
and inputs (i.e. $n_{\ell}=2n_{x}+2n_{u}=14$), and for a reasonable choice of
probabilistic parameters, i.e. $\varepsilon=0.05$ and $\delta=10^{-6}$, we get
$N_{LT}^{\textsc{smpc}}=114,530$, which in turn corresponds to more than $1.6$
million linear inequalities. For this reason, in [11] a post-processing step
was proposed to remove redundant constraints. While it is indeed true that all
the cumbersome computations may be performed offline, it is still the case
that, in applications with stringent requirements on the solution time, the
final number of inequalities may easily become unbearable.

Remark 2 motivates the approach presented in the next section, which builds
upon the results presented in [13]. We show how the probabilistic scaling
approach directly leads to approximations of user-chosen complexity, which can
be directly used in applications instead of creating the need for a post-
processing step to reduce the complexity of the sampled set.

##  4 The Probabilistic Scaling Approach

We propose a novel sample-based approach, alternative to the randomized
procedures proposed so far, which allows to maintain the nice probabilistic
features of these techniques, while at the same time providing the designer
with a way of tuning the complexity of the approximation.

The main idea behind this approach consists of first obtaining a simple
initial approximation of the shape of the probabilistic set
$\mathbb{X}_{\varepsilon}$ by exploiting scalable simple approximating sets
(Scalable SAS) of the form

${\mathbb{S}}(\gamma)=\theta_{c}\oplus\gamma{\mathbb{S}}.$ (22)

These sets are described by a center point $\theta_{c}$ and a low-complexity
shape set ${\mathbb{S}}$. The center $\theta_{c}$ and the shape ${\mathbb{S}}$
constitute the design parameters of the proposed approach. By appropriately
selecting the shape ${\mathbb{S}}$, the designer can control the complexity of
the approximating set.

Note that we do not ask this initial set to have any guarantee of
probabilistic nature. What we ask is that this set is being able to “capture”
somehow the shape of the set $\mathbb{X}_{\varepsilon}$. Recipes on a possible
procedure for constructing this initial set are provided in section 5. The set
${\mathbb{S}}$ constitutes the starting point of a scaling procedure, which
allows to derive a probabilistic guaranteed approximation of the
$\varepsilon$-CCS, as detailed in the next section. In particular, we show how
an optimal scaling factor $\gamma$ can be derived so that the set (22) is
guaranteed to be an inner approximation of $\mathbb{X}_{\varepsilon}$ with the
desired confidence level $\delta$. We refer to the set ${\mathbb{S}}(\gamma)$
as Scalable SAS.

###  4.1 Probabilistic Scaling

In this section, we address the problem of how to scale the set
${\mathbb{S}}(\gamma)$ around its center $\theta_{c}$ to guarantee, with
confidence level $\delta\in(0,1)$, the inclusion of the scaled set into
$\mathbb{X}_{\varepsilon}$. Within this sample-based procedure we assume that
$N_{\gamma}$ iid samples $\\{w^{(1)},\ldots,w^{(N_{\gamma})}\\}$ are obtained
from $\mathsf{Pr}_{\mathbb{W}}$ and based on these, we show how to obtain a
scalar $\bar{\gamma}>0$ such that

$\mathsf{Pr}_{\mathbb{W}^{N_{\gamma}}}\\{{\mathbb{S}}(\bar{\gamma})\subseteq\mathbb{X}_{\varepsilon}\\}\geq
1-\delta.$

To this end, we first define the scaling factor associated to a given
realisation of the uncertainty.

######  Definition 3 (Scaling factor).

Given a Scalable SAS ${\mathbb{S}}(\gamma)$, with given center $\theta_{c}$
and shape ${\mathbb{S}}\subset\Theta$, and a realization $w\in\mathbb{W}$, we
define the scaling factor of ${\mathbb{S}}(\gamma)$ relative to $w$ as

$\gamma(w)\doteq\left\\{\begin{array}[]{cc}0&\,\,\,\mbox{if}\;\theta_{c}\not\in\mathbb{X}(w)\\\
\max\limits_{{\mathbb{S}}(\gamma)\subseteq\mathbb{X}(w)}\gamma&\,\,\,\mbox{otherwise}.\end{array}\right.$

with $\mathbb{X}(w)$ defined as in (16).

That is $\gamma(w)$ represents the maximal scaling that can be applied to
${\mathbb{S}}(\gamma)=\theta_{c}\oplus\gamma{\mathbb{S}}$ around the center
$\theta_{c}$ so that ${\mathbb{S}}(\gamma)\subseteq\mathbb{X}(w)$. The
following theorem states how to obtain, by means of sampling, a scaling factor
$\bar{\gamma}$ that guarantees, with high probability, that
${\mathbb{S}}(\bar{\gamma})\subseteq\mathbb{X}_{\varepsilon}$.

######  Theorem 1 (Probabilistic scaling).

Given a candidate Scalable SAS ${\mathbb{S}}(\gamma)$, with
$\theta_{c}\in\mathbb{X}_{\varepsilon}$, accuracy parameter
$\varepsilon\in(0,1)$, confidence level $\delta\in(0,1)$, and a discarding
integer parameter $r\geq 0$, let $N_{\gamma}$ be chosen such that

$\mathbf{B}(r;N_{\gamma},\varepsilon)\leq\delta.$ (23)

Draw $N_{\gamma}$ iid samples $\\{w^{(1)},w^{(2)},\ldots,w^{(N_{\gamma})}\\}$
from distribution $\mathsf{Pr}_{\mathbb{W}}$, compute the corresponding
scaling factor

$\gamma_{i}\doteq\gamma(w^{(i)}),$ (24)

for $i\in[N_{\gamma}]$ according to Definition 3, and let
$\bar{\gamma}=\gamma_{1+r:N_{\gamma}}$. Then, with probability no smaller than
$1-\delta$,

${\mathbb{S}}(\bar{\gamma})=\theta_{c}\oplus\bar{\gamma}{\mathbb{S}}\subseteq\mathbb{X}_{\varepsilon}.$

Proof: If $\bar{\gamma}=0$, then we have
${\mathbb{S}}(\bar{\gamma})\equiv\theta_{c}\in\mathbb{X}_{\varepsilon}$.
Hence, consider $\bar{\gamma}>0$. From Property 1 in Appendix A.2, we have
that $\bar{\gamma}0$ satisfies, with probability no smaller than $1-\delta$,
that
$\mathsf{Pr}_{\mathbb{W}}\\{{\mathbb{S}}(\gamma)\not\subseteq\mathbb{X}(w)\\}\leq\varepsilon$.
Equivalently,
$\mathsf{Pr}_{\mathbb{W}}\\{{\mathbb{S}}(\gamma)\subseteq\mathbb{X}(w)\\}>1-\varepsilon.$
This can be rewritten as $\mathsf{Pr}_{\mathbb{W}}\\{F(w)^{\top}\theta\leq
g(w),\;\;\forall\theta\in{\mathbb{S}}(\gamma)\\}>1-\varepsilon,$ and it
implies that the probability of violation in
$\theta_{c}\oplus\bar{\gamma}{\mathbb{S}}$ is no larger than $\varepsilon$,
with probability no smaller than $1-\delta$. ∎

In the light of the theorem above, from now on we will assume that the
Scalable SAS is such that $\theta_{c}\in\mathbb{X}_{\varepsilon}$. The above
result leads to the following simple algorithm, in which we summarise the main
steps for constructing the scaled set, and we provide an explicit way of
determining the discarding parameter $r$.

Algorithm 1 Probabilistic SAS Scaling

1:Given a candidate Scalable SAS ${\mathbb{S}}(\gamma)$, and probability
levels $\varepsilon$ and $\delta$, choose
$N_{\gamma}\geq\frac{7.47}{\varepsilon}\ln\frac{1}{\delta}\quad\text{ and
}\quad r=\left\lfloor\frac{\varepsilon N_{\gamma}}{2}\right\rfloor.$ (25)

2:Draw $N_{\gamma}$ samples of the uncertainty
$w^{(1)},\ldots,w^{(N_{\gamma})}$

3:for $i=1$ to $N_{\gamma}$ do

4: Solve the optimization problem  $\displaystyle\gamma_{i}\doteq$
$\displaystyle\max_{{\mathbb{S}}(\gamma)\subseteq\mathbb{X}(w^{(i)})}\gamma$
(26)

5:end for

6:Return $\bar{\gamma}=\gamma_{1+r:N_{\gamma}}$, the $(1+r)$-th smallest value
of $\gamma_{i}$.

A few comments are in order regarding the algorithm above. In step 4, for each
uncertainty sample $w^{(i)}$ one has to solve an optimization problem, which
amounts to finding the largest value of $\gamma$ such that
${\mathbb{S}}(\gamma)$ is contained in the set $\mathbb{X}(w^{(i)})$ defined
in (16). If the SAS is chosen accurately, we can show that this problem is
convex and computationally very efficient: this is discussed in Section 5.
Then, in step 6, one has to re-order the set
$\\{\gamma_{1},\gamma_{2},\ldots,\gamma_{N_{\gamma}}\\}$ so that the first
element is the smallest one, the second element is the second smallest one,
and so on and so fort, and then return the $r+1$-th element of the reordered
sequence. The following Corollary applies to Algorithm 1.

######  Corollary 1.

Given a candidate SAS set in the form
${\mathbb{S}}(\gamma)=\theta_{c}\oplus\gamma{\mathbb{S}}$, assume that
$\theta_{c}\in\mathbb{X}_{\varepsilon}$. Then, Algorithm 1 guarantees that
${\mathbb{S}}(\bar{\gamma})\subseteq\mathbb{X}_{\varepsilon}$ with probability
at least $1-\delta$.

Proof: The result is a direct consequence of Theorem 1, which guarantees that,
for given $r\geq 0$,
$\mathsf{Pr}\\{{\mathbb{S}}(\gamma)\subseteq\mathbb{X}_{\varepsilon}\\}$ is
guaranteed if the scaling is performed on a number of samples satisfying (23).
From [42, Corollary 1]) it follows that, in order to satisfy (23) it suffices
to take $N_{\gamma}$ such that

$N_{\gamma}\geq\frac{1}{\varepsilon}\left(r+\ln\frac{1}{\delta}+\sqrt{2r\ln\frac{1}{\delta}}\right).$
(27)

Since $r=\lfloor\frac{\varepsilon N}{2}\rfloor$, we have that
$r\leq\frac{\varepsilon N}{2}$. Thus, inequality (27) is satisfied if

$\displaystyle N_{\gamma}$ $\displaystyle\geq$
$\displaystyle\frac{1}{\varepsilon}\left(\frac{\varepsilon
N_{\gamma}}{2}+\ln\frac{1}{\delta}+\sqrt{\varepsilon
N_{\gamma}\ln\frac{1}{\delta}}\right)$ $\displaystyle=$
$\displaystyle\frac{N_{\gamma}}{2}+\frac{1}{\varepsilon}\ln\frac{1}{\delta}+\sqrt{N_{\gamma}\frac{1}{\varepsilon}\ln\frac{1}{\delta}}.$

Letting $\nabla\doteq\sqrt{N_{\gamma}}$ and
$\alpha\doteq\sqrt{\frac{1}{\varepsilon}\ln\frac{1}{\delta}}$222Note that both
quantities under square root are positive., the above inequality rewrites
$\nabla^{2}-2\alpha\nabla-2\alpha^{2}\geq 0,$ which has unique positive
solution $\nabla\geq(1+\sqrt{3})\alpha$. In turn, this rewrites as

$N_{\gamma}\geq\frac{(1+\sqrt{3})^{2}}{\varepsilon}\ln\frac{1}{\delta}.$

The formula (25) follows by observing that $(1+\sqrt{3})^{2}<~{}7.47$. ∎

In the next sections, we provide a “library” of possible candidates SAS
shapes. We remind that these sets need to comply to two main requirements: i)
being a simple and low-complexity representation; and ii) being able to
capture the original shape of the $\varepsilon$-CCS. Moreover, in the light of
the discussion after Algorithm 1, we also ask these sets to be convex.

##  5 Candidate SAS: Sampled-polytope

First, we note that the most straightforward way to design a candidate SAS is
again to recur to a sample-based procedure: we draw a fixed number $N_{S}$ of
“design” uncertainty samples333These samples are denoted with a tilde to
distinguish them from the samples used in the probabilistic scaling procedure.
$\\{\tilde{w}^{(1)},\ldots,\tilde{w}^{(N_{S})}\\}$, and construct an initial
sampled approximation by introducing the following sampled-polytope SAS

${\mathbb{S}}_{N_{S}}=\bigcap_{j=1}^{N_{S}}\mathbb{X}(\tilde{w}^{(j)}).$ (28)

Note that the sampled polytope ${\mathbb{S}}_{N_{S}}$, by construction, is
given by the intersection of $n_{\ell}N_{S}$ half-spaces. Hence, we observe
that this approach provides very precise control on the final complexity of
the approximation, through the choice of the number of samples $N_{S}$.
However, it is also clear that a choice for which $N_{S}<<N_{LT}$ implies that
the probabilistic properties of ${\mathbb{S}}_{N_{S}}$ before scaling will be
very bad. However, we emphasize again that this initial geometry doesn’t have
nor require any probabilistic guarantees, which are instead provided by the
probabilistic scaling discussed in Section 4.1. It should be also remarked
that this is only one possible heuristic. For instance, along this line one
could as well draw many samples and then apply a clustering algorithm to boil
it down to a desired number of samples.

We remark that, in order to apply the scaling procedure, we need to define a
center around which to apply the scaling procedure. To this end, we could
compute the so-called Chebyshev center, defined as the center of largest ball
inscribed in ${\mathbb{S}}_{N_{S}}$, i.e.
$\theta_{c}=\mathsf{Cheb}({\mathbb{S}}_{N_{S}})$. We note that computing the
Chebyshev center of a given polytope is an easy convex optimization problem,
for which efficient algorithms exist, see e.g. [43]. A possible alternative
would be the analytic center of ${\mathbb{S}}_{N_{S}}$, whose computation is
even easier (see [43] for further details). Once the center $\theta_{c}$ has
been determined, the scaling procedure can be applied to the set
${\mathbb{S}}_{N_{S}}(\gamma)\doteq\theta_{c}\oplus\gamma\\{{\mathbb{S}}_{N_{S}}\ominus\theta_{c}\\}$.
Note that the center needs to be inside $\mathbb{X}_{\varepsilon}$. Aside for
that, the choice of $\theta_{c}$ only affects the goodness of the shape, but
we can never know a priori if the analytic center is a better choice than any
random center in $\mathbb{X}_{\varepsilon}$.

(a) ${\mathbb{S}}_{N_{S}}$ with $N_{S}=100$. $\rightarrow$ $\gamma=0.8954$

(b) ${\mathbb{S}}_{N_{S}}$ with $N_{S}=1,000$. $\rightarrow$ $\gamma=1.2389$

(c) LT-based (Lemma 1). $N_{LT}=52,044$

Figure 3: (a-b) Probabilistic scaling approximations of the $\varepsilon$-CCS.
Scaling procedure applied to a sampled-polytope with $N_{S}=100$ (a) and
$N_{S}=1,000$ (b). The initial sets are depicted in red, the scaled ones in
green. (c) Approximation obtained by direct application of Lemma 1. Note that,
in this latter case, to plot the set without out-of-memory errors a pruning
procedure [44] of the $52,044$ linear inequalities was necessary.

######  Example 3 (Sample-based approximations).

To illustrate how the proposed scaling procedure works in practice in the case
of sampled-polytope SAS, we revisit Example 2. To this end, a pre-fixed number
$N_{S}$ of uncertainty samples were drawn, and the set inequalities

$F(\tilde{w}^{(j)})\theta\leq g(\tilde{w}^{(j)}),\quad j\in[N_{S}],$

with $F(w),g(w)$ defined in (7), were constructed, leading to the candidate
set ${\mathbb{S}}_{N_{S}}$. Then, the corresponding Chebyshev center was
computed, and Algorithm 1 was applied with $\varepsilon=0.05$,
$\delta=10^{-6}$, leading to $N_{\gamma}=2,120$.

We note that, in this case, the solution of the optimization problem in (26)
may be obtained by bisection on $\gamma$. Indeed, for given $\gamma$, checking
if ${\mathbb{S}}_{N_{S}}(\gamma)\subseteq\mathbb{X}(w^{(i)})$ amounts to
solving some simple linear programs.

Two different situations were considered: a case where the number of
inequalities is rather small $N_{S}=100$, and a case where the complexity of
the SAS is higher, i.e. $N_{S}=1,000$. The outcome procedure is illustrated in
Figure 3. We can observe that, for a small $N_{S}$ – Fig. 3(a) – the initial
approximation is rather large (although it is contained in
$\mathbb{X}_{\varepsilon}$, we remark that we do not have any guarantee that
this will happen). In this case, the probabilistic scaling returns
$\gamma=0.8954$ which is less than one. This means that, in order to obtain a
set fulfilling the desired probabilistic guarantees, we need to shrink it
around its center. In the second case, for a larger number of sampled
inequalities – Fig. 3(b) \- the initial set (the red one) is much smaller, and
the scaling procedure inflates the set by returning a value of $\gamma$
greater than one, i.e. $\gamma=1.2389$. Note that choosing a larger number of
samples for the computation of the initial set does not imply that the final
set will be a better approximation of the $\varepsilon$-CCS.

Finally, we compare this approach to the scenario-like ones discussed in
Subsection 3.3. To this end, we also draw the approximation obtained by
directly applying the Learning Theory bound (20). Note that in this case,
since $n_{\theta}=3$ and $n_{\ell}=4$, we need to take $N_{LT}=13,011$
samples, corresponding to $52,044$ linear inequalities. The resulting set is
represented in Fig. 3(c). We point out that using this approximation i) the
set is much more complex, since the number of involved inequalities is much
larger, ii) the set is much smaller, hence providing a much more conservative
approximation of the $\varepsilon$-CCS. Hence, the ensuing chance-constrained
optimization problem will be computationally harder, and lead to a solution
with a larger cost or even to an infeasible problem, in cases where the
approximating set is too small.

##  6 Candidate SAS: Norm-based SAS

In this section, we propose a procedure in which the shape of the scalable SAS
may be selected a-priori. This corresponds to situations where the designer
wants to have full control in the final shape in terms of structure and
complexity. The main idea is to define so-called norm-based SAS of the form

${\mathbb{S}_{p}}(\gamma)\doteq\theta_{c}\oplus\gamma H\mathbb{B}_{p}^{s}$
(29)

where $\mathbb{B}_{p}^{s}$ is a $\ell_{p}$-ball in $\mathbb{R}^{s}$,
$H\in\mathbb{R}^{n_{\theta},s}$, with $s\geq n_{\theta}$, is a design matrix
(not necessarily square), and $\gamma$ is the scaling parameter. Note that
when the matrix $H$ is square (i.e. $s=n_{\theta}$) and positive definite
these sets belong to the class of $\ell_{p}$-norm based sets originally
introduced in [45]. In particular, in case of $\ell_{2}$ norm, the sets are
ellipsoids. This particular choice is the one studied in [14]. Here, we extend
this approach to a much more general family of sets, which encompasses for
instance zonotopes, obtained by letting $p=\infty$ and $s\geq n_{\theta}$.
Zonotopes have been widely studied in geometry, and have found several
applications in systems and control, in particular for problems of state
estimation and robust Model Predictive Control, see e.g. [46].

###  6.1  Scaling factor computation for norm-bases SAS

We recall that the scaling factor $\gamma(w)$ is defined as $0$ if
$\theta_{c}\not\in\mathbb{X}(w)$ and as the largest value $\gamma$ for which
${\mathbb{S}_{p}}(\gamma)\subseteq\mathbb{X}(w)$ otherwise. The following
theorem, whose proof is reported in Appendix A.3, provides a direct and simple
way to compute in closed form the scaling factor for a given candidate norm-
based SAS.

######  Theorem 2 (Scaling factor for norm-based SAS).

Given a norm-based SAS ${\mathbb{S}}(\gamma)$ as in (29), and a realization
$w\in\mathbb{W}$, the scaling factor $\gamma(w)$ can be computed as

$\gamma(w)=\min_{\ell\in[n_{\ell}]}\;\gamma_{\ell}(w),$

with $\gamma_{\ell}(w)$, $\ell\in[n_{\ell}]$, given by

$\gamma_{\ell}(w)=\left\\{\begin{array}[]{ccl}0&\mbox{if
}&\tau_{\ell}(w)<0,\\\ \infty&\mbox{if}&\tau_{\ell}(w)\geq 0\mbox{ and
}\rho_{\ell}(w)=0,\\\
{\displaystyle{\frac{\tau_{\ell}(w)}{\rho_{\ell}(w)}}}&\mbox{if}&\tau_{\ell}(w)\geq
0\mbox{ and }\rho_{\ell}(w)>0,\end{array}\right.$ (30)

where $\tau_{\ell}(w)\doteq g_{\ell}(w)-f_{\ell}^{T}(w)\theta_{c}$ and
$\rho_{\ell}(w)\doteq\|H^{T}f_{\ell}(w)\|_{p^{*}}$, with $\|\cdot\|_{p}^{*}$
being the dual norm of $\|\cdot\|_{p}$.

Note that $\gamma(w)$ is equal to zero if and only if $\theta_{c}$ is not
included in the interior of $\mathbb{X}(w)$.

###  6.2 Construction of a candidate norm-based set

Similarly to Section 5, we first draw a fixed number $N_{S}$ of “design”
uncertainty samples $\\{\tilde{w}^{(1)},\ldots,\tilde{w}^{(N_{S})}\\},$ and
construct an initial sampled approximation by introducing the following
sampled-polytope SAS ${\mathbb{S}}_{N_{S}}$ as defined in
$\eqref{eq:sampledSAS}$. Again, we consider the Chebyshev center of
${\mathbb{S}}_{N_{S}}$, or its analytical center as a possible center
$\theta_{c}$ for our approach.

Given ${\mathbb{S}}_{N_{S}}$, $s\geq n_{\theta}$ and $p\in\\{1,2,\infty\\}$,
the objective is to compute the largest set $\theta_{c}\oplus
H\mathbb{B}^{s}_{p}$ included in ${\mathbb{S}}_{N_{S}}$. To this end, we
assume that we have a function $\mathsf{Vol}_{p}(H)$ that provides a measure
of the size of $H\mathbb{B}^{s}_{p}$. That is, larger values of
$\mathsf{Vol}_{p}(H)$ are obtained for increasing sizes of
$H\mathbb{B}^{s}_{p}$.

######  Remark 3 (On the volume function).

The function $\mathsf{Vol}_{p}(H)$ may be seen as a generalization of the
classical concept of Lebesgue volume of the set ${\mathbb{S}}_{N_{S}}$.
Indeed, when $H$ is a square positive definite matrix, some possibilities are
$\mathsf{Vol}_{p}(H)=\log\,\det(H)$ – which is directly proportional to the
classical volume definition, or $\mathsf{Vol}_{p}(H)=\rm{tr}\,H$ – which for
$p=2$ becomes the well known sum of ellipsoid semiaxes (see [47] and [43,
Chapter 8]). These measures can be easily generalized to non square matrices.
It suffices to compute the singular value decomposition. If $H=U\Sigma V^{T}$,
we could use the measures $\mathsf{Vol}_{p}(H)=\rm{tr}\,\Sigma$ or
$\mathsf{Vol}_{p}(H)=\log\,\det(\Sigma)$.
For non square matrices $H$, specific results for particular values of $p$ are
known. For example, we remind that if $p=\infty$ and
$H\in\mathbb{R}^{n_{\theta}\times s}$, $s\geq n_{\theta}$, then
$\theta_{c}\oplus H\mathbb{B}^{s}_{\infty}$ is a zonotope. Then, if we denote
as generator each of the columns of $H$, the volume of a zonotope can be
computed by means of a sum of terms (one for each different way of selecting
$n_{\theta}$ generators out of the $s$ generators of $H$); see [48], [49].
Another possible measure of the size of a zonotope $\theta_{c}\oplus
H\mathbb{B}^{s}_{\infty}$ is the Frobenious norm of $H$ [48].

Given an initial design set ${\mathbb{S}}_{N_{S}}$, we elect as our candidate
Scalable SAS the largest “volume” norm-based SAS contained in
${\mathbb{S}}_{N_{S}}$. Formally, this rewrites as the following optimization
problem

$\displaystyle\max\limits_{\theta_{c},H}~{}\mathsf{Vol}_{p}(H)$
$\displaystyle\text{subject to }\theta_{c}\oplus
H\mathbb{B}_{p}^{s}\subseteq{\mathbb{S}}_{N_{S}}$

As it has been shown, this problem is equivalent to

$\displaystyle\min\limits_{\theta_{c},H}$ $\displaystyle-\mathsf{Vol}_{p}(H)$
s.t. $\displaystyle
f_{\ell}^{T}(\tilde{w}^{(j)})\theta_{c}+\|H^{T}f_{\ell}(w^{(j)})\|_{p^{*}}-g_{\ell}(w^{(j)})\leq
0,$ $\displaystyle\qquad\qquad\qquad\ell\in[n_{\ell}],\;j\in[N_{S}],$

where we have replaced the maximization of $\mathsf{Vol}_{p}(H)$ with the
minimization of -$\mathsf{Vol}_{p}(H)$.

We notice that the constraints are convex on the decision variables; also, the
functional to minimize is convex under particular assumptions. For example
when $H$ is assumed to be square and positive definite and
$\mathsf{Vol}_{p}(H)=\log\det(H)$. For non square matrices, the constraints
remain convex, but the convexity of the functional to be minimized is often
lost. In this case, local optimization algorithms should be employed to obtain
a possibly sub-optimal solution.

(a) $\gamma=0.9701$

(b) $\gamma=1.5995$

(c) $\gamma=0.9696$

(d) $\gamma=1.5736$

Figure 4: Scaling procedure applied to (a) ${\mathbb{S}}_{1}$-SAS with
$N_{S}=100$, (b) ${\mathbb{S}}_{1}$-SAS with $N_{S}=1,000$ (b),
${\mathbb{S}}_{\infty}$-SAS with $N_{S}=100$ (c), and $\ell_{\infty}$-poly
with $N_{S}=1,000$ (d). The initial set is depicted in red, the final one in
green. The sampled design polytope ${\mathbb{S}}_{N_{S}}$ is represented in
black.

######  Example 4 (Norm-based SAS).

We revisit again Example 2 to show the use of norm-based SAS. We note that, in
this case, the designer can control the approximation outcome by acting upon
the number of design samples $N_{S}$ used for constructing the set
${\mathbb{S}}_{N_{S}}$. In Figure 4 we report two different norm-based SAS,
respectively with $p=1$ and $p=\infty$, and for each of them we consider two
different values of $N_{S}$, respectively $N_{S}=100$ and $N_{S}=1,000$.
Similarly to what observed for the sampled-polys, we see that for larger
$N_{S}$, the ensuing initial set becomes smaller. Consequently, we have an
inflating process for small $N_{S}$ and a shrinkage one for large $N_{S}$
However, we observe that in this case, the final number of inequalities is
independent on $N_{S}$, being equal to $3n_{\theta}+1=10$ for
${\mathbb{S}}_{1}$ and $2n_{\theta}$ for ${\mathbb{S}}_{\infty}$.

####  6.2.1 Relaxed computation

It is worth remarking that that the minimization problem of the previous
subsection might be infeasible. In order to guarantee the feasibility of the
problem, a soft-constrained optimization problem is proposed. With a relaxed
formulation, $\theta_{c}$ is not guaranteed to satisfy all the sampled
constraints. However $\theta_{c}\in{\mathbb{S}}_{N_{S}}$ is not necessary to
obtain an $\varepsilon$-CSS (in many practical applications, every element of
$\Theta$ has a non zero probability of violation and ${\mathbb{S}}_{N_{S}}$ is
empty with non-zero probability). Moreover, a relaxed formulation is necessary
to address problems in which there is no element of $\Theta$ with probability
of violation equal to zero (or significantly smaller than $\varepsilon$). Not
considering the possibility of violations is an issue especially when $N_{S}$
is large, because the probability of obtaining an empty sampled set
${\mathbb{S}}_{N_{S}}$ grows with the number of samples $N_{S}$.

Given $\xi>0$ the relaxed optimization problem is

$\displaystyle\min\limits_{\theta_{c},H,\tau_{1},\ldots,\tau_{N_{S}}}~{}-\mathsf{Vol}_{p}(H)+\xi\sum\limits_{j=1}^{N_{S}}\max\\{\tau_{j},0\\}$
(31) $\displaystyle\text{s.t.
}\;f_{\ell}^{T}(w^{(j)})\theta_{c}+\|H^{T}f_{\ell}(w^{(j)})\|_{p^{*}}-g_{\ell}(w^{(j)})\leq\tau_{j},$
$\displaystyle\qquad\qquad\qquad\ell\in[n_{\ell}],\;j\in[N_{S}].$

The parameter $\xi$ serves to provide an appropriate trade off between
satisfaction of the sampled constraints and the size of the obtained region. A
possibility to choose $\xi$ would be to choose it in such a way that the
fraction of violations $n_{viol}/N_{S}$ (where $n_{viol}$ is the number of
elements $\tau_{j}$ larger than zero) is smaller than $\varepsilon/2$.

##  7 Numerical example: Probabilistic set membership estimation

We now present a numerical example in which the results of the paper are
applied to the probabilistic set membership estimation problem, introduced in
subSection 2.3. We consider the universal approximation functions given by
Gaussian radial basis function networks (RBFN) [50].

Given the nodes $[x_{1},x_{2},\ldots,x_{M}]$ and the variance parameter $c$,
the corresponding Gaussian radial basis function network is defined as

${\rm{RBFN}}(x,\theta)=\theta^{T}\varphi(x),$

where
$\theta=\left[\begin{array}[]{ccc}\theta_{1}&\ldots&\theta_{M}\end{array}\right]^{T}$
represents the weights and

$\varphi(x)=\left[\begin{array}[]{ccc}\exp\left(\frac{-\|x-x_{1}\|^{2}}{c}\right)&\ldots&\exp\left(\frac{-\|x-x_{M}\|^{2}}{c}\right)\end{array}\right]^{T}$

is the regressor function. Given $\delta\in(0,1)$ and $\varepsilon\in(0,1)$,
the objective is to obtain, with probability no smaller than $1-\delta$, an
inner approximation of the probabilistic feasible parameter set
${\mathsf{FPS}}_{\varepsilon}$, which is the set of parameters
$\theta\in\mathbb{R}^{M}$ that satisfies

$\mathsf{Pr}_{\mathbb{W}}\\{|y-\theta^{T}\varphi(x)|\leq\rho\\}\geq
1-\varepsilon,$ (32)

where $x$ is a random scalar with uniform distribution in $[-5,5]$ and

$y=\sin(3x)+\sigma,$

where $\sigma$ is a random scalar with a normal distribution with mean $5$ and
variance 1.

We use the procedure detailed in Sections 4, 5 and 6 to obtain an SAS of
${\mathsf{FPS}}_{\varepsilon}$. We have taken a grid of $M=20$ points in the
interval $[-5,5]$ to serve as nodes for the RBFN, and a variance parameter of
$c=0.15$. We have taken $N_{S}=350$ random samples $w=(x,y)$ to compute the
initial geometry, which has been chosen to be an $\ell_{\infty}$ norm-based
SAS of dimension 20 with a relaxation parameter of $\xi=1$ (see (31)). The
chosen initial geometry is $\theta_{c}\oplus H\mathbb{B}^{20}_{\infty}$, where
$H$ is constrained to be a diagonal matrix.

When the initial geometry is obtained, we scale it around its center by means
of probabilistic scaling with Algorithm 1. The number of samples required for
the scaling phase to achieve $\varepsilon=0.05$ and $\delta=10^{-6}$ is
$N_{\gamma}=2065$ and the resulting scaling factor is $\gamma=0.3803$. The
scaled geometry $\theta_{c}\oplus\gamma H\mathbb{B}^{20}_{\infty}$ is, with a
probability no smaller than $1-\delta$, an inner approximation of
${\mathsf{FPS}}_{\varepsilon}$ which we will refer to as
${\mathsf{FPS}}_{\varepsilon}^{\delta}$. Since it is a transformation of an
$\ell_{\infty}$ norm ball with a diagonal matrix $H$, we can write it as

${\mathsf{FPS}}_{\varepsilon}^{\delta}=\\{\theta:\theta^{-}\leq\theta\leq\theta^{+}\\},$

where the extreme values $\theta^{-},\theta^{+}\in\mathbb{R}^{20}$ are
represented in Figure 5 [51], along with the central value
$\theta_{c}\in\mathbb{R}^{20}$.

Figure 5: Representation of the extreme values $\theta^{+}$ and $\theta^{-}$
and the central value $\theta_{c}$ of the
${\mathsf{FPS}}_{\varepsilon}^{\delta}$.

Once the ${\mathsf{FPS}}_{\varepsilon}^{\delta}$ has been computed, we can use
its center $\theta_{c}$ to make the point estimation
$y\approx\theta_{c}^{T}\varphi(x)$. We can also obtain probabilistic upper and
lower bounds of $y$ by means of equation (32). That is, every point in
${\mathsf{FPS}}_{\varepsilon}^{\delta}$ satisfies, with confidence $1-\delta$:

$\displaystyle\mathsf{Pr}_{\mathbb{W}}\\{y\leq\theta^{T}\varphi(x)+\rho\\}\geq
1-\varepsilon,$ (33)
$\displaystyle\mathsf{Pr}_{\mathbb{W}}\\{y\geq\theta^{T}\varphi(x)-\rho\\}\geq
1-\varepsilon.$

We notice that the tightest probabilistic bounds are obtained with
$\theta^{+}$ for the lower bound and $\theta^{-}$ for the upper one. That is,
we finally obtain that, with confidence $1-\delta$:

$\displaystyle\mathsf{Pr}_{\mathbb{W}}\\{y\leq{\theta^{-}}^{T}\varphi(x)+\rho\\}\geq
1-\varepsilon,$ (34)
$\displaystyle\mathsf{Pr}_{\mathbb{W}}\\{y\geq{\theta^{+}}^{T}\varphi(x)-\rho\\}\geq
1-\varepsilon.$

Figure 6 shows the results of both the point estimation and the probabilistic
interval estimation.

Figure 6: Real values of $y$ vs central estimation (blue) and interval
prediction bounds (red).

##  8 Conclusions, extensions, and future directions

In this paper, we proposed a general approach to construct probabilistically
guaranteed inner approximations of the chance-constraint set
$\mathbb{X}_{\varepsilon}$. The approach is very general and flexible.

First, we remark that the proposed scaling approach is not limited to sets
defined by linear inequalities, but immediately extends to more general sets.
Indeed, we may consider a generic binary performance function
$\phi:\Theta\times\mathbb{W}\to\\{0,\,1\\}$ defined as 444Clearly, this
formulation encompasses the setup discussed, obtained by simply setting
$\phi(\theta,w)=\left\\{\begin{array}[]{ll}0&\text{if $F(w)\theta\leq
g(w)$}\\\ 1&\text{otherwise.}\end{array}\right.$

$\phi(\theta,q)=\left\\{\begin{array}[]{ll}0&\text{if $\theta$ meets design
specifications for $w$}\\\ 1&\text{otherwise.}\end{array}\right.$ (35)

In this case, the violation probability may be written as
$\mathsf{Viol}(\theta)\doteq\mathsf{Pr}_{\mathbb{W}}\,\\{\,\psi(\theta,w)=1\,\\}=\mathbb{E}(\theta)$,
and we can still define the set $\mathbb{X}_{\varepsilon}$ as in (5). Then,
given an initial SAS candidate, Algorithm 1 still provides a valid
approximation. However, it should be remarked that, even if we choose a “nice”
SAS as those previously introduced, the nonconvexity of $\phi$ will most
probably render step 4 of the algorithm intractable. To further elaborate on
this point, let us focus on the case when the design specification may be
expressed as a (nonlinear) inequality of the form

$\psi(\theta,q)\leq 0.$

Then, step 4 consist in solving the following nonconvex optimization problem

$\displaystyle\gamma_{i}\doteq$ $\displaystyle\arg\max\gamma$ (36)
$\displaystyle\text{s.t.}\quad{\mathbb{S}}(\gamma)\subseteq\mathbb{X}(w^{(i)})=\Bigl{\\{}\theta\in\Theta\;|\;\psi(\theta,w^{(i)})\leq
0\Bigr{\\}}.$

We note that this is general a possibly hard problem. However, there are cases
when this problem is still solvable. For instance, whenever $\psi(\theta,q)$
is a convex function of $\theta$ for fixed $w$ and the set ${\mathbb{S}}$ is
also convex, the above optimization problem may be formulated as a convex
program by application of Finsler lemma. We remark that, in such situations,
the approach proposed here is still completely viable, since all the
derivations continue to hold.

Second, we remark that the paper open the way to the design of other families
of Scaling SAS. For instance, we are currently working on using the family of
sets defined in the form of polynomial superlevel sets (PSS) proposed in [52].

##  Appendix A Appendix

###  A.1 Proof of Lemma 1

To prove the lemma, we first recall the following definition from [41].

######  Definition 4 ($(\alpha,k)$-Boolean Function).

The function $h:\Theta\times\mathbb{W}\to\mathbb{R}$ is an
$(\alpha,k)$-Boolean function if for fixed $w$ it can be written as an
expression consisting of Boolean operators involving $k$ polynomials
$p_{1}(\theta),p_{2}(\theta),\ldots,p_{k}(\theta),$ in the components
$\theta_{i}$, $i\in[n_{\theta}]$ and the degree with respect to $\theta_{i}$
of all these polynomials is no larger than $\alpha$.

Let us now define the binary functions

$h_{\ell}(\theta,w)\doteq\left\\{\begin{array}[]{rl}0&\mbox{ if
}f_{\ell}(w)\theta\leq g_{\ell}(w)\\\ 1&\mbox{
otherwise}\end{array}\right.,\;\ell\in[n_{\ell}].$

Introducing the function
$h(\theta,w)\doteq\max\limits_{\ell=1,\ldots,n_{\ell}}h_{\ell}(\theta,w),$ we
see that the violation probability can be alternatively written as
$\mathsf{Viol}(\theta)\doteq\mathsf{Pr}_{\mathbb{W}}\,\\{\,h(\theta,w)=1\,\\}.$
The proof immediately follows by observing that $h(\theta,w)$ is an
$(1,n_{\ell})$-Boolean function, since it can be expressed as a function of
$n_{\ell}$ Boolean functions, each of them involving a polynomial of degree 1.
Indeed, it is proven in [41, Theorem 8], that, if
$h:\Theta\times\mathbb{W}\to\mathbb{R}$ is an $(\alpha,k)$-Boolean function
then, for $\varepsilon\in(0,0.14)$, with probability greater than $1-\delta$
we have $\mathsf{Pr}_{\mathbb{W}}\,\\{\,h(\theta,w)=1\,\\}\leq\varepsilon$ if
$N$ is chosen such that

$N\geq\frac{4.1}{\varepsilon}\Big{(}\ln\frac{21.64}{\delta}+4.39n_{\theta}\,\log_{2}\Big{(}\frac{8e\alpha
k}{\varepsilon}\Big{)}\Big{)}.$

###  A.2 Property 1

######  Property 1.

Given $\varepsilon\in(0,1)$, $\delta\in(0,1)$, and $0\leq r\leq N$, let $N$ be
such that $\mathbf{B}(r;N,\varepsilon)\leq\delta$. Draw $N$ iid sample-sets
$\\{\mathbb{X}^{(1)},\mathbb{X}^{(2)},\ldots,\mathbb{X}^{(N)}\\}$ from a
distribution $\mathsf{Pr}_{\mathbb{X}}$. For $i\in[N]$, let
$\gamma_{i}\doteq\gamma(\mathbb{X}^{(i)})$, with $\gamma(\cdot)$ as in
Definition 3, and suppose that $\bar{\gamma}=\gamma_{1+r:N}>0$. Then, with
probability no smaller than $1-\delta$, it holds that
$\mathsf{Pr}_{\mathbb{X}}\\{\theta_{c}\oplus\bar{\gamma}{\mathbb{S}}\not\subseteq\mathbb{X}\\}\leq\varepsilon$.

Proof: It has been proven in [38, 39] that if one discards no more than $r$
constraints on a convex problem with $N$ random constraints, then the
probability of violating the constraints with the solution obtained from the
random convex problem is no larger than $\varepsilon\in(0,1)$, with
probability no smaller than $1-\delta$, where

$\delta=\left(\begin{array}[]{c}d+r-1\\\ d-1\\\
\end{array}\right)\sum\limits_{i=0}^{d+r-1}\left(\begin{array}[]{c}N\\\ i\\\
\end{array}\right)\varepsilon^{i}(1-\varepsilon)^{N-i},$

and $d$ is the number of decision variables. We apply this result to the
following optimization problem

$\max\limits_{\gamma}\gamma\text{ subject to
}\theta_{c}\oplus\gamma{\mathbb{S}}\subseteq\mathbb{X}^{(i)},\;\;i\in[N].$

From Definition 3, we could rewrite this optimization problem as

$\max\limits_{\gamma}\gamma\text{ subject to
}\gamma\leq\gamma(\mathbb{X}^{(i)}),\;i\in[N].$

We first notice that the problem under consideration is convex and has a
unique scalar decision variable $\gamma$. That is, $d=1$. Also, the non-
degeneracy and uniqueness assumption required in the application of the
results of [38] and [39] are satisfied. Hence, if we allow $r$ violations in
the above minimization problem, we have that with probability no smaller than
$1-\delta$, where

$\delta=\left(\begin{array}[]{c}r\\\ 0\\\
\end{array}\right)\sum\limits_{i=0}^{r}\left(\begin{array}[]{c}N\\\ i\\\
\end{array}\right)\varepsilon^{i}(1-\varepsilon)^{N-i}=\mathbf{B}(r;N,\varepsilon),$

the solution $\bar{\gamma}$ of problem (A.2) satisfies
$\mathsf{Pr}_{\mathbb{X}}\\{\bar{\gamma}>\gamma(\mathbb{X})\\}\leq\varepsilon.$
We conclude from this, and Definition 3, that with probability no smaller than
$1-\delta$,

$\mathsf{Pr}_{\mathbb{X}}\\{\theta_{c}\oplus\bar{\gamma}{\mathbb{S}}\not\subseteq\mathbb{X}\\}\leq\varepsilon.$

Finally, note that the optimization problem under consideration can be solved
directly by ordering the values $\gamma_{i}=\gamma(\mathbb{X}^{(i)})$. It is
clear that if $r\geq 0$ violations are allowed, then the optimal value for
$\gamma$ is $\bar{\gamma}=\gamma_{r+1:N}$. ∎

###  A.3 Proof of Theorem 2

Note that, by definition, the condition $\theta_{c}\oplus\gamma
H\mathbb{B}^{s}_{p}\subseteq\mathbb{X}(w)$ is equivalent to

$\max\limits_{z\in\mathbb{B}^{s}_{p}}f_{\ell}^{T}(w)(\theta_{c}+\gamma
Hz)-g_{\ell}(w)\leq 0,\;\ell\in[n_{\ell}].$

Equivalently, from the dual norm definition, we have

$f_{\ell}^{T}(w)\theta_{c}+\gamma\|H^{T}f_{\ell}(w)\|_{p^{*}}-g_{\ell}(w)\leq
0,\;\ell\in[n_{\ell}].$

Denote by $\gamma_{\ell}$ the scaling factor $\gamma_{\ell}$ corresponding to
the $\ell$-th constraint

$f_{\ell}^{T}(w)\theta_{c}+\gamma_{\ell}\|H^{T}f_{\ell}(w)\|_{p^{*}}-g_{\ell}(w)\leq
0.$

With the notation introduced in the Lemma, this constraint rewrites as

$\gamma_{\ell}\rho_{\ell}(w)\leq\tau_{\ell}(w).$

The result follows noting that the corresponding scaling factor
$\gamma_{\ell}(w)$ can be computed as

$\gamma_{\ell}(w)=\max_{\gamma_{\ell}\rho_{\ell}(w)\leq\tau_{\ell}(w)}\gamma_{\ell},$

and that the value for $\gamma(w)$ is obtained from the most restrictive one.
∎

## References

  * [1] A. Prékopa, _Stochastic Programming_. Springer Science & Business Media, 2013.
  * [2] N. V. Sahinidis, “Optimization under uncertainty: state-of-the-art and opportunities,” _Computers & Chemical Engineering_, vol. 28, no. 6-7, pp. 971–983, 2004.
  * [3] A. Ben-Tal and A. Nemirovski, “Robust convex optimization,” _Mathematics of Operations Research_ , vol. 23, pp. 769–805, 1998.
  * [4] G. Calafiore, F. Dabbene, and R. Tempo, “Research on probabilistic methods for control system design,” _Automatica_ , vol. 47, pp. 1279–1293, 2011.
  * [5] R. Tempo, G. Calafiore, and F. Dabbene, _Randomized Algorithms for Analysis and Control of Uncertain Systems: with Applications_. Springer Science & Business Media, 2012.
  * [6] M. Mammarella, E. Capello, F. Dabbene, and G. Guglieri, “Sample-based SMPC for tracking control of fixed-wing UAV,” _IEEE Control Systems Letters_ , vol. 2, no. 4, pp. 611–616, 2018.
  * [7] J. Li, W. Zhan, Y. Hu, and M. Tomizuka, “Generic tracking and probabilistic prediction framework and its application in autonomous driving,” _IEEE Transactions on Intelligent Transportation Systems_ , 2019.
  * [8] M. Chamanbaz, F. Dabbene, and C. Lagoa, _Algorithms for Optimal AC Power Flow in the Presence of Renewable Sources_. Wiley Encyclopedia of Electrical and Electronics Engineering, 2020, pp. 1–13.
  * [9] M. Chamanbaz, F. Dabbene, and C. M. Lagoa, “Probabilistically robust AC optimal power flow,” _IEEE Transactions on Control of Network Systems_ , vol. 6, no. 3, pp. 1135–1147, 2019.
  * [10] X. Geng and L. Xie, “Data-driven decision making in power systems with probabilistic guarantees: Theory and applications of chance-constrained optimization,” _Annual Reviews in Control_ , vol. 47, pp. 341–363, 2019\.
  * [11] M. Lorenzen, F. Dabbene, R. Tempo, and F. Allgöwer, “Stochastic MPC with offline uncertainty sampling,” _Automatica_ , vol. 81, no. 1, pp. 176–183, 2017.
  * [12] M. Mammarella, M. Lorenzen, E. Capello, H. Park, F. Dabbene, G. Guglieri, M. Romano, and F. Allgöwer, “An offline-sampling SMPC framework with application to autonomous space maneuvers,” _IEEE Transactions on Control Systems Technology_ , pp. 1–15, 2018.
  * [13] T. Alamo, V. Mirasierra, F. Dabbene, and M. Lorenzen, “Safe approximations of chance constrained sets by probabilistic scaling,” in _2019 18th European Control Conference (ECC)_. IEEE, 2019, pp. 1380–1385.
  * [14] M. Mammarella, T. Alamo, F. Dabbene, and M. Lorenzen, “Computationally efficient stochastic mpc: a probabilistic scaling approach,” in _Proc. of 4th IEEE Conference on Control Technology and Applications_ , 2020.
  * [15] M. Ahsanullah, V. Nevzorov, and M. Shakil, _An introduction to Order Statistics_. Paris: Atlantis Press, 2013\.
  * [16] B. Miller and H. Wagner, “Chance constrained programming with joint constraints,” _Operations Research_ , vol. 13, pp. 930–945, 1965.
  * [17] L. Khachiyan, “The problem of calculating the volume of a polyhedron is enumerably hard,” _Russian Mathematical Surveys_ , 1989.
  * [18] A. Shapiro, D. Dentcheva, and A. Ruszczyński, _Lectures on stochastic programming: modeling and theory_. SIAM, 2014.
  * [19] W. van Ackooij, “Eventual convexity of chance constrained feasible sets,” _Optimization_ , vol. 64, no. 5, pp. 1263–1284, 2015.
  * [20] A. Prékopa, T. Rapcsák, and I. Zsuffa, “Serially linked reservoir system design using stochastic programing,” _Water Resources Research_ , vol. 14, no. 4, 1978.
  * [21] D. Dentcheva, B. Lai, and A. Ruszczyński, “Dual methods for probabilistic optimization problems*,” _Mathematical Methods of Operations Research_ , vol. 60, no. 2, pp. 331–346, 2004.
  * [22] M. Lorenzen, F. Dabbene, R. Tempo, and F. Allgöwer, “Constraint-tightening and stability in stochastic model predictive control,” _IEEE Transactions on Automatic Control_ , vol. 62, no. 7, pp. 3165–3177, 2017.
  * [23] A. Vicino and G. Zappa, “Sequential approximation of feasible parameter sets for identification with set membership uncertainty,” _IEEE Transactions on Automatic Control_ , vol. 41, no. 6, pp. 774–785, 1996.
  * [24] J. M. Bravo, T. Alamo, and E. F. Camacho, “Bounded error identification of systems with time-varying parameters,” _IEEE Transactions on Automatic Control_ , vol. 51, no. 7, pp. 1144–1150, 2006.
  * [25] S. Kataoka, “A stochastic programming model,” _Econometrica: Journal of the Econometric Society_ , pp. 181–196, 1963.
  * [26] A. Prékopa, “Logarithmic concave measures with application to stochastic programming,” _Acta Scientiarum Mathematicarum_ , pp. 301–316, 1971.
  * [27] C. M. Lagoa, “On the convexity of probabilistically constrained linear programs,” in _Proceedings of the 38th IEEE Conference on Decision and Control (Cat. No.99CH36304)_ , vol. 1, 1999, pp. 516–521 vol.1.
  * [28] G. C. Calafiore and L. E. Ghaoui, “On distributionally robust chance-constrained linear programs,” _Journal of Optimization Theory and Applications_ , vol. 130, no. 1, pp. 1–22, 2006.
  * [29] R. Henrion and C. Strugarek, “Convexity of chance constraints with independent random variables,” _Computational Optimization and Applications_ , vol. 41, no. 2, pp. 263–276, 2008.
  * [30] L. Hewing and M. N. Zeilinger, “Stochastic model predictive control for linear systems using probabilistic reachable sets,” in _2018 IEEE Conference on Decision and Control (CDC)_ , 2018, pp. 5182–5188.
  * [31] S. Yan, P. Goulart, and M. Cannon, “Stochastic model predictive control with discounted probabilistic constraints,” in _2018 European Control Conference (ECC)_. IEEE, 2018, pp. 1003–1008.
  * [32] A. Nemirovski and A. Shapiro, “Convex approximations of chance constrained programs,” _SIAM Journal on Optimization_ , vol. 17, no. 4, pp. 969–996, 2006.
  * [33] W. Chen, M. Sim, J. Sun, and C.-P. Teo, “From CVaR to uncertainty set: Implications in joint chance-constrained optimization,” _Operations Research_ , vol. 58, no. 2, pp. 470–485, 2010.
  * [34] A. Jasour, N. S. Aybat, and C. M. Lagoa, “Semidefinite programming for chance constrained optimization over semialgebraic sets,” _SIAM Journal on Optimization_ , vol. 25, no. 3, pp. 1411–1440, 2015.
  * [35] J. B. Lasserre, “Representation of chance-constraints with strong asymptotic guarantees,” _IEEE Control Systems Letters_ , vol. 1, no. 1, pp. 50–55, 2017\.
  * [36] G. Calafiore and M. Campi, “The scenario approach to robust control design,” _IEEE Transactions on Automatic Control_ , vol. 51, no. 5, pp. 742–753, 2006\.
  * [37] M. Campi and S. Garatti, “The exact feasibility of randomized solutions of robust convex programs,” _SIAM Journal of Optimization_ , vol. 19, pp. 1211—1230, 2008.
  * [38] G. Calafiore, “Random convex programs,” _SIAM Journal of Optimization_ , vol. 20, pp. 3427–3464, 2010.
  * [39] M. Campi and S. Garatti, “A sampling-and-discarding approach to chance-constrained optimization: feasibility and optimality,” _Journal of Optimization Theory and Applications_ , vol. 148, pp. 257–280, 2011.
  * [40] V. Vapnik, _Statistical Learning Theory_. New York: John Wiley and Sons, 1998.
  * [41] T. Alamo, R. Tempo, and E. F. Camacho, “Randomized strategies for probabilistic solutions of uncertain feasibility and optimization problems,” _IEEE Transactions on Automatic Control_ , vol. 54, no. 11, pp. 2545–2559, 2009.
  * [42] T. Alamo, R. Tempo, A. Luque, and D. Ramirez, “Randomized methods for design of uncertain systems: Sample complexity and sequential algorithms,” _Automatica_ , vol. 52, pp. 160–172, 2015.
  * [43] S. Boyd and L. Vandenberghe, _Convex Optimization_. Cambridge University Press, 2004.
  * [44] M. Herceg, M. Kvasnica, C. N. Jones, and M. Morari, “Multi-parametric toolbox 3.0,” in _2013 European control conference (ECC)_. IEEE, 2013, pp. 502–510.
  * [45] F. Dabbene, C. Lagoa, and P. Shcherbakov, “On the complexity of randomized approximations of nonconvex sets,” in _2010 IEEE International Symposium on Computer-Aided Control System Design_. IEEE, 2010, pp. 1564–1569.
  * [46] V. T. H. Le, C. Stoica, T. Alamo, E. F. Camacho, and D. Dumur, _Zonotopes: From Guaranteed State-estimation to Control_. Wiley, 2013.
  * [47] F. Dabbene, D. Henrion, C. Lagoa, and P. Shcherbakov, “Randomized approximations of the image set of nonlinear mappings with applications to filtering,” _IFAC-PapersOnLine_ , vol. 48, no. 14, pp. 37–42, 2015.
  * [48] T. Alamo, J. M. Bravo, and E. F. Camacho, “Guaranteed state estimation by zonotopes,” _Automatica_ , vol. 41, no. 6, pp. 1035–1043, 2005.
  * [49] E. Gover and N. Krikorian, “Determinants and the volumes of parallelotopes and zonotopes,” _Linear Algebra and its Applications_ , vol. 433, no. 1, pp. 28–40, 2010.
  * [50] M. D. Buhmann, “Radial basis functions,” _Acta numerica_ , vol. 9, pp. 1–38, 2000.
  * [51] L. J, “Plotrix: a package in the red light district of r,” _R-News_ , vol. 6, no. 4, pp. 8–12, 2006.
  * [52] F. Dabbene, D. Henrion, and C. M. Lagoa, “Simple approximations of semialgebraic sets and their applications to control,” _Automatica_ , vol. 78, pp. 110 – 118, 2017.