{'fraction_non_alphanumeric': 0.0393981894683157, 'fraction_numerical': 0.008882655446470314, 'mean_word_length': 4.458130364184644, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 4, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Randomized Approximations of the Image Set of Nonlinear Mappings with Applications to Filtering

The aim of this paper is twofold: In the first part, we leverage recent results on scenario design to develop randomized algorithmsfor approximating the image set of a nonlinear mapping, that is, a (possibly noisy) mapping of a set via a nonlinear function.We introduce minimum-volume approximations which have the characteristic of guaranteeing a low probability of violation, i.e.,we admit for a probability that some points in the image set are not contained in the approximating set,but this probability is kept below a pre-specified threshold.In the second part of the paper, this idea is then exploited to develop a new family of randomized prediction-corrector filters.These filters represent a natural extension and rapprochement of Gaussian and set-valued filters,and bear similarities with modern tools such as particle filters.


Introduction
In recent years, randomized algorithms have gained increasing popularity in the field of control of uncertain systems; e.g., see Calafiore et al. (2011); Tempo et al. (2013), due to their ability of dealing with large and complex uncertainty structures, thus extending the applicability of the robust control methods. This is obtained by shifting from the robustness paradigm, where one looks for guaranteed performance which should hold for every possible instance of the uncertainty, to an approach where probabilistic guarantees are accepted, i.e., performance is guaranteed only within a given level of probability ε > 0. The main technical tool that permits to obtain computationally tractable solutions are randomized algorithms, which could be seen as extensions of the classical Monte Carlo method, and are based on the extraction of random samples of the uncertainty.
In this paper, we exploit these ideas for finding reliable approximations of the image of a set through a nonlinear mapping, which we refer to as the image set. This set is in general nonconvex (possibly not connected), so that classical approximations may be rather difficult to compute and in general may turn out to be very conservative.
The first part of the paper adapts and significantly extends recent results presented by some of the authors in Dabbene et al. (2010), where a new definition of "goodness" of approximation was provided in probabilistic terms. Namely, an approximating set A of a set X is deemed to be "good" if it contains "most" of the points in X with high probability. Contrary to classical approximating sets, which are generally either outer or inner, the ensuing approximating set "optimally describes" the set without neither containing nor being contained in it. This new concept allows to obtain generally tighter approximations, providing a probabilistic characterization of the set, which is particularly appealing in many contexts (even if it may not be desired in others, such e.g. as safety analysis). We recall that outer bounding sets, that is, approximations which are guaranteed to contain the set X , have been very popular in the set-membership approach, and have been used in designing set-theoretic state estimators for uncertain discrete-time nonlinear dynamic systems; e.g., see El Ghaoui and Calafiore (2001); Alamo et al. (2008). Inner ellipsoidal approximations have been introduced for instance in the context of nonlinear programming problems (Nesterov and Nemirovski, 1994) and tolerance design (Wojciechowski and Vlach, 1993).
In Sections II and III, the results of Dabbene et al. (2010) are particularized to the specific problem at hand, and are also generalized by considering a new family of approximating sets which are based on the construction of minimum volume polynomial superlevel sets, recently introduced in Dabbene and Henrion (2013). With these approximating sets, the original convexity requirement can be relaxed, since they can be nonconvex/nonconnected, thus allowing for better descriptions.
The second part of the paper extends these ideas to the design of probabilistic filters for nonlinear discrete time dynamical systems subject to random uncertainty. In particular, the image set approximation is used to design a probabilistically optimal predictor filter. The rationale behind this approach is to "describe" the state position at step k by a set where the state lies with high probability. The prediction step is then combined with a correction step where the propagated set is trimmed based on the available measurements at time k and on the measurement noise assumptions. A detailed discussion on the proposed randomized prediction correction filter is given in Section IV, where a discussion on the related literature is also reported. Numerical examples conclude the paper.
Notation: For a symmetric matrix P , P 0 means that P is positive definite. We denote by B p the unit ball in the p norm: B p . = {z ∈ R n : z p ≤ 1} . The volume or, more precisely, the Lebesgue measure of a compact set X is denoted by vol X . = X dx. The uniform measure over X is denoted by λ X , i.e. λ X is such that, for any set Y ⊆ X , λ X (Y) = vol Y/vol X . The set of all polynomials of order less than or equal to σ is denoted by P σ . The monomial basis for this set is represented by the (column) vector π σ (x) and any polynomial q ∈ P σ can be expressed in the form q(x) = π σ (x)q = π σ/2 (x)Qπ σ/2 (x) where q is a vector and Q is a symmetric matrix of appropriate dimensions, referred to as the Gram matrix.

Minimum volume approximations
Consider the following nonlinear mapping: which represents the one-step evolution of a discrete-time dynamical system, with x, x + ∈ R n representing the current and future states, respectively, and w ∈ R nw describing a process noise vector. We assume that the current state is confined within a compact set X ⊂ R n and that the noise w also belongs to a given compact set W. The problem we are interested in is to find a good approximation to the set of points that can be obtained from (1) by starting from x ∈ X and accounting for all possible values of the noise w ∈ W, that is, find an approximation of the image set defined as Motivated by the computational complexity issues and the fact that deterministic formulations of the approximation problem may not be suitable in practical situations, in Dabbene et al. (2010) the authors proposed an original approach for tackling the problem, based on a probabilistic viewpoint. To this end, a probabilistic information over X and W is assumed to be known (which is actually the case in many applications), and an approximation A of X + is deemed to be "good" if it contains all points in X + with high probability. More formally, we assume that the sets X and W are endowed with probability measures µ X and µ W . Then, for a given reliability level ε, the following concept of ε-probabilistic approximation is introduced.
Definition 1 (ε-probabilistic approximation of X + ) The set A is an ε-probabilistic approximation of the set X + if Note that the probability in (2) is measured with respect to the underlying measures µ X and µ W . The left-hand side of the equation is referred to as the violation probability of the set A.
The main characteristic of this approach is that an ε-probabilistic approximating set has neither to cover nor to be fully contained in X ; it just has to guarantee that the violation probability of the set A is bounded by ε. Clearly, we are interested in finding the smallest among such sets. In the sequel, we first define the two families of approximating sets considered in this paper.

Convex approximating sets: Ellipsoids, parallelotopes, and hyperrectangles
The following general description of the norm-based approximating set (NAS) was introduced in Dabbene et al. (2010): Note that the family of sets above is parameterized by the positive-definite shape matrix P , and by the center c ∈ R n ; it represents a generalization of the classical ellipsoidal set for norms different from the Euclidean one.
Indeed, for p = 2, we obtain the ellipsoid and for p = ∞ we get a so-called elementary parallelotope, a special-type polytope with parallel faces (these can be viewed as a particular class of zonotopes (Alamo et al., 2005) with positive definite generator matrix). In particular, if P is chosen to be diagonal, we obtain a classical hyper-rectangle.
It follows that, in general, the problem of finding the minimum volume NAS containing X + can be rewritten in the form of the following robust convex problem: (NAS-robust)

Nonconvex approximations: Polynomial superlevel sets
In Dabbene and Henrion (2013), nonconvex set approximations based on the superlevel set of a multidimensional polynomial have been introduced, and shown to represent a simple and efficient way for describing complex shaped sets. Formally, assume we are given a compact semialgebraic set being given polynomials (the set S is usually a hyperrectangle). Then, given a polynomial of degree σ > 0, i.e., q ∈ P σ , its polynomial superlevel set U(q) is defined as In the sequel, we will refer to this family of sets as polynomial-based approximating sets (PAS). The goal is then to find a polynomial q ∈ P σ whose corresponding PAS U(q) contains X + and has minimal volume. In Dabbene and Henrion (2013), the following problem was formulated: and the L 1 measure above was shown to be a good surrogate to the volume of the set U(q). Moreover, theoretically rigorous convergence results were provided in Dabbene and Henrion (2013), showing that, if X + is semialgebraic, a hierarchy of outer approximations obtained by solving the above problem for increasing values of σ converges in volume (or, equivalently, almost uniformly) to the set X + .
We now make the important observation that this result still holds for generic compact sets, not necessarily of semialgebraic form. Hence, we propose the following family of PASs, where we further relax the positivity requirement of q over S by requiring that the polynomial can be expressed as a sum-of-squares (SOS) over S: (PAS-robust) Note that the requirement that q is SOS over S can be formulated using Putinar's Positivstellensatz and a hierarchy of finite-dimensional convex LMI relaxations which are linear in the coefficients of q. More specifically, we write q = r 0 + m i r i b i , where r 0 , r 1 , . . . , r m are polynomial sum-of-squares of given degree σ r , to be found. For simplicity, we choose σ r = σ. For each fixed degree, the problem of finding such polynomials is an LMI; e.g., see (Lasserre, 2015, Section 3.2). Note again that problem (PAS-robust) is a robust convex program (that is, for every fixed x + , we get a convex program).

Randomized approximations
In this section, we propose a simple randomized algorithm to construct, with arbitrarily high probability, ε-probabilistic approximations of the image set X + . The suggested randomized procedure is rather straightforward: We draw N samples of x ∈ X according to the measure µ X , obtaining the multisample and, similarly, we draw N instances of the random noise w ∈ W, drawn according to the (given) noise measure µ W Based on these samples, we compute the points Then, we construct the approximating set A as the minimum-size bounding set containing the random points x + | N 1 thus obtained. This is doable by solving one of the following two scenario problems 1 : and 1 Note the slight abuse of notation adopted: A NAS denotes the NAS set defined by the optimal values of P and c obtained by solving problem (NAS-random). Same for A PAS . (PAS-random) The idea at the basis of the approach is depicted in Fig. 2, where random samples of x and w are drawn and mapped via the nonlinear function f (·, ·). Then, an approximation of X + is obtained by constructing a minimum volume approximating set around these points. Note that both problems above represent simple convex optimization problems which can be efficiently solved by available SDP solvers. In particular: (a) problem (NAS-random) reduces to a linear program (LP) for p = 1 and diagonal P ; (b) the constraints q(x (i) + ) ≥ 1 in problem (PAS-random) can also be expressed as an LP in the coefficients of q(x). The procedure is summarized in the following algorithm.
The key features of the algorithm above are as follows: (i) there is no need in knowing the probability measures µ X and µ W explicitly, only the capability of obtaining random samples from these are mandatory for running the algorithm; (ii) formal results on the probabilistic quality of the obtained approximations are available. In particular, the properties of the algorithm follow from the well-known results of Calafiore and Campi (2006) on the so-called scenario optimization, which allow to bound a priori (before the algorithm is run) the violation probability of the solution.
Theorem 1 Let A (i.e., A NAS or A PAS ) be the outcome of Algorithm 1. Fix a risk level ε > 0 and let d be the "size" of the optimization problem NAS-random (or PAS-random). Then we have Moreover, the above bound is tight, see Campi and Garatti (2008); Calafiore (2010).
Note that the result above makes explicit use of the size d of the convex problem. This represents the number of free variables in the convex optimization problem; see Table 1. In order to obtain an expression for the desired sample size N (the number of samples ensuring that the designed approximation is ε-optimal with probability at least δ), we need to invert equation (7). The best bound currently available is proposed in Alamo et al. (2015), where it is shown that, given ε, δ ∈ (0, 1), if N is chosen as then A is an ε-probabilistic approximation of X + with probability at least (1 − δ).

Numerical Example
To show the approximating sets obtained by Algorithm 1, we consider the nonlinear mapping x + (1) = sin x(2) + 3 cos x(2) + w(1), Figure 3 depicts N = 200 randomly generated points in X + , and the corresponding NAS and PAS probabilistic approximations. The shaded area represents the "true" set X + obtained by generating one million points in it. These approximations are indeed seen to be neither inner nor outer ones.

Randomized Filtering
The results of the previous sections are here reinterpreted in a filtering setting. The idea is rather intuitive: One can apply the image set description in an iterative way and use the approximation as a description of the state position at time k. For the sake of simplicity, we limit our exposition to the case where NAS descriptions are employed.
Consider the following discrete-time (autonomous) nonlinear system Equation (9) represents the systems dynamics, and x k , x k+1 ∈ R n denote the state at time instants k and k + 1, respectively, while w k ∈ R nw is a noise vector, usually referred to as process noise. Equation (10) represents a nonlinear measurement mapping associated with the system equations, where y k ∈ R ny is the measurement vector and v k ∈ R nv is the measurement noise vector.
Various approaches have been adopted in the literature to address the filtering problem for this type of systems. In the Extended Kalman filtering approach, and, in general, in all the subsequently developed Gaussian Filtering techniques, (e.g., see Ito. and Xiong (2000)), probabilistic information on the states and noise is considered. The main approximation made in this setup is that all information about the state and noise is collected in the first two moments of this distribution; in other words, the state and (update and measurement) noise can be described by a Gaussian distribution. Note that, even if at k = 0 everything could indeed be Gaussian, this is clearly not true after the first step, since Gaussian distributions are generally not preserved under nonlinear mappings. Approaches of this line are also those based on Unscented Kalman Filtering, e.g., see Sarkka (2007).
A philosophically different approach is the one based on set-valued filtering, in which no probabilistic assumptions are made on the nature of the noise and on the a priori information on x; instead, a so-called unknown-but-bounded approach is used. In those works (e.g., see Calafiore (2005)), the a priori information is given in terms of a bounded set where the initial state is guaranteed to lie, and at each time step both process and measurement noise are assumed to take values from given bounded sets. Then, all information on the predicted state is itself captured by a set (usually an ellipsoid) which is guaranteed to contain the state at time k + 1.
In this work, we propose a rapprochement between those two models and assume that the process and measurement noise are random, with given distributions over bounded support sets. That is, we assume that, at each time instant k, w k ∈ W and v k ∈ V, where W, V are bounded sets 2 . Moreover, we assume that w k has stochastic nature, and we denote by P W the probability measure on W. Note that no stochasticity assumption is made on v k .
The next approximation is crucial in our setup, and, from the philosophical point of view, it plays the same role of the Gaussian approximation commonly made in Gaussian filtering, in which the first two moments of the distribution are taken as representatives of the whole pdf of the propagated state.
Approximation 1 (State description) At time k, all state information is captured by a NAS set A k = A(P k , c k ) (e.g., an ellipsoid) with center c k and shape matrix P k . Moreover, we assume that x k has stochastic nature, and it is uniformly distributed over A k .
Note again that this constitutes one of the main approximations introduced in our setup: At each step we 'reset' the density information about the state and make the implicit assumption that x k has a uniform distribution λ with ellipsoidal support A k .

Probabilistically Guaranteed Simulation
In this section, we concentrate on the state dynamics equation only and propose a prediction filter for simulating the behavior of system (9). To this end, assume that at step k, the state x k is guaranteed (in probability) to lie inside a NAS set, that is x k ∈ A k . = A(x k , P k ). Then, we can apply the set image approximation algorithm; namely, given a probabilistic parameter ε, we choose N so as to guarantee a desired approximation level, and we start by drawing uniform random samples x Based on these samples, we compute the points x . . , N , and solve the NAS-random problem. The proposed procedure is summarized in Algorithm 2.

Algorithm 2 One-step Randomized Prediction Filter
Remark 1 (m-step ahead prediction) The iterative one-step prediction scheme discussed above has been introduced having in mind the prediction-correction scheme that we formulate in the next section. However, there are cases when it could be convenient to directly perform an m-step ahead prediction, in one shot. To do this, a simple solution is to "propagate" the state samples m-times ahead in the following way: Generate N samples x

Measurement update
If a measurement is available, it can be immediately used to reject incorrect predictions, as discussed in Algorithm 3 (RPCF) below.
Remark 2 (Similarities with particle filtering) In the so called particle filtering, particles are generated and fed through the dynamical system. Then, they are kept/discarded based on Bayesian-type reasonings. Our setup has strong similarities. Indeed, one can construct a "particle-like" implementation of the predictor-corrector randomized algorithm, where at each step not all samples are generated again, but instead, the non-discarded ones are fed through. This algorithm seems to provide a good practical behavior. However, it still needs a thorough theoretical analysis, which will be the subject of future work.

Numerical Example
As a filtering example, we tested the specific two-dimensional nonlinear system with scalar measurements, borrowed from Alamo et al. (2008), which was considered in the context of set filtering by means of zonotopic approximation sets: where the propagation noise is |w k (1)| ≤ 0.1, |w k (2)| ≤ 0.1; the measurement error is |v k | ≤ 0.2, k > 0; the initial state belongs to the box X 0 = [−3, 3] × [−3, 3]. Note that the measurements are linear in both x and v, and the dynamics is linear in w.
To shape the experiment, we first generated a "realized" trajectory x sim k and the respective sequence of measurements y k , k = 1, . . . , K. Then, given A k at step k, we sample N points in A k and propagate them via (11)-(12) to obtain a setX k of candidate next-step points x (i) k+1 . We call x (i) k+1 ∈X k good if it fits the measurement, i.e.
see (13). We call this account for the measurements as correction step or rejection (of bad points). Because of the linearity of the measurements, the good points are those who fall in the slab S k+1 defined by (14).
Next, we build a minimum volume set A k+1 (of the shape adopted) that contains only the good points and adopt it as a probabilistic approximation set for x k+1 . Figure 4 represents the evolution of the sets A k in the case of ellipsoidal shapes (the circle represents the ellipse E([0.6, 0.07] , 6.8I) of initial uncertainty). The stars correspond to the reference trajectory x sim k .
For the same realization we plot log vol(A k ) , see Fig. 5 (left). We remark that this quantity is of the same order of magnitude (or better) than that for the zonotopes in Alamo et al. (2008), but of course the computations are way simpler, and this advantage will get stronger as the dimension of x grows.
The filtering goal was to estimate x k (1), the first component of x k . In Fig. 5 (right), we plot the x(1)-span of our ellipses, representing the estimation error, (bold dashed curves), while with the solid line we plot the worst-case bounds on x(1) obtained from the reference trajectory: x sim k (1) ± 0.1 (see definition (11)-(12) keeping in mind that the bounded propagation noise |w k (1, 2)| ≤ 0.1 enters additively in the dynamics).
Notice that sometimes the inner corridor goes outside the bounds obtained by our method. This is exactly what illustrates the nature of our approach: We do neglect low-probability events.
Remark 3 (Comments on re-sampling and re-use) In the RPCF Algorithm, many samples could be rejected due to the measurement equation. Experiments seem to show that in such cases it is preferable to re-sample the samples rejected in the measurement phase, so as to guarantee that the new set A k+1 is constructed based on a sufficiently large number of samples. Another point worth noticing is that, on the contrary, having built the A k at step k, we already possess N good ≤ N points inside it, so that they can be re-used (propagated), and just the rest N − N good ones are to be newly sampled; this is called re-use.

Conclusions
We developed an efficient algorithm for constructing ε-probabilistic approximations of the image set of a noisy nonlinear mapping. The recursive application of such technique, combined with an ad-hoc rejection procedure, leads to the design of a new family of randomized prediction-correction filters, which showed themselves pretty competitive with both classical set-theoretical and Gaussian filtering methods. The derivations of the RCPF were limited to the NAS setting, but we stress that the use of PAS descriptions for the construction of prediction/correction filters is also possible in principle. However, we admit that in that case, the computations may become too cumbersome for online implementations.