# Generating Formal Safety Assurances for High-Dimensional Reachability

Albert Lin1 and Somil Bansal2 1Author is with CS at Princeton:
<EMAIL_ADDRESS>2Author is with ECE at University of Southern
California<EMAIL_ADDRESS>Project Website: https://sia-lab-
git.github.io/DeepReach˙Verification/. This research is supported in part by
the NVIDIA Academic Hardware Grant Program and the NSF REU Program at USC in
robotics and autonomous systems.

###### Abstract

Providing formal safety and performance guarantees for autonomous systems is
becoming increasingly important. Hamilton-Jacobi (HJ) reachability analysis is
a popular formal verification tool for providing these guarantees, since it
can handle general nonlinear system dynamics, bounded adversarial system
disturbances, and state and input constraints. However, it involves solving a
PDE, whose computational and memory complexity scales exponentially with
respect to the state dimensionality, making its direct use on large-scale
systems intractable. A recently proposed method called DeepReach overcomes
this challenge by leveraging a sinusoidal neural PDE solver for high-
dimensional reachability problems, whose computational requirements scale with
the complexity of the underlying reachable tube rather than the state space
dimension. Unfortunately, neural networks can make errors and thus the
computed solution may not be safe, which falls short of achieving our
overarching goal to provide formal safety assurances. In this work, we propose
a method to compute an error bound for the DeepReach solution. This error
bound can then be used for reachable tube correction, resulting in a safe
approximation of the true reachable tube. We also propose a scenario-based
optimization approach to compute a probabilistic bound on this error
correction for general nonlinear dynamical systems. We demonstrate the
efficacy of the proposed approach in obtaining probabilistically safe
reachable tubes for high-dimensional rocket-landing and multi-vehicle
collision-avoidance problems.

##  I Introduction

It is becoming increasingly important that we can design provably safe
controllers for autonomous systems. Hamilton-Jacobi (HJ) reachability analysis
provides a powerful framework to design such controllers for general nonlinear
dynamical systems [1, 2]. In reachability analysis, the safe states for the
system are characterized through the Backward Reachable Tube (BRT) of the
system. This is the set of states from which trajectories will eventually
reach some given target set despite the best control effort. Thus, if the
target set represents the set of undesirable states, the BRT represents unsafe
states for the system and should be avoided. Along with the BRT, reachability
analysis also provides a safety controller to keep the system outside the BRT.

Traditionally, the BRT computation in HJ reachability is formulated as an
optimal control problem. The BRT can then be obtained as a sub-zero level
solution of the corresponding value function. Obtaining the value function
requires solving a partial differential equation (PDE) over a state-space
grid, resulting in an exponentially scaling computation complexity with the
number of states [3]. To overcome this challenge, a variety of solutions have
been proposed that trade off between the class of dynamics they can handle,
the approximation quality of the BRT, and the required computation. These
include specialized methods for linear and affine dynamics [4, 5, 6, 7, 8, 9,
10, 11, 12], polynominal dynamics [13, 14, 15, 16], monotonic dynamics [17],
and convex dynamics [18] (see [3, 19] for a survey). Owing to the success of
deep learning, there has also been a surge of interest in approximating high-
dimensional BRTs [20, 21, 22, 23, 24] and optimal controllers [25] through
deep neural networks (DNN). Building upon this line of work, the authors in
[19] have proposed DeepReach – a toolbox that leverages recent advances in
neural implicit representations and neural PDE solvers to compute a value
function and a safety controller for high-dimensional systems. Compared to the
aforementioned methods, DeepReach can handle general nonlinear dynamics, the
presence of exogenous disturbances, as well as state and input constraints
during the BRT computation. However, even if we can now obtain an approximate
BRT for high-dimensional systems, it is limited in usefulness by its
approximate nature. In particular, the learned BRT can be overly optimistic;
thus, no formal safety guarantees can be provided on the obtained BRT and
controller.

In this work, our goal is to compute a BRT and a safety controller for high-
dimensional systems with provable safety guarantees. Specifically, we build
upon DeepReach and propose a verification method to provide safety assurances
on the DeepReach solution. The key insight of our work is to rely on the
consistency between the learned value function and the implicit safety
controller induced by this value function to compute a uniform correction
bound for the value function. Essentially, if a particular state is outside
the BRT (i.e., “safe” as per the learned value function), then starting from
this state, the corresponding safety controller should keep the system
trajectory outside the target set at all times (i.e., the state should also be
“safe” under the prescribed controller). Otherwise, we cannot ensure safety of
this state, and the value function should be corrected such that this state is
inside the BRT. Once all such states have been added to the learned BRT, we
have obtained a provably safe approximation of the true BRT.

We show that the computation of the correction bound can be posed as an
optimization problem. However, in general, it is challenging to tractably
compute this bound, since the learned value function can be highly nonlinear
and hard to optimize. We propose a scenario-based optimization method to
compute this correction. Scenario optimization is a sampling-based method to
solve semi-infinite optimization problems, and has been widely used for system
and control design [26, 27, 28, 29]. The proposed method is not restricted to
a specific class of system dynamics or value functions. Given that the value
function correction is obtained via randomization and, hence, is a random
quantity, we provide probabilistic guarantees on the safety of the recovered
BRT. However, this confidence is a design parameter and can be chosen as close
to 1 as desired (within a simulation budget).

To summarize, the key contributions of this paper are:

  * •

an error correction mechanism for DeepReach that results in a provably safe
BRT and safety controller for general dynamical systems;

  * •

a practical method to compute a probabilistic bound on this error correction
that is not restricted to a specific class of systems, resulting in a
tractable computation of probabilistically safe reachable tubes; and

  * •

a demonstration of the proposed approach for various dynamical systems,
inspired by rocket landing and multi-vehicle collision avoidance problems.

##  II Problem Setup

Consider a dynamical system with state $x\in X\subseteq\mathbb{R}^{n}$,
control $u\in\mathcal{U}$, and dynamics $\dot{x}=f(x,u)$ governing how $x$
evolves over time until a final time $T$. Let $\xi_{x,t}^{u}(\tau)$ denote the
state achieved at time $\tau\in[t,T]$ by starting at initial state $x$ and
time $t$ and applying control $u(\cdot)$ over $[t,\tau]$. Let $\mathcal{L}$
represent a target set of states that the agent wants to either reach (e.g.
goal states) or entirely avoid (e.g. obstacles).

Running example: Dubins3D Avoid. As a running example, consider a simple low-
dimensional system in the literature known as Dubins3D. It involves a car with
position $(p_{x},p_{y})$, heading $\theta$, velocity $v$, and steering control
$u_{1}\in[u_{\min},u_{\max}]$. The car’s state $x$ evolves according to

$\displaystyle\dot{p_{x}}=v\cos{\theta},\quad\dot{p_{y}}=v\sin{\theta},\quad\dot{\theta}=u_{1}$

It wants to avoid a circle with radius $R$ centered at the origin. Thus, we
define the target set $\mathcal{L}$ of this system to be:

$\displaystyle\mathcal{L}=\\{x:\sqrt{p_{x}^{2}+p_{y}^{2}}\leq R\\}$

We use Dubins3D as a benchmark because it is tractable for traditional
reachability methods to compute a solution for, and we can thus compare our
results with a ground truth. In Section V, we present high-dimensional
problems which traditional methods struggle with.

In this setting, we are interested in computing the system’s initial-time
Backward Reachable Tube, which we denote as BRT. We define BRT as the set of
all initial states in $X$ from which the agent will eventually reach
$\mathcal{L}$ within the time horizon $[0,T]$, despite best control efforts:

$\displaystyle\small{\text{BRT}=\\{x:x\in X,\forall
u(\cdot),\exists\tau\in[0,T],\xi_{x,0}^{u}(\tau)\in\mathcal{L}\\}}$ (1)

When $\mathcal{L}$ represents unsafe states for the system, staying outside of
BRT is desirable. When $\mathcal{L}$ instead represents the states that the
agent wants to reach (e.g., a goal set), BRT is defined as the set of all
initial states in $X$ from which the agent, acting optimally, can eventually
reach $\mathcal{L}$ within $[0,T]$. Thus, staying within BRT is desirable.

Our goal in this work is to compute a provable approximation of the safe set.
Specifically, we want to compute an approximation $\hat{\mathcal{S}}$ such
that $\hat{\mathcal{S}}\subseteq\text{BRT}^{C}$
($\hat{\mathcal{S}}\subseteq\text{BRT}$ when $\mathcal{L}$ represents goal
states). We are especially interested in settings where the system is high-
dimensional, with which current state-of-the-art reachability methods
struggle.

##  III Background: Hamilton-Jacobi (HJ) Reachability and DeepReach

In this work, we build upon Hamilton-Jacobi reachability analysis to compute
an approximation of the safe set. Here, we provide a quick overview of
Hamilton-Jacobi reachability analysis and a specific toolbox to compute high-
dimensional reachable sets, DeepReach.

###  III-A Hamilton-Jacobi (HJ) Reachability.

In HJ reachabilty, computing BRT is formulated as an optimal control problem.
We will explain it in the context of $\mathcal{L}$ being a set of undesirable
states. In the end, we will comment on the case where $\mathcal{L}$ is a set
of desirable states.

To compute BRT, we first define a target function $l(x)$ such that the sub-
zero level of $l(x)$ yields $\mathcal{L}$:

$\small{\mathcal{L}=\\{x:l(x)\leq 0\\}}\vspace{-0.25em}$ (2)

$l(x)$ is commonly a signed distance function to $\mathcal{L}$. For example,
we can choose $l(x)=\sqrt{p_{x}^{2}+p_{y}^{2}}-R$ for our Dubins3D Avoid
running example.

Next, we define the cost function of a state corresponding to some policy
$u(\cdot)$ to be the minimum of $l(x)$ over its trajectory:

$\small{J_{u(\cdot)}(x,t)=\min_{\tau\in[t,T]}l(\xi_{x,t}^{u}(\tau)).}\vspace{-0.25em}$
(3)

Since the system wants to avoid $\mathcal{L}$, our goal is to maximize
$J_{u(\cdot)}(x,t)$. Thus, the value function corresponding to this optimal
control problem is:

$\small{V(x,t)=\sup_{u(\cdot)}J_{u(\cdot)}(x,t).}\vspace{-0.25em}$ (4)

By defining our optimal control problem in this way, we can recover BRT using
the value function. In particular, the value function being sub-zero implies
that the target function is sub-zero somewhere along the optimal trajectory,
or in other words, that the system has reached $\mathcal{L}$. Thus, BRT is
given as the sub-zero level set of the value function at the initial time:

$\text{BRT}=\\{x:x\in X,V(x,0)\leq 0\\}$ (5)

The value function in Equation (4) can be computed using dynamic programming,
resulting in the following final value Hamilton-Jacobi-Bellman Variational
Inequality (HJB-VI):

$\displaystyle\min\Big{\\{}D_{t}V(x,t)+H(x,t),l(x)-V(x,t)\Big{\\}}=0,$ (6)

with the terminal value function $V(x,T)=l(x)$. $D_{t}$ and $\nabla$ represent
the time and spatial gradients of the value function. $H$ is the Hamiltonian
that encodes the role of dynamics and the optimal control.

$\displaystyle H(x,t)=\max_{u}$ $\displaystyle\langle\nabla
V(x,t),f(x,u)\rangle.$ (7)

The value function in Equation (4) induces the optimal safety controller:

$\displaystyle u^{*}(x,t)=\underset{u}{\arg\max}$ $\displaystyle\langle\nabla
V(x,t),f(x,u)\rangle.$ (8)

Intuitively, the safety controller aligns the system dynamics in the direction
of the value function’s gradient, thus steering the system towards higher-
value states.

We have just explained the case where $\mathcal{L}$ represents a set of
undesirable states. When the system instead wants to reach $\mathcal{L}$, an
infimum is used instead of a supremum in Equation $\eqref{eq:valuefunc}$. The
control wants to reach $\mathcal{L}$, hence there is a minimum instead of a
maximum in Equation (7) and Equation (8).

Traditionally, the value function is computed by solving the HJB-VI over a
discretized grid in the state space. Unfortunately, doing so involves
computation whose memory and time complexity scales exponentially with respect
to the system dimensionality, making these methods practically intractable for
high-dimensional systems, such as those beyond 5D. Fortunately, a recent deep
learning approach, DeepReach, has been proposed to enable Hamilton-Jacobi
reachability for high-dimensional systems.

###  III-B DeepReach Approximate Solutions.

Instead of solving the HJB-VI over a grid, DeepReach represents the value
function as a sinusoidal deep neural network (DNN) to learn a parameterized
approximation of the value function [19]. Thus, memory and complexity
requirements for training scale with the value function complexity rather than
the grid resolution. To train the DNN, DeepReach uses self-supervision on the
HJB-VI itself. Ultimately, it takes as input a state $x$ and time $t$, and it
outputs a learned value $\tilde{V}(x,t)$. $\tilde{V}(x,t)$ also induces a
corresponding policy $\tilde{\pi}(x,t)$. We refer interested readers to [19]
for further details.

The learned $\tilde{V}(x,t)$ from DeepReach can be used to obtain a BRT as in
Equation (5), but it will only be as accurate as $\tilde{V}(x,t)$.
Unfortunately, like any learning method, $\tilde{V}(x,t)$ can deviate
substantially from the true $V(x,t)$.

Our goal is to recover the biggest safe set $\hat{\mathcal{S}}$ using
$\tilde{V}(x,t)$ that we can (probabilistically) guarantee to be fully
contained within the true safe set. Although we work with DeepReach solutions
in particular for this problem setup, our proposed approach in Section IV can
verify any general $\tilde{V}(x,t)$ and $\tilde{\pi}(x,t)$, regardless of
whether DeepReach, a level-set method, or some other tool is used to obtain
them.

##  IV Approach

Our approach is to apply a minimal error correction to $\tilde{V}(x,t)$ such
that the corrected value function can be used to extract a safe set. We will
explain our approach in the context of $\mathcal{L}$ being a set of
undesirable states. In the end of each subsection, we will comment on the case
where $\mathcal{L}$ represents a set of desirable states. In Section IV-A, we
propose an error metric for correcting the value function. In Section IV-B, we
propose a practical method to compute a probabilistic bound on this error
metric through scenario optimization.

###  IV-A Error Metric

We propose a new error metric $\delta_{\tilde{V},\tilde{\pi}}$ for the value
function correction. It is defined as the maximum learned value of an
empirically unsafe initial state under the induced policy $\tilde{\pi}$:

$\displaystyle\small{\delta_{\tilde{V},\tilde{\pi}}\coloneqq\max_{x\in
X}\\{\tilde{V}(x,0):J_{\tilde{\pi}}(x,0)\leq 0\\}},$ (9)

where $J_{\tilde{\pi}}(x,0)$ is the cost function associated with the
trajectory obtained by using the policy $\tilde{\pi}(x,t)$ from an initial
state $x$ and initial time $t=0$ until $T$ (see Equation (3)). Here on, we use
$\delta$ as a shorthand for $\delta_{\tilde{V},\tilde{\pi}}$ for brevity
purposes.

Intuitively, $\delta$ finds the tightest level of the learned value function
that separates the states that are safe under the induced policy
$\tilde{\pi}(x,t)$ from the ones that are not. Thus, any initial state within
the super-$\delta$ level set of $\tilde{V}(x,0)$ is guaranteed to be safe
under the (possibly sub-optimal) policy $\tilde{\pi}(x,t)$. Lemma 1 formalizes
this claim.

###### Lemma 1

If we compute $\tilde{\mathcal{S}}$ as:

$\displaystyle\small{\tilde{\mathcal{S}}=\\{x\in X:\tilde{V}(x,0)>\delta\\}}$
(10)

then $\tilde{\mathcal{S}}\subseteq\text{BRT}^{C}$.

The proof for Lemma 1 is in Appendix-A of the extended version of this article
[30]. Intuitively, Lemma 1 states that if we could exactly compute $\delta$,
then we can recover a $\tilde{\mathcal{S}}$ that is guaranteed to be a subset
of the true safe set. Furthermore, it is clear that, by definition, $\delta$
is the smallest (uniform) value adjustment required on $\tilde{V}(x,0)$ to
guarantee its safety. Thus, $\tilde{\mathcal{S}}$ is the largest safe set we
can recover from the learned value function by such a uniform error correction
procedure.

Unfortunately, computing $\delta$ is a challenging optimization problem, since
both the value function and the cost function in the metric definition (9) are
typically non-convex functions of $x$. In the next section, we propose a
scenario-based optimization approach to compute a high-confidence bound on
$\delta$ and thereby generate a high-confidence safe set.

###### Remark 1

We assume there exists an unsafe state for the metric definition (9).
Otherwise, define $\delta_{\tilde{V},\tilde{\pi}}$ trivially as $-\infty$.

So far, we have discussed an error correction metric for the case where
$\mathcal{L}$ represents a set of undesirable states. When $\mathcal{L}$
represents a set of desirable states, we use a minimum instead of a maximum
and flip the cost inequality in the metric definition (9). To recover the safe
set, we extract the sub-$\delta$ level set of the value function in Equation
(10).

###  IV-B Scenario Optimization Verification Method

We will compute a high-confidence bound on $\delta$ by utilizing a random
sampling procedure referred to as scenario optimization in the systems and
control design literature [29]. The Scenario Optimization Verification Method
to compute an approximation $\hat{\delta}$ is summarized in Algorithm 1.

0: $X$, $N$, $M$, $\tilde{V}(x,0)$, $J_{\tilde{\pi}}(x,0)$

1: $\delta_{0}\leftarrow-\infty$

2: for $i=0,1,\ldots,M-1$ do

3: $\mathcal{D}_{i}\leftarrow$ Sample $N$ states IID from $\\{x:x\in
X,\tilde{V}(x,0)>\delta_{i}\\}$

4: if $\exists x\in\mathcal{D}_{i}:J_{\tilde{\pi}}(x,0)\leq 0$ then

5:
$\delta_{i}\leftarrow\max_{x\in\mathcal{D}_{i}}\\{\tilde{V}(x,0):J_{\tilde{\pi}}(x,0)\leq
0\\}$

6: else

7: break

8: end if

9: end for

10: return $\hat{\delta}:=\delta_{i}$

Algorithm 1 Scenario Optimization Verification

At a high-level, Algorithm 1 computes a converging sequence of $\delta_{i}$
that approximates $\delta$ via random sampling until we no longer find any
safety violations or reach a maximum number of iterations, $M$. Specifically,
at each iteration $i$, we randomly sample $N$ initial states within the
super-$\delta_{i}$ level set of $\tilde{V}$ (Line 3) using rejection sampling
and compute the costs of associated trajectories $J_{\tilde{\pi}}(x,0)$ under
$\tilde{\pi}(x,t)$ as the controller (Equation (3)). We next compute the
maximum learned value among the states that violate safety and use that as the
new estimate of $\delta$ (Line 5)), essentially discarding a level region of
$\tilde{V}(x,0)$ that is empirically unsafe under $\tilde{\pi}$. Thus, with
each iteration, we obtain a tighter approximation of $\delta$ and terminate
when no more safety violations are found or a maximum number of iterations is
achieved.

While Algorithm 1 is simple to understand and execute, scenario optimization
provides a formal guarantee for bounding $\delta$ with $\hat{\delta}$ when the
algorithm terminates prior to reaching the maximum number of iterations. Thus,
we can use the approximated $\hat{\delta}$ instead of $\delta$ to compute an
approximate safe set $\hat{\mathcal{S}}$ similar to the one defined in Lemma
1:

$\displaystyle\small{\hat{\mathcal{S}}=\\{x:x\in
X,\tilde{V}(x,0)>\hat{\delta}\\}}$ (11)

Crucially, we can make a formal probabilistic safety guarantee for
$\hat{\mathcal{S}}$. This is summarized in Theorem 1 and proven in Appendix-C
of the extended version of this article [30].

######  Theorem 1 (Scenario Optimization Verification Theorem)

Select a violation parameter $\epsilon\in(0,1)$ and a confidence parameter
$\beta\in(0,1)$. Pick $N$ such that

$\displaystyle\small{N\geq\frac{2}{\epsilon}\left(\ln{\frac{1}{\beta}}+1\right)}$
(12)

Suppose Algorithm 1 converges, then with probability at least $1-\beta$, the
recovered safe set $\hat{\mathcal{S}}$ in Equation (11) satisfies:

$\begin{aligned}
\small{\underset{x\in\hat{\mathcal{S}}}{\mathbb{P}}(V(x,0)\leq
0)\leq\epsilon}\end{aligned}\vspace{0.50em}\hbox to0.0pt{$\qquad\Box$\hss}$
(13)

In practice, we find that the algorithm converges within 3-4 iterations.
Intuitively, the higher the quality of the value function approximation, the
faster the convergence of $\delta_{i}$.

Disregarding the confidence parameter $\beta$ for a moment, Theorem 1 states
that the volume of the unsafe states within $\hat{\mathcal{S}}$ is smaller
than or equal to the prescribed $\epsilon$ value, as long as we sample enough
states to satisfy Inequality (12) during the computation of $\hat{\delta}$. As
$\epsilon$ approaches $0$, the number of safety violations in the recovered
safe set also approaches $0$. In turn, the simulation effort grows unbounded
since $N$ is inversely proportional to $\epsilon$.

To interpret the confidence parameter $\beta$, note that $\hat{\mathcal{S}}$
is a random variable that depends on a randomly sampled set of initial states.
It may be the case that we just happen to draw a poorly representative sample,
in which case the $\epsilon$ bound does not hold. $\beta$ controls the
probability of this adverse event happening, which regards the correctness of
the final guarantee in Equation (13). Fortunately, $N$ only grows
logarithmically with $\frac{1}{\beta}$, so $\beta$ can be chosen to be an
extremely small value such as $10^{-12}$. $1-\beta$ should then be so close to
$1$ that it does not have any practical importance.

When $\mathcal{L}$ represents a set of desirable states, we initialize
$\delta_{0}$ to be $\infty$ instead of $-\infty$, flip the inequalities, and
take a minimum instead of a maximum in Algorithm 1. We flip the value
inequalities in Equation (11) and Equation (13).

###### Remark 2

Note that when Algorithm 1 does not converge, scenario optimization can still
be used to bound $\delta$. However, the safety guarantees are more involved.
We defer a detailed investigation of these guarantees to future work.

Running example: Dubins3D Avoid. We now demonstrate the Scenario Optimization
Verification method on the Dubins3D Avoid running example introduced in
Section II. We choose system parameters
$v=0.6m/s,u_{\min}=-1.1rad/s,u_{\max}=1.1rad/s,R=0.25m$, and use DeepReach to
learn the system’s value function. The used neural network parameters and
architecture are the same as in the DeepReach paper [19], i.e., a DNN with 3
hidden layers and 512 neurons in each layer.

For this first example, we purposefully sabotage the training of DeepReach by
permitting it to train for only 8,000 epochs, far less than $100$K epochs used
in the DeepReach paper. We do this to illustrate our method’s utility even
with low-quality learned solutions. We apply our method to a well-trained
DeepReach solution in Section V-A.

The overall training took less than 30 minutes on an NVIDIA GeForce RTX
3090Ti. We execute our verification approach choosing
$\beta=10^{-16},\epsilon=10^{-3}$, resulting in an $N=75683$ as per Theorem 1.
Thus, we will be $1-10^{-16}$ confident that at least $1-10^{-3}$ of our
recovered set will be provably safe. Algorithm 1 converges in 4 iterations in
this case, resulting in a $\hat{\delta}=0.3728$. Slices of the trained BRT
(yellow), the recovered BRT (blue), and the ground truth BRT boundary (black)
are shown in Figure 1. The ground truth is computed by a state-of-the-art PDE
solver, Level Set Toolbox (LST) [31], which computes the value function over a
discrete grid of size $101\times 101\times 101$.

Figure 1: Dubins3D Avoid: Slices of the trained, recovered, and ground truth
border BRTs for three values of $\theta$ from a sabotaged DeepReach solution.
The recovered BRT completely encompasses the ground truth BRT and provides a
probabilistically safe approximation of the safe set.

As evident from the figure, the recovered BRT slices completely encompass the
true unsafe states (the blue region encompasses the states within the black
boundary). Hence, the complement of the recovered BRT is safe, as we expect
based on Theorem 1. In the rightmost slice, observe that the trained BRT is
smaller than the ground truth BRT and is thus unsafe as it is. Our
verification approach correctly expands the trained BRT past the ground truth
boundary, ensuring safety. Validation by sampling 1M states within the
recovered safe set reveals a violation rate of $3\times 10^{-6}<<\epsilon$.

Note that the trained solution’s safety violations are concentrated near
$\theta=\pi$. It is likely these errors which are responsible for the large
expansion of the recovered BRT even for slices that are almost accurate (left
and middle slices in Figure 1). This is because our method expands the BRT
uniformly. However, too much expansion is undesirable, since we want to
recover as much of the safe set as possible.

One way to overcome this challenge is to apply refined error correction to
$\tilde{V}(x,t)$ by independently verifying separate regions of the states-
space (i.e., separately consider each region as our $X$). This allows us to
recover sets with the same safety guarantee but expanded differently depending
on their own region’s error. However, this method’s efficacy depends heavily
on the choice of sub-regions, which are hard to prespecify, especially for
high-dimensional systems. Instead, we can resort to a data-driven approach.

Specifically, we trained a simple multilayer perceptron (MLP) on 1M randomly
sampled initial states and their empirical safety violation costs (the
training took an hour on a standard GPU). Even though the predictor may be
inaccurate (after all, predicting the safe set is the original challenging
problem), we hypothesize that it will have learned something useful enough
about the general distribution of errors throughout the state space. We
divided the output of the MLP in 10 bins, each corresponding to a different
range of empirical cost. Thus, we expect the last bin to correspond to the
states with maximal safety violations (in this case, that will correspond to
states with $\theta=\pi$). We next ran Algorithm 1 independently for each bin
to compute an error correction. The corresponding results are shown in Figure
2. As evident, we get a much tighter approximation of the safe set near the
slices that have small errors, allowing us to recover a much bigger safe set
overall. In this case, the recovered safe set volume increased by 35%.

Figure 2: Dubins3D Avoid: Slices of the trained, recovered, and ground truth
border BRTs for three values of $\theta$ from a sabotaged DeepReach solution
when recovery is refined by MLP binning. The recovered BRT completely
encompasses the ground truth BRT and provides a probabilistically safe
approximation of the safe set that is larger than when binning is not used.

As discussed earlier, the performance improvement of this approach depends
heavily on the selection of subregions. We defer this to future work, and for
the rest of this paper, we only focus on uniform value function correction.

##  V Case Studies

We will evaluate our approach on various reachability problems. First, we show
results for the avoid version and the reach version of the low-dimensional
Dubins3D system for which we have a ground truth BRT available for comparison.
Next, we recover safe sets for a 9D multivehicle collision avoidance problem
and for a 6D rocket landing problem, with which traditional methods struggle.
Like our running example, we use $\beta=10^{-16},\epsilon=10^{-3},N=75683$.

###  V-A Dubins3D Avoid

In Figure 3, we show the results of applying our method to a well-trained
DeepReach solution for the same Dubins3D Avoid system as in the running
example in Section IV. The solution is trained to consider the periodicity of
$\theta$ and for 100K, instead of 8K, epochs. The proposed algorithm results
in a very small $\hat{\delta}=-0.0016$. Interestingly, the fact that
$\hat{\delta}<0$ indicates that the learned value function is conservative.
Our method shrinks the trained BRT to recover a larger safe set. Validation by
sampling 1M states within the recovered safe set reveals a violation rate of
$3.2\times 10^{-5}<<\epsilon$.

The fact that a much smaller error correction is found for the trained
solution shown in Figure 3 than the one in Figure 1 indicates it is of a much
higher quality. This illustrates how the error correction metric $\delta$ can
also be used to evaluate the quality of different approximate value function
solutions. For example, we can evaluate the relative performance of different
DNN hyperparameters, which is especially helpful when the ground truth BRT is
not available for comparison.

Figure 3: Dubins3D Avoid: Slices of the trained, recovered, and ground truth
border BRTs for three values of $\theta$ from a well-trained DeepReach
solution. The trained BRT completely encompasses the recovered BRT, so we plot
the recovered BRT border on top for visualization purposes.

###  V-B Dubins3D Reach

In Figure 4, we show the results of applying our method to a DeepReach
solution trained with the same scheme as for the Dubins3D Avoid solution in
Section V-A, but for the reach version. That is, $\mathcal{L}$ is now a set of
desirable states, so a safe approximate BRT should be fully contained by the
true BRT. The trained BRT is just barely crossing the ground truth boundary,
and the proposed method cuts it down to just behind the boundary, ensuring
system safety.

Figure 4: Dubins3D Reach: Slices of the trained, recovered, and ground truth
border BRTs for three values of $\theta$ from a well-trained DeepReach
solution. The recovered BRT is correctly a subset of the ground truth BRT.

###  V-C Multivehicle Collision Avoidance

We now consider a 9D collision avoidance system involving 3 independent
Dubins3D cars. The $i$th car has position $(p_{xi},p_{yi})$, heading
$\theta_{i}$, velocity $v$, and steering control
$u_{i}\in[u_{\min},u_{\max}]$. The dynamics of vehicle $i$ are given as:

$\displaystyle\dot{p}_{xi}=v\cos{\theta_{i}},\quad\dot{p}_{yi}=v\sin{\theta_{i}},\quad\dot{\theta_{i}}=u_{i}$

The (undesirable) target set is given by the states where any of the vehicle
pairs is in collision:

$\mathcal{L}=\\{x:\min\\{d(Q_{1},Q_{2}),d(Q_{1},Q_{3}),d(Q_{2},Q_{3})\\}\leq
R\\}$

where $d(Q_{i},Q_{j})$ is the distance between cars $i,j$. We choose
$v=0.6,u_{\min}=-1.1,u_{\max}=1.1,R=0.25$ for our case study. The results of
applying our method to a DeepReach solution are shown in Figure 5.

Figure 5: Multivehicle Collision Avoidance: Slices of the trained and
recovered BRTs for three values of $\theta_{1}$ from a DeepReach solution.

This high-dimensional example is typically difficult to compute with
traditional methods, yet here we demonstrate successful recovery of a safe set
with formal safety guarantees. While the recovery preserves much of the safe
set that is learned, a significant amount is pruned off. The percent reduction
in safe set size (calculated from 1M samples) is much larger in this example
(26%) than for the Dubins3D solutions in Section V-A (0%) and Section V-B
(10%), suggesting that the learned value function is more inaccurate for
higher dimensional systems.

To validate safety of the recovered BRT, in Figure 6, we plot the minimum
pairwise distance between vehicles along trajectories spawning from states
sampled within the trained and recovered safe sets. Note that the trained safe
set contains states which result in a collision, as shown by the mass to the
left of the dotted line representing the collision distance. After
verification, the recovered safe set shows no such mass. Validation by
sampling 1M samples in the recovered safe set reveals a violation rate of
$5\times 10^{-6}<<\epsilon$.

Figure 6: Multivehicle Collision Avoidance: Histogram of the minimum pairwise
spare distance between vehicles along trajectories which spawn from states
sampled within the trained and recovered safe sets. The collision distance is
indicated by the dotted line.

###  V-D Rocket Landing

Consider a 6D rocket landing system with position $(p_{x},p_{y})$, heading
$\theta$, velocity $(v_{x},v_{y})$, angular velocity $\omega$, and torque
controls $\tau_{1},\tau_{2}\in[-250,250]$.

The dynamics are:

$\displaystyle\dot{p_{x}}=v_{x},~{}\dot{p_{y}}=v_{y},~{}\dot{\theta}=\omega,~{}\dot{\omega}=0.3\tau_{1},$
$\displaystyle\dot{v_{x}}=\tau_{1}\cos{\theta}-\tau_{2}\sin{\theta},~{}\dot{v_{y}}=\tau_{1}\sin{\theta}+\tau_{2}\cos{\theta}-g,$

where $g=9.81$ is acceleration due to gravity. The target set is the set of
states where the rocket reaches a rectangular landing zone of side length 20
centered at the origin:

$\displaystyle\mathcal{L}=\\{x:|p_{x}|<20.0,p_{y}<20.0\\}$

In Figure 7, we show the results of applying our method to a DeepReach
solution. We also plot the trajectories emanating from 20 randomly sampled
initial states within the recovered safe set, demonstrating that they do
indeed reach the target set (the green region) by following the induced
policy. Validation by sampling 1M states within the recovered safe set reveals
a violation rate of $5\times 10^{-6}<<\epsilon$.

Figure 7: Rocket Landing: Slices of the trained and recovered BRTs for three
values of $v_{y}$ from a DeepReach solution. The system trajectories safely
reach the target set (green) from the initial states within the recovered BRT.

##  VI Discussion and Future Work

In this work, we present an approach to compute an error bound for DeepReach
solutions to recover a provably safe approximation of the true reachable tube.
We also propose a practical method to compute a probabilistic bound on this
error correction that is not restricted to a specific class of systems. This
allows us to utilize the power of learning-based reachability methods to
provide probabilistic safety assurances for high-dimensional dynamical
systems. We apply our method to obtain probabilistically safe reachable tubes
for high-dimensional rocket-landing and multi-vehicle collision-avoidance
problems which traditional methods struggle with.

In the future, we will explore a more refined approach to error correction (as
opposed to a uniform correction), as discussed briefly in the end of Section
IV. Other directions include considering worst-case disturbances in the system
dynamics and using the verification approaches proposed here for a targeted
refinement of DeepReach solutions.

## References

  * [1] J. Lygeros, “On reachability and minimum cost optimal control,” _Automatica_ , vol. 40, no. 6, pp. 917–927, 2004.
  * [2] I. Mitchell, A. Bayen, and C. J. Tomlin, “A time-dependent Hamilton-Jacobi formulation of reachable sets for continuous dynamic games,” _IEEE Transactions on Automatic Control (TAC)_ , vol. 50, no. 7, pp. 947–957, 2005.
  * [3] S. Bansal, M. Chen, S. Herbert, and C. J. Tomlin, “Hamilton-Jacobi Reachability: A brief overview and recent advances,” in _IEEE Conference on Decision and Control (CDC)_ , 2017.
  * [4] M. R. Greenstreet and I. Mitchell, “Integrating projections,” in _Hybrid Systems: Computation and Control_ , T. A. Henzinger and S. Sastry, Eds. Berlin, Heidelberg: Springer Berlin Heidelberg, 1998, pp. 159–174.
  * [5] G. Frehse, C. Le Guernic, A. Donzé, S. Cotton, R. Ray, O. Lebeltel, R. Ripado, A. Girard, T. Dang, and O. Maler, “SpaceEx: Scalable verification of hybrid systems,” in _International Conference Computer Aided Verification_ , 2011.
  * [6] A. B. Kurzhanski and P. Varaiya, “Ellipsoidal techniques for reachability analysis: internal approximation,” _Systems & Control Letters_, 2000.
  * [7] A. Kurzhanski and P. Varaiya, “On ellipsoidal techniques for reachability analysis. part ii: Internal approximations box-valued constraints,” _Optimization Methods and Software_ , vol. 17, pp. 207–237, 01 2002.
  * [8] J. N. Maidens, S. Kaynama, I. M. Mitchell, M. M. Oishi, and G. A. Dumont, “Lagrangian methods for approximating the viability kernel in high-dimensional systems,” _Automatica_ , 2013.
  * [9] A. Girard, “Reachability of uncertain linear systems using zonotopes,” in _International Workshop on Hybrid Systems: Computation and Control_ , 2005, pp. 291–305.
  * [10] M. Althoff, O. Stursberg, and M. Buss, “Computing reachable sets of hybrid systems using a combination of zonotopes and polytopes,” _Nonlinear analysis: hybrid systems_ , vol. 4, no. 2, pp. 233–249, 2010.
  * [11] S. Bak, H.-D. Tran, and T. T. Johnson, “Numerical verification of affine systems with up to a billion dimensions,” in _International Conference on Hybrid Systems: Computation and Control_ , 2019, pp. 23–32.
  * [12] P. Nilsson and N. Ozay, “Synthesis of separable controlled invariant sets for modular local control design,” in _American Control Conference_ , 2016, pp. 5656–5663.
  * [13] A. Majumdar, R. Vasudevan, M. M. Tobenkin, and R. Tedrake, “Convex optimization of nonlinear feedback controllers via occupation measures,” _The International Journal of Robotics Research_ , vol. 33, no. 9, pp. 1209–1230, 2014. [Online]. Available: https://doi.org/10.1177/0278364914528059
  * [14] A. Majumdar and R. Tedrake, “Funnel libraries for real-time robust feedback motion planning,” _The International Journal of Robotics Research_ , vol. 36, no. 8, pp. 947–982, 2017.
  * [15] T. Dreossi, T. Dang, and C. Piazza, “Parallelotope bundles for polynomial reachability,” in _International Conference on Hybrid Systems: Computation and Control_ , 2016.
  * [16] D. Henrion and M. Korda, “Convex computation of the region of attraction of polynomial control systems,” _IEEE Transactions on Automatic Control_ , vol. 59, no. 2, pp. 297–312, 2014.
  * [17] S. Coogan and M. Arcak, “Efficient finite abstraction of mixed monotone systems,” in _Proceedings of the 18th International Conference on Hybrid Systems: Computation and Control_ , 2015, pp. 58–67.
  * [18] Y. T. Chow, J. Darbon, S. Osher, and W. Yin, “Algorithm for overcoming the curse of dimensionality for time-dependent non-convex hamilton–jacobi equations arising from optimal control and differential games problems,” _Journal of Scientific Computing_ , vol. 73, no. 2-3, pp. 617–643, 2017.
  * [19] S. Bansal and C. J. Tomlin, “DeepReach: A deep learning approach to high-dimensional reachability,” in _IEEE International Conference on Robotics and Automation (ICRA)_ , 2021.
  * [20] V. Rubies-Royo, D. Fridovich-Keil, S. Herbert, and C. J. Tomlin, “A classification-based approach for approximate reachability,” in _International Conference on Robotics and Automation_. IEEE, 2019, pp. 7697–7704.
  * [21] J. F. Fisac, N. F. Lugovoy, V. Rubies-Royo, S. Ghosh, and C. J. Tomlin, “Bridging Hamilton-Jacobi Safety Analysis and Reinforcement Learning,” _International Conference on Robotics and Automation_ , 2019.
  * [22] B. Djeridane and J. Lygeros, “Neural approximation of pde solutions: An application to reachability computations,” in _Conference on Decision and Control_ , 2006, pp. 3034–3039.
  * [23] K. Niarchos and J. Lygeros, “A neural approximation to continuous time reachability computations,” in _Conference on Decision and Control_ , 2006, pp. 6313–6318.
  * [24] J. Darbon, G. P. Langlois, and T. Meng, “Overcoming the curse of dimensionality for some hamilton–jacobi partial differential equations via neural network architectures,” _Research in the Mathematical Sciences_ , vol. 7, no. 3, pp. 1–50, 2020.
  * [25] D. Onken, L. Nurbekyan, X. Li, S. W. Fung, S. Osher, and L. Ruthotto, “A neural network approach for high-dimensional optimal control applied to multiagent path finding,” _IEEE Transactions on Control Systems Technology_ , 2022.
  * [26] S. Ghosh, S. Bansal, A. Sangiovanni-Vincentelli, S. A. Seshia, and C. Tomlin, “A new simulation metric to determine safe environments and controllers for systems with unknown dynamics,” in _Proceedings of the 22nd ACM International Conference on Hybrid Systems: Computation and Control_ , 2019, pp. 185–196.
  * [27] M. C. Campi and S. Garatti, “A sampling-and-discarding approach to chance-constrained optimization: feasibility and optimality,” _Journal of Optimization Theory and Applications_ , 2011.
  * [28] G. C. Calafiore and M. C. Campi, “The scenario approach to robust control design,” _IEEE Transactions on Automatic Control_ , 2006.
  * [29] M. C. Campi, S. Garatti, and M. Prandini, “The scenario approach for systems and control design,” _Annual Reviews in Control_ , 2009.
  * [30] A. Lin and S. Bansal, “Generating formal safety assurances for high-dimensional reachability,” _arXiv preprint arXiv:2209.12336_ , 2022\.
  * [31] I. Mitchell, _A Toolbox of Level Set Methods_ , 2009, http://people.cs.ubc.ca/~mitchell/ToolboxLS/index.html.

## Appendix

###  -A Proof of Lemma 1

###### Proof:

Assume for a contradiction that we compute $\tilde{\mathcal{S}}$ as in
Equation (10) and that

$\displaystyle\small{\tilde{\mathcal{S}}\not\subseteq\text{BRT}^{C}}$ (14)

which implies

$\displaystyle\small{\exists
x^{{}^{\prime}}:x^{{}^{\prime}}\in\tilde{\mathcal{S}},x^{{}^{\prime}}\not\in\text{BRT}^{C}}$
(15)

$x^{{}^{\prime}}\in\tilde{\mathcal{S}}$ implies that

$\displaystyle\small{\tilde{V}(x^{{}^{\prime}},0)>\delta}$ (16)

by Equation (10).

On the other hand, $x^{{}^{\prime}}\not\in\text{BRT}^{C}$ implies that
$V(x^{{}^{\prime}},0)\leq 0$. We know that
$\forall(x,t),J_{\tilde{\pi}}(x,t)\leq V(x,t)$ from Equation (4). Thus,
$J_{\tilde{\pi}}(x^{{}^{\prime}},0)\leq V(x^{{}^{\prime}},0)\leq 0$. Then

$\displaystyle\small{\tilde{V}(x^{{}^{\prime}},0)\leq\delta}$ (17)

by Equation (9).

We have found the contradiction that $\tilde{V}(x^{{}^{\prime}},0)>\delta$ and
$\tilde{V}(x^{{}^{\prime}},0)\leq\delta$. Thus, our assumption must be false,
and it must be true that if we compute $\tilde{\mathcal{S}}$ as in Equation
(10), then $\tilde{\mathcal{S}}\subseteq\text{BRT}^{C}$.

To prove the equivalent of Lemma 1 but in the case when $\mathcal{L}$
represents a set of desirable states, we use BRT instead of $\text{BRT}^{C}$
and flip the value and cost inequalities. ∎

###  -B Lemma 2

###### Lemma 2

Select a violation parameter $\epsilon\in(0,1)$ and a confidence parameter
$\beta\in(0,1)$. Pick $N$ such that

$\displaystyle\small{N\geq\frac{2}{\epsilon}\left(\ln{\frac{1}{\beta}}+1\right)}$
(18)

Then, with probability at least $1-\beta$, the solution $\hat{\delta}$ to
Algorithm 1 executed with $N$ satisfies the following condition.

$\displaystyle\small{\underset{\\{x:x\in
X,\tilde{V}(x,0)>\hat{\delta}\\}}{\mathbb{P}}\left(\left(\tilde{V}(x,0):J_{\tilde{\pi}}(x,0)\leq
0\right)>\hat{\delta}\right)\leq\epsilon}$ (19)

###### Proof:

Lemma 2 follows from Theorem 1 in Scenario Optimization [29]. To use the
theorem, we need to prove the following conditions:

  1. 1.

computing $\delta$ can be converted into a standard Scenario Optimization
problem

  2. 2.

Algorithm 1 obtains $\hat{\delta}$ by sampling i.i.d from the same space that
$\delta$ is optimized over, $X$, but with a probability distribution that is
uniform over $\\{x:x\in X,\tilde{V}(x,0)>\hat{\delta}\\}$

$\delta$ can be formalized as the following optimization problem:

$\displaystyle\min{g}$ $\displaystyle\text{s.t. }\forall x\in
X,\left(\tilde{V}(x,0):J_{\tilde{\pi}}(x,0)\leq 0\right)\leq g$

This is a semi-infinite optimization problem where the constraints are linear,
and thus convex, in the optimization variable $g$ for any given $x$. Thus,
Lemma 2 follows from Theorem 1 in Scenario Optimization [29] by replacing $c$
by $1$, $\gamma$ by $g$, $\Delta$ by the state space of interest $X$, and $f$
by $\left(\tilde{V}(x,0):J_{\tilde{\pi}}(x,0)\leq 0\right)-g$. To use Theorem
1, we also require that i.i.d samples are chosen according to the uniform
distribution over $\\{x:x\in X,\tilde{V}(x,t)>\hat{\delta}\\}$. This can be
proven by observing that, in Algorithm 1, the last iteration of the while loop
samples $x$ randomly and independently from $\\{x:x\in
X,\tilde{V}(x,t)>\hat{\delta}\\}$ where $\hat{\delta}$ is exactly the final
metric returned. Importantly, the algorithm does not update $\hat{\delta}$ in
the final loop (since the loop subsequently breaks).

To prove the equivalent of Lemma 2 but in the case when $\mathcal{L}$
represents a set of desirable states, we flip the value and cost inequalities.
We also take a maximum instead of a minimum and flip the inequality with $g$
when we formalize $\delta$ as an optimization problem. ∎

###  -C Proof of Theorem 1

###### Proof:

We rewrite the LHS of Inequality (19) of Lemma 2 in Appendix-B by substituting
in $\hat{\mathcal{S}}$ by Equation (11):

$\displaystyle\underset{\\{x:x\in
X,\tilde{V}(x,0)>\hat{\delta}\\}}{\mathbb{P}}$
$\displaystyle\left(\left(\tilde{V}(x,0):J_{\tilde{\pi}}(x,0)\leq
0\right)>\hat{\delta}\right)$ $\displaystyle=\underset{\\{x:x\in
X,\tilde{V}(x,0)>\hat{\delta}\\}}{\mathbb{P}}$
$\displaystyle\left(\left(\tilde{V}(x,0)>\hat{\delta}\right)\cap\left(J_{\tilde{\pi}}(x,0)\leq
0\right)\right)$ $\displaystyle=\underset{\\{x:x\in
X,\tilde{V}(x,0)>\hat{\delta}\\}}{\mathbb{P}}$
$\displaystyle\left(J_{\tilde{\pi}}(x,0)\leq 0\right)$
$\displaystyle=\underset{x\in\hat{\mathcal{S}}}{\mathbb{P}}$
$\displaystyle\left(J_{\tilde{\pi}}(x,0)\leq 0\right)$

where the second-to-last line follows because the probability is over the set
of states that are already given to have a learned value greater than
$\hat{\delta}$.

We know that $\forall(x,t),J_{\tilde{\pi}}(x,t)\leq V(x,t)$ from Equation (4).
Thus, it follows that:

$\displaystyle\small{\underset{x\in\hat{\mathcal{S}}}{\mathbb{P}}\left(V(x,0)\leq
0\right)\leq\underset{x\in\hat{\mathcal{S}}}{\mathbb{P}}\left(J_{\tilde{\pi}}(x,0)\leq
0\right)\leq\epsilon}$ (20)

To prove the equivalent of Theorem 1 but in the case when $\mathcal{L}$
represents a set of desirable states, we flip the value and cost inequalities.
∎


# DeepLens: Towards a Visual Data Management System

Sanjay Krishnan  Adam Dziedzic  Aaron J. Elmore

University of Chicago  {skr, ady<EMAIL_ADDRESS>
###### Abstract

Advances in deep learning have greatly widened the scope of automatic computer
vision algorithms and enable users to ask questions directly about the content
in images and video. This paper explores the necessary steps towards a future
Visual Data Management System (VDMS), where the predictions of such deep
learning models are stored, managed, queried, and indexed. We propose a query
and data model that disentangles the neural network models used, the query
workload, and the data source semantics from the query processing layer. Our
system, DeepLens, is based on dataflow query processing systems and this
research prototype presents initial experiments to elicit important open
research questions in visual analytics systems. One of our main conclusions is
that any future “declarative” VDMS will have to revisit query optimization and
automated physical design from a unified perspective of performance and
accuracy tradeoffs. Physical design and query optimization choices can not
only change performance by orders of magnitude, they can potentially affect
the accuracy of results.

##  1 Introduction

Recent advances in deep learning have enabled new forms of analysis on images
and videos [19]. A typical analytics task is to find all images in a corpus
that contain a particular object (e.g., detecting a red car in CCTV footage)
using neural network-based object detection models. Answering such queries at
scale is a formidable computer systems challenge with significant recent
interest from both industry and academia [22, 21, 15, 3, 13, 4, 6, 27, 10,
11].

Recent work on visual analytics is clearly inspired by ideas from relational
database systems [21, 13, 14]: there is a high-level query language and an
optimized execution layer that processes the queries. _But, are these new
systems truly declarative like an RDMBS?_ We argue that existing systems
neglect some of the critical challenges in the space; giving an illusion of a
declarative interface for simple query processing tasks. This paper explores
the execution tradeoffs in visual analytics and illustrates a complex
interplay storage, latency, and accuracy considerations. Abstracting these
tradeoffs away from an end-user will be the primary technical challenge in a
future Visual Data Mangement System (VDMS).

An obvious first omission from existing proposals is a principled study of
input data formats. Today’s systems interface with video data _as independent
frame-ordered sequences of raw images_. While this format is convenient for
implementation, it may not be the optimal physical layout for all use cases.
Video data is often stored and transferred in an encoded formats, such as
H.264, which can significantly reduce the storage size of high-resolution
video (sometimes by multiple orders of magnitude). On the other hand, due to
the sequential nature of most video encoding schemes (where decoding requires
scanning preceding frames), encoding precludes random access to particular
frame numbers or times in the video–affecting possible execution paths to
answer a given query. Furthermore, most encoding schemes are lossy, which
leads to downstream accuracy considerations after decoding. One could also
imagine hybrid physical layouts that segment the video into short encoded
clips and then bucket them by time. The complexity of this space means that
any VDMS has to allow for a multitude of encoding and physical layout schemes,
and decouple these layouts from the programming model.

Beyond just physical layout, answering complex queries in such systems is
challenging. Consider a variant of the task before: given two videos find a
particular object that appears in both videos. For example, we might be
interested in determining if the same car appeared in two different CCTV
feeds. To answer this query, one has to first find potential target objects in
both videos and then match them against each other. Even systems that support
_join queries_ across videos [21], have difficulty answering questions when
the join predicate considers pixel data. To be able to answer such a query
efficiently, there are implementation decisions about indexing (can a
multidimensional index be used for faster matching), optimization (if we do
index, which video to scan and what type of an index to use), and compression
(can the matching be performed on low-dimensional features instead of raw
frames). Indexing and query optimization are significant missing components in
current proposals, which are necessary to scale such systems to more complex
analytics tasks.

To go beyond these limitations, this paper explores the missing pieces towards
a full-featured VDMS. For the programming abstraction, we start with a analog
to the relational model. All visual corpora can be described as an unordered
collection of subimages (called patches) and an associated key-value
dictionary storing information about them (e.g., neural net classifications,
their provenance). These patches are abstract data types (stored in raw pixels
or pre-compressed to features) that can be generated from object detection
models that crop identified objects, classical visual segmentation methods
that crop based on color, or even whole images. The query model, which allows
for Selects and Joins over patch collections, has set semantics as in SQL and
makes no structural assumptions about the data source, e.g., any timing data
in video or geometric data from camera correspondences are stored as
additional attributes. In other words, every operator defines closed-algebra:
collection of patches in and collection of patches out. Optimizations for
specific workloads are introduced through physical design: data layout, and
both single-attribute and multi-dimensional indexing on patch collections.

Figure 1: DeepLens has a dataflow-like architecture for processing visual
analytics queries. All analysis is cast as relational queries on relations of
image patches. Intermediate results can be materialized and indexed.

This logical layer of indirection disentangles the process that generates the
patches (decoding, neural network inference, segmentation) with the downstream
query processing of the collection. The query operators are agnostic to
whether those patches represent the output of a neural network, black-box
geometric computer vision pipeline, or even the output from crowdsourcing.

We built a research prototype system, DeepLens, that implements this model. To
illustrate the benefits of physical and logical separation, we design a
benchmark of six analytics tasks and vary what is indexed, how the operators
are implemented, the underlying algorithms, and the properties of the
underlying hardware. A substantial number of queries of interest require
matching predictions across multiple videos or matching outputs across
multiple prediction pipelines, and our central insights is that _indexing_ and
_lineage_ are critical missing pieces from the current discourse on visual
analytics.

### Key Insights

Encoding:  Compressing video data with standard codecs can save up to 50x in
storage, but requires additional computation and scan costs to process.
Segmenting video into clips and then compressing those clips provides a
flexible way to optimize this tradeoff between storage and computation.

Indexing:  We discovered that index usage is crucial and improves query
performance by up-to 612x. However, particularly for the geometric indices,
the size and dimensionality of the data plays an essential role in how
beneficial an index will be.

Lineage:  We found that many tasks require relating processed results back to
the base data. Maintaining and indexing tuple-level lineage led to a 60x
improvement in one benchmark query.

CPU vs. GPU:  Cost-models that accurately balance CPU and GPU utilization for
query optimization will be a significant challenge. Some queries benefited by
nearly two orders of magnitude, while others are slower when processing tasks
are offloaded to a GPU.

Managing Uncertainty:  Unlike in relational databases, queries in a VDMS are
approximate by nature. Traditional optimizations, like filter pushdown, may
additionally create different accuracy profiles.

One of our main conclusions is that any future “declarative” VDMS will have to
_revisit query optimization and automated physical design from a unified
perspective of performance and accuracy tradeoffs._

##  2 Background and API

As a direct consequence of the rapid progress in computer vision, visual
analytics workloads are becoming increasingly complex with requirements. One
can easily make brittle assumptions about ordering of images or the structure
of any associated metadata. Or, one might use ad-hoc data structures that will
not scale past memory limits. This echoes the challenges of the pre-relational
days of database research before the idea of data independence where the query
processing, the data, and the storage format are all intertwined. We want a
unified model that covers much of visual analytics that can allow us to
leverage logical-physical separation principles that the database community
has pioneered.

###  2.1 The “Narrow Waist” of Visual Analytics

DeepLens is based on a logical model where data are unordered collections of
image patches (featurized sub-images) with their associated derivation
information and metadata. The user writes programs against such a logical
model and a query optimizer selects an appropriate execution path. The patches
and metadata are typed to allow the system to validate logical equivalence.
This architecture abstracts lower level concerns about the physical layout,
while simultaneously ensuring that the system can reason about the application
logic. Physical concerns include video encoding, temporal segmentation, and
cross-video indexing. Application-level concerns include data lineage, whether
a component (e.g., an object detector) can be replaced with another one, and
latency-accuracy tradeoffs. We envision that such a narrow waist design is
crucial as we scale visual analytics to harder problems, such as in robotics
where it is increasingly common to integrate predictions from multiple neural
networks [9].

###  2.2 Iterator API for Patches

After surveying a number of use cases from robotics to traffic camera
analysis, we arrived at a simple query model, where visual analytics queries
are relational queries over collections of subimages (called patches) and
associated metadata about how those subimages were generated. The query
processing engine is agnostic to how those patches are generated. They can be
whole images, smaller tiled subimages, or even subimages extracted by an
object detection neural networks. In pseudocode, a Patch object contains a
pointer to the image that generated it, the data contained in the patch, and a
key-value metadata dictionary:

⬇

Patch(ImgRef, Data, MetaData)

We are purposefully vague about the structure of Data and MetaData. We only
assume that Data is an n-dimensional dense vector (can represent pixel content
or image features) and MetaData has a key-value dictionary of attributes of
this data. Inspired by dataflow systems, operators in the system implement
iterators over tuples of Patch objects:

⬇

Operator(Iterator<Tuple<Patch>> in,

Iterator<Tuple<Patch>> out)

Lineage is maintained as every operator is required to update the ImgRef
attribute to retain a lineage chain back to the original image.

####  2.2.1 Example 1.

Let us first consider an example studied in existing work, such as [21, 22,
13]. Consider a CCTV feed of a parking lot that collects and stores video. We
want to evaluate the parking lot’s utilization so we want to count the number
of cars in each frame of the video. We first run all of the frames through the
Single-Shot Detector (SSD) object detection network [20], which returns
bounding boxes and labels for common objects detected in the image. Each of
these bounding boxes defines a patch:

⬇

SSDPatch(Frame, Bbox,

{’label’: L, ’frameno’: F })

Over these SSD patch objects, the query of interest is well expressed in
relational algebra over the metadata dictionary (a filter over labels and an
aggregation over frame numbers).

####  2.2.2 Example 2

Now, let us consider a marginally more difficult example. Suppose we are given
two videos from two different cameras, we want to find all cars that appear in
both videos. The first step is exactly the same as the previous example, we
break the frames up into SSDPatch objects. Over the two sets of SSDPatch, we
need to compare all of the bounding boxes and return those that are
sufficiently similar in terms of image content. What is different about this
query is that it leverages both the metadata and the pixel data in the
bounding boxes.

Interestingly, there are a number of unresolved questions about how to
actually process this query efficiently. Naively, one could compare all pairs
of bounding boxes and then return those of sufficient image similarity. Most
image matching algorithms use lower dimensional features to match, so another
option is to pre-compute the relevant features and build a multidimensional
index over one of the sets of SSDPatch objects, e.g., a KD-Tree over a set of
color histograms.

###  2.3 Related Work

The new era of “deep” multimedia systems are not actually multimedia database
systems in the traditional sense [5]. In current systems, the queries issued
to pixel data are very restricted: either binary neural network predicates
[15, 27, 2, 11] or highly structured [21, 22]. Our work is inspired by prior
work on multimedia databases, which have long acknowledged the importance of
indexing and query processing strategies [26, 5]. DeepLens revisits the idea
of a multimedia database in the era of deep learning, where the content–both
structured and unstructured–is populated by the outputs of a neural network
inference pipeline. Recent work can be best summarized as filter optimization
[15, 27, 2, 11]; how to evaluate a neural network predicate as quickly as
possible while satisfying an accuracy constraint. The Optasia system [21] goes
beyond simple predicates and does consider joins; however, as noted earlier
the joins do not efficiently handle predicates over the image data (not just
the metadata). We project into the future that the community will soon move
past filters and also consider joins and more complex query operators.
Addressing the more complex visual analytics queries will require leveraging
ideas from classical work on indexing strategies [5]. However, problems in
this new setting are higher dimensional and have to manage uncertainty.
Another important design decision was to natively manage image transformation
and provenance code in our system. We chose this architecture since recent
work manipulates neural network structure, e.g., by cascading models to
improve accuracy or performance [15, 2, 10]. Another recent paper that we
would like to highlight is a data management system for Augmented Reality [8],
and many of the ideas relating to multidimensional indexing will be very
relevant. This paper evaluates all of these design trade-offs and proposes an
initial architecture for a modern VDMS. One of our main conclusions is that
any successful VDMS system will have to feature a sophisticated and robust
query optimizer and some level of automated physical design. We envision that
new results in learning-based query optimization [12, 16] and ideas inspired
by the new automated tools for physical design such as [25, 23] will be
crucial.

##  3 Storage Layer

Our architecture is designed as a middle ground between the pure stateless
dataflow approach and the RDBMS approach. DeepLens provides dataflow operators
to shape, transform, and manipulate image corpora at scale, but also supports
materialization and indexing of intermediate results. Figure 1 illustrates the
basic architecture with three main layers: (1) Persistent Storage, (2) ETL,
and (3) Query Processing. We first describe the storage layer in DeepLens.

###  3.1 Raw Data

In our initial prototype, we only consider videos at rest and not streaming
videos. Our loading API is as follows:

⬇

Load(filename, filter=True)

Videos are loaded into the system returning an iterator that returns a patch
collection where each patch is a full video frame: 3-D dense arrays
representing width and height of images, as well as corresponding RGB
channels. The loader can take a filter as an optional argument and it only
returns those frames that satisfy the filter condition. The loader abstracts
the encoding scheme of the underlying video from the user. DeepLens supports
the following storage formats:

Frame File:  In the most basic format, we treat each frame of a video as a
single record. This frame can either be stored in its raw pixel values or in
some featurized format. The _Frame File_ is implemented with BerkeleyDB
111Implemented with bsddb3 (Python binding for BerkeleyDB). The image and
feature data is serialized in a binary format before insertion. By default,
they are stored in a sorted file by frame number for videos and wall clock
time. This is because many queries of interest examine particular time
segments of videos. The sorted file allows for quick retrieval of temporal
predicates. The advantage of the Frame File is a temporal filter push down;
the disadvantage is that it can require significantly more storage for video.
For example, the H.264 codec often reduces the storage footprint of video by
20+ times. Storing video as a sequence of image frames precludes the use of
sequential encoding techniques.

Encoded File:  Therefore, we also support storing video in common encoded
formats (e.g., OGG, MPEG4). The tradeoff is that encoding precludes pushing
down temporal predicates since many encoding formats require a sequential
decoding procedure. Furthermore, the cost of decoding the video must also be
factored in.

Segmented File:  As a hybrid between the Frame File and the Encoded File, we
have the Segmented File. This storage format segments the video into short
clips and stores the encoded clips in BerkeleyDB. We can benefit from coarse-
grained temporal filter push down, while having some benefits of encoding.

###  3.2 Derived Data

Any of the intermediate results in DeepLens can be materialized as well. We
also support the construction of indexes on the materialized data. The
challenge is that every data type requires a specialized index structure. Over
string valued or discrete metadata, the index choices are straight-forward. We
support hash tables and B+ Trees over any key (both of which are implemented
with BerkeleyDB).

For the multidimensional data (patch segmentation parameters or image
features), indexing is a little more nuanced and dependent on the workload. As
a concrete example, let us consider patches that are parametrized by “bounding
boxes” (x1,y1,x2,y2). Even though this data is multidimensional, we might
effectively be interested in single dimensional queries. For example, find
bounding boxes left of a certain point in the image. To support such queries,
we found that it was far more efficient to use a B+ Tree rather than an R-Tree
multidimensional index. On the other hand, if we are interested in containment
and intersection queries, then we need a true multidimensional index. We
provide an interface to a disk-based R-Tree implemented with
libspatialindex222https://libspatialindex.github.io/.

However, we found that existing R-Tree implementations are optimized for
geospatial problems in 2D. They could not be efficiently modified for higher
dimensional data, for instance, for image matching queries, where we compare
features of two images and threshold the similarity. For this class of queries
and data, a data structure called a Ball-Tree was the most effective at
answering Euclidean threshold queries in high-dimensional spaces [17].

In summary, DeepLens supports a large number of single-dimensional (Hash, B+
Tree, Sorted Files) and multi-dimensional indexes (R-Tree and Ball-Tree). The
physical design problems are interesting as even the same attributes might be
indexed in different ways depending on the workload.

### Future Work: Storage Advisor

In this paper, we prototyped all of the different storage options, and
manually selected the appropriate ones per use case. In the future, we would
like to have a storage advisor that can analyze a workload or an SLO and
return an optimized storage scheme. We see an analog with the work that has
been done in the database community on materialization and storage format
tuning [1].

##  4 Visual ETL

Once loaded, the image dataflow is passed into patch generators that turns the
image into a set of patches. Associated with each patch is a key-value
dictionary storing metadata about how the patch was generated and any property
of the patch. The patches are then fed into a composition of transformers that
featurize, compress, or otherwise store detected properties of the patches
into the dictionary. Over tuples of transformed patches, users can build a
directed computation graph with operators (e.g., select, join). Tuple-level
lineage is automatically maintained by the system allowing any downstream
patch to be associated with its base data. Any of the intermediate results can
be materialized.

###  4.1 Visual ETL

In a sense, image data is unstructured data. Semantics from image data have to
be first extracted with computer vision algorithms before structured queries
can be executed. The ETL layer defines the generation of Patch objects and
their manipulation.

Patch Generators:  We provide a library of Patch Generators. These generators
take as input an iterator over raw images and return an iterator over Patch
objects. Our experiments consider three instantiations of these generators:
object detection (ones that segment objects out of the image return
segmentation masks), optical character recognition (one that identifies text
in an image and returns a segmentation as well as the recognized strings) and
whole-image patches (returns the whole image).

Transformers:  The main content of a patch is still raw pixel data, which is
often not very useful on its own. To be able to compare or manipulate the
patches, we need a featurized representation. This leads to the next module of
the ETL layer that defines Transformers. A transformer takes as input an
iterator over Patch objects and returns an iterator over transformed Patch
objects. In our experiments, we consider two transformers: color histogram
features for image matching and a depth prediction neural network which
predicts the 3D geometry of a patch.

Materialize:  DeepLens allows any stage of this pipeline to be materialized
and persisted to disk or memory.

###  4.2 Validation

The entire API in DeepLens is typed, which allows us to validate pipelines.
The computer vision community is rapidly changing, and we foresee new visual
processing primitives will constantly be added to any VDMS. Beyond the
standard int, float, string types, our type system maintains the resolution
and dimensions of each patch (almost all neural networks used today require
fixed input resolutions). We also include the domains of any discrete metadata
created when available. For example, object detection networks have a closed-
world of labels that they could predict. We include this information in the
type system. Any downstream operator (e.g., filter) that consumes those labels
can be validated to see if that label is plausibly produced by the pipeline.

### Future Work: Pipeline Synthesis

Long term, we would like a system that could declaratively synthesize a
pipeline given a library of generators and transformers and latency/accuracy
constraints. Concretely, suppose we had a library of general purpose pre-
trained object detection models and some special case programmed models (e.g.,
a car detector). We envision a system that scores each model with a
precision/recall profile for a desired dataset, and can choose the model that
is most appropriate for a query. This can only be done if there is an
appropriate type system that understands what labels can be predicted and when
detectors are interchangeable.

##  5 Query Processing

Our query processing engine is designed like a dataflow query processing
system [7]. In our initial prototype, we implement Select and Join operators.
The design of the Join operators are the most interesting so we highlight it
here.

Nested Loop Join:  If no indexes are available, the most generic operator is a
nested loop join operator. This join operator can execute arbitrary
$\theta$-joins on the data. It compares all pairs of patches from two
collections and returns those that satisfy a predicate.

Index Joins:  If a multi-dimensional or single dimensional index is available,
we can use that index to enable equality joins, range joins, or similarity
joins.

On-The-Fly Index Similarity Join:  We found that for image matching queries
where one of the relations was relatively small, the index could be
constructed on-the-fly. We load the smaller relation into an in-memory Ball-
Tree. Then, probe using the other collection of patches.

###  5.1 Lineage

Many visual analytics tasks of interest relate processed results back to the
base data. For example, we might process a single set of images in two
different ways, e.g., segmenting the image with an object detector and using a
depth prediction model to determine relative distance between pixels. To
relate these results, we have to run a _backtracing_ query (select all raw
images that contributed to a patch). This is similar to queries studied in
recent lineage systems [24].

DeepLens natively tracks tuple-level lineage. Every Patch object maintains a
descriptor how it was generated from either a raw image or another patch. Its
relationship to the base data is maintained as a sequence of pointers. This
information is stored as attributes in the metadata key-value dictionary so
indexes and queries can be natively supported on them.

### Future Work: Visual Query Optimizer

We believe that the key to a true VDMS will be a query optimizer than can
manage and abstract all of the complexity of the physical operations.
Classical database operators are mostly I/O bound, but on the other hand, many
VDMS operators are compute-bound. The caveat is that device placement
(CPU/GPU) matters for these operators. Cost-models that accurately balance CPU
and GPU utilization for query optimization will be a significant challenge. In
our experiments, some queries benefited by nearly two orders of magnitude by
GPU optimizations, while others are actually slower when the processing tasks
are offloaded to a GPU. We hope to explore new strategies for joint multi-
query optimization and optimal physical design. Doing so might well require
significant machine learning to navigate a noisy and analytically complex cost
model. Image queries are always approximate, and managing uncertainty and
quantifying the accuracy effects of a certain plan will be a serious
challenge.

##  6 Benchmarks

We briefly discuss how we designed the benchmark datasets and corresponding
workload.

###  6.1 Datasets

From publicly available sources, we constructed datasets to simulate
envisioned use cases for DeepLens. An important consideration was images of
varying format and size. We decided that it was important to evaluate a
diversity of tasks and the robustness of the injestion pipeline to different
formats, sizes, and content.

PC. This dataset is designed to simulate a dataset of images found on a
personal computer. It consists of 779 photographs, screenshots, and document
scans.

TrafficCam. This dataset consists of 24 mins and 30 secs of high-definition
(1080p) traffic camera video (35280 frames).

Football. This dataset consists of 15 low-definition (720p) videos of American
football clips of the same team ranging from 30 secs to 1 mins (15244 total
images).

###  6.2 Queries

We propose a benchmark workload of 6 queries on these datasets to evaluate
DeepLens. These queries are inspired by problems considered in prior work and
consist of tasks that involve querying pixel data, the results of neural
network predictions, relating results back to base data, and combinations of
these tasks.

q1. _Find all near-duplicates in the PC dataset_. This query is inspired by
classical multimedia information retrieval problems such as reverse image
search (find the closest image to a query image).

q2. _Count all of the frames with at least one vehicle present in the
TrafficCam dataset_. This query is inspired by recent work that uses neural
networks to analyze traffic and movement patterns [15]. This is a simple query
that takes the output of a neural network that identifies objects in the frame
and simply queries the output.

q3. _Track one player’s trajectory in every play in the Football dataset_.
Given segmentation output that identifies a player in frame and OCR output
that identifies a number if one is visible, we have to relate that sequence of
bounding boxes back to the original image.

q4. _Count all distinct pedestrians in the TrafficCam dataset_. This query is
a variant of q2. The distinct qualifier makes this query significantly more
challenging as it requires deduplicating candidate pedestrians detected in the
video.

q5. _Lookup the presence of a string in the PC dataset_. Apply OCR to all of
the images, and store a collection of strings discovered. We run a query to
identify the first image with a target string.

q6. _Find all tuples of pedestrians (p1,p2) where p1 is behind p2 in the
TrafficCam dataset_. This query is inspired from applications in robotics and
navigation, where one has to estimate how far a given object is located from a
camera. This problem, called depth prediction, has recently been a subject of
research interest in computer vision [18]. We leverage the published
code333https://github.com/iro-cp/FCRN-DepthPrediction and the pre-trained
parameters to annotate all detected pedestrians in the TrafficCam dataset with
depth predictions and find such pairs.

##  7 Experiments

Recent work on video and image analytics systems has largely focused on
optimizing neural network inference. For queries that are simply filters or
aggregates, neural network inference clearly dominates the processing cost. In
particular, we use our benchmark to understand the execution profiles of more
complex queries where the subsequent processing time before or after neural
network inference might dominate.

Figure 2: Encoding a video with a sequential codec can reduce storage costs by
over 50x without loss of accuracy.

###  7.1 Data Encoding

In the first experiment, we illustrate some of the interesting challenges with
data encoding. We compare different encoding formats for the video in Q2 in
terms of storage cost as well as the overall accuracy of the pipeline (Figure
2). We estimated the accuracy by manual annotation of the video. We compare
three levels of lossy encoding: High, Medium, Low.

Figure 3: Hybrid storage formats can support coarse-grained filter push down
as well as take advantage of sequential compression.

The RAW encoding (where every frame is an image) rests at about 107 GB on
disk. In contrast, a state-of-the-art H.264 encoding444Implemented with
https://github.com/cisco/openh264 only takes up 2.5 GB on disk. These
encodings are all lossy, but there is negligible impact on downstream accuracy
for the "high quality" encodings. For larger compression ratios, we do see a
degradation in downstream performance.

Of course, there is a latency consideration as well. We evaluate a subset of
the input formats in terms of their end-to-end latency (including decoding) in
Figure 3. We add a temporal filter to Q2 (selecting a small subset of frames
to analyze). The H.264 encoding cannot support a true filter push down as the
codec algorithm is sequential and requires a scan of the preceding video until
that point. On the other hand, the JPEG and RAW formats can trivially support
the push down optimization. A hybrid scheme that breaks the video into short
clips and then encodes them, can allow for coarse-grained filter push down. We
manually tuned this granularity for best performance.

These encoding experiments illustrate the complexity of visual analytics
physical layout problems. VDMSes must support a variety of these formats and
structures to be effective. We hope this is also an opportunity to re-visit
automatic storage tuning and declarative interfaces for storage configuration.

Figure 4: DeepLens significantly speeds up “query time” by using indexes. The
queries that match multidimensional features can be sped up by up-to 600x.

###  7.2 The Power of Indexes

Next we study indexing and materializing intermediate results. In these
experiments we separate “Query time” (processing patches) from “ETL time”
(generating patches). Many ETL processes can be reused for multiple queries,
and those costs can be amortized.

Figure 4 plots the query times with and without indexing. Our baseline is the
same query processing engine with no indexes. We compare this baseline to a
hand-tuned version where we manually select the best physical design for a
query.

The queries that benefit the most from the indexes are ones that require image
matching, and are up-to 612x faster for q4 and 59x faster for q1. Note that
these matching costs are incurred after neural network inference, so
techniques such as [15, 3, 13] would have no effect on this step. Queries that
have to relate results back to the base images also see improvements since
they can leverage the lineage information and do not need to rescan the base
data. q3 requires a backtracing query to match the bounding boxes and the OCR
output in pixels on the original image and has a 41x improvement. Similarly,
q6 runs 2.5x faster. q5 is illustrative of a query whose predicate does not
benefit from any of the available indexes.

Figure 5: We evaluate the pipeline runtime, including ETL and "on-the-fly"
index creation, for an optimized DeepLens (DL) vs. the baseline (BL). On many
queries it is beneficial to materialize intermediate results and build indexes
to speed up future performance. Indexing has a relatively small overhead given
the compute-intensive nature of the queries.  Figure 6: Building
multidimensional indexes can be very costly and initial experiments indicate
that construction time scales poorly with the increase of data size.

###  7.3 Overhead For Indexing

Existing systems in this domain do not make a distinction between Query time
and ETL time. We believe this separation is justified because indexes can be
reused by other queries and their construction costs amortize. Regardless,
several of the queries execute faster even if the indexes are built “on-the-
fly” (Figure 5). For example, q1 executes nearly 5 times faster than the
baseline and q4 executes 3.5 times faster than the baseline. The index
significantly reduces the number of image matching operations that dominate
the runtime, and thus, the cost of building and persisting the index is
offset.

The indexes with largest overhead to build are the spatial indexes. Figure 6
plots the construction time of the multidimensional and single dimensional
indexes supported in DeepLens as a function of the number of tuples indexed.
The R-Tree is nearly 20x slower to construct than a B+ Tree. However, there
are some mitigating factors that are interesting to consider in future work.
Since visual analytics is approximate by nature, perhaps exact
multidimensional indexing is unnecessary. For some workloads may suffice to
apply single dimensional indices and merge results with independence
assumptions. For others, locality sensitive hashing or similar approximations
may suffice.

Figure 7: We evaluate the execution time of a Ball-Tree join as function of
the size of the indexed relation in the high-dimensional and low-dimensional
case. As the data structure is increasingly filled the execution time grows
non-linearly.

###  7.4 Subtleties in Query Optimization

For queries that join multiple image collections, neural network inference
time does not necessarily dominate. Optimizing these spatial joins turned out
to be more challenging than we thought.

####  7.4.1 Nonlinear join costs

To avoid a large library of hand-written rules, the natural question to ask is
whether it is possible to build a cost-based query optimizer for this system.
Figure 7 illustrates the difficulties with the execution time of a Ball-Tree
join as function of the size of the indexed relation in the high-dimensional
and low-dimensional case. As the data structure is increasingly filled the
execution time grows non-linearly. The non-linearity is also data-dependent
and is more extreme in higher dimensional data. Accurately modeling the
relationship between input relation size and operator cost is crucical for
cost-based query optimization. Non-linearities are also known to affect
standard join ordering heuristics [16].

####  7.4.2 Balancing CPU vs. GPU

Processing these complex queries requires a nuanced strategy in balancing CPU
vs. GPU usage. For the ETL phase of these queries, which is dominated by
neural network inference time, using the GPU is almost universally better.
Figure 8 (left) plots the processing time on each of the six benchmark queries
for a vanilla CPU implementation (CPU), a vectorized execution (AVX), and a
GPU implementation (GPU). Just by changing the underlying execution
architecture there were up-to 12x changes in execution time. On the other
hand, the results were more mixed for the query time. For the two image
matching queries, q1 and q4, we implemented an all pairs matching comparison
with vectorization and on the GPU. For the larger query (q4) there is a
significant performance benefit from using the GPU (34% faster). For the
smaller query (q1), the overhead of using the GPU outweighs the costs. Any
cost-based optimizer has to weigh overheads before selecting a plan.

Figure 8: The execution architecture has a considerable impact on both ETL and
Query time: vanilla CPU implementation (CPU), a vectorized execution (AVX),
and GPU implementation (GPU). GPUs are significantly more efficient for the
neural network dominated ETL time, but have more mixed results for the
subsequent query processing on two image matching queries (q1, q4).

####  7.4.3 Accuracy implications of different query plans

Unlike in relational databases, queries in a VDMS are approximate by nature.
Cascades of approximate operators can have correlations in their errors that
are hard to reason about. These errors can be affected by query optimizer
choices. Let us take q4 as an example. To process this query, the system first
identifies patches using an object detector, then filters those patches to
those that are people, and then matches all pairs of patches to deduplicate.
In another approach, the system would first identify patches using an object
detector, then match all pairs of patches (regardless if they are “people”) to
deduplicate, and finally filter on those pairs that have at least one person
label. The second approach goes against typical query optimization principles
of filter pushdown–but we see that it is actually a more accurate strategy.

Execution method for q4 | Recall | Precision | Runtime
---|---|---|---
Patch, Filter, Match | 0.73 | 0.97 | 34.56
Patch, Match, Filter | 0.82 | 0.98 | 62.11
Table 1: Accuracy vs. runtime for different execution methods of query q4.

##  8 Conclusion and Future Work

The growing maturity of neural networks for image classification and
segmentation problems has made visual analytics an attractive and inter-
discplinary field of study. We explore the intersection of these neural
network models and data management by designing a query processing engine
called DeepLens. Our long-term goal is develop a scalable and performant
visual data management system. While recent work has been very focused on
speeding up neural network inference [15, 27, 2, 11], we find that as we start
processing increasingly complex queries, the neural network inference time no
longer dominates and our experiments illustrate other bottlenecks in the
system. Addressing these complex queries requires a unified data and query
model and the logical-physical separation seen in a traditional RDBMS.

Perhaps, one of the most compelling reasons to study the design of VDMS is
that it brings together many hard challenges in database research. In
particular, we believe the query optimization challenge here is significant.
We hope to explore new strategies for joint multi-query optimization and
optimal physical design. Doing so might well require significant machine
learning to navigate a noisy and analytically complex cost model. Image
queries are always approximate, and managing uncertainty and quantifying the
accuracy effects of a certain plan will be a serious challenge. We plan on
considering improved techniques for image compression and approximate query
processing. The natural numerical representation of images provides us with
more structure than we typically have in normal database workloads.

Acknowledgements: We thank Eugene Wu and Michael Maire for their helpful
feedback on early drafts. We also acknowledge NVIDIA for a generous equipment
grant. This work is supported in part by the Center For Unstoppable Computing
(CERES) at the University of Chicago.

## References

  * [1] S. Agrawal, E. Chu, and V. Narasayya.  Automatic physical design tuning: workload as a sequence.  In Proceedings of the 2006 ACM SIGMOD international conference on Management of data, pages 683–694. ACM, 2006.
  * [2] M. R. Anderson, M. Cafarella, G. Ros, and T. F. Wenisch.  Physical representation-based predicate optimization for a visual analytics database.  arXiv preprint arXiv:1806.04226, 2018.
  * [3] M. R. Anderson, M. Cafarella, T. F. Wenisch, and G. Ros.  Predicate optimization for a visual analytics database.
  * [4] S. Chetlur, C. Woolley, P. Vandermersch, J. Cohen, J. Tran, B. Catanzaro, and E. Shelhamer.  cudnn: Efficient primitives for deep learning.  arXiv preprint arXiv:1410.0759, 2014.
  * [5] C. Faloutsos.  Searching multimedia databases by content, volume 3.  Springer Science & Business Media, 2012.
  * [6] Z. Feng, J. Wang, J. Harkes, P. Pillai, and M. Satyanarayanan.  Eva: An efficient system for exploratory video analysis.
  * [7] G. Graefe.  Volcano – an extensible and parallel query evaluation system.  IEEE Transactions on Knowledge and Data Engineering, 6(1):120–135, 1994.
  * [8] B. Haynes, A. Mazumdar, A. Alaghi, M. Balazinska, L. Ceze, and A. Cheung.  Lightdb: a dbms for virtual reality video.  Proceedings of the VLDB Endowment, 11(10):1192–1205, 2018.
  * [9] R. Hodson.  How robots are grasping the art of gripping?, 2018.
  * [10] A. H. Jiang, D. L.-K. Wong, C. Canel, L. Tang, I. Misra, M. Kaminsky, M. A. Kozuch, P. Pillai, D. G. Andersen, and G. R. Ganger.  Mainstream: Dynamic stem-sharing for multi-tenant video processing.  In 2018 USENIX Annual Technical Conference. USENIX Association, 2018.
  * [11] J. Jiang, G. Ananthanarayanan, P. Bodik, S. Sen, and I. Stoica.  Chameleon: scalable adaptation of video analytics.  In Proceedings of the 2018 Conference of the ACM Special Interest Group on Data Communication, pages 253–266. ACM, 2018.
  * [12] T. Kaftan, M. Balazinska, A. Cheung, and J. Gehrke.  Cuttlefish: A lightweight primitive for adaptive query processing.  arXiv preprint arXiv:1802.09180, 2018.
  * [13] D. Kang, P. Bailis, and M. Zaharia.  BlazeIt: Fast Exploratory Video Queries using Neural Networks.  arXiv preprint arXiv:1805.01046, 2018.
  * [14] D. Kang, P. Bailis, and M. Zaharia.  Challenges and Opportunities in DNN-Based Video Analytics: A Demonstration of the BlazeIt Video Query Engine.  CIDR, 2019.
  * [15] D. Kang, J. Emmons, F. Abuzaid, P. Bailis, and M. Zaharia.  NoScope: optimizing neural network queries over video at scale.  Proceedings of the VLDB Endowment, 10(11):1586–1597, 2017.
  * [16] S. Krishnan, Z. Yang, K. Goldberg, J. Hellerstein, and I. Stoica.  Learning to optimize join queries with deep reinforcement learning.  arXiv preprint arXiv:1808.03196, 2018.
  * [17] N. Kumar, L. Zhang, and S. Nayar.  What is a good nearest neighbors algorithm for finding similar patches in images?  In European conference on computer vision, pages 364–378. Springer, 2008.
  * [18] I. Laina, C. Rupprecht, V. Belagiannis, F. Tombari, and N. Navab.  Deeper depth prediction with fully convolutional residual networks.  2016 Fourth International Conference on 3D Vision (3DV), pages 239–248, 2016.
  * [19] Y. LeCun, Y. Bengio, and G. Hinton.  Deep learning.  nature, 521(7553):436, 2015.
  * [20] W. Liu, D. Anguelov, D. Erhan, C. Szegedy, S. Reed, C.-Y. Fu, and A. C. Berg.  Ssd: Single shot multibox detector.  In European conference on computer vision, pages 21–37. Springer, 2016.
  * [21] Y. Lu, A. Chowdhery, and S. Kandula.  Optasia: A relational platform for efficient large-scale video analytics.  In Proceedings of the Seventh ACM Symposium on Cloud Computing, pages 57–70. ACM, 2016.
  * [22] Y. Lu, A. Chowdhery, S. Kandula, and S. Chaudhuri.  Accelerating machine learning inference with probabilistic predicates.  In Proceedings of the 2018 International Conference on Management of Data, pages 1493–1508. ACM, 2018.
  * [23] A. Pavlo, G. Angulo, J. Arulraj, H. Lin, J. Lin, L. Ma, P. Menon, T. C. Mowry, M. Perron, I. Quah, et al.  Self-driving database management systems.
  * [24] F. Psallidas and E. Wu.  Smoke: Fine-grained lineage at interactive speed.  Proceedings of the VLDB Endowment, 11(6):719–732, 2018.
  * [25] A. Sharma, F. M. Schuhknecht, and J. Dittrich.  The case for automatic database administration using deep reinforcement learning.  arXiv preprint arXiv:1801.05643, 2018.
  * [26] A. Yoshitaka and T. Ichikawa.  A survey on content-based retrieval for multimedia databases.  IEEE Transactions on Knowledge and Data Engineering, 11(1):81–93, 1999.
  * [27] C. Zhang, Q. Cao, H. Jiang, W. Zhang, J. Li, and J. Yao.  Ffs-va: A fast filtering system for large-scale video analytics.  In Proceedings of the 47th International Conference on Parallel Processing, page 85. ACM, 2018.