Scaling the Convex Barrier with Active Sets
Authors: Alessandro De Palma
Harkirat Behl
Rudy R Bunel
Philip Torr
M. Pawan Kumar
Venue: ICLR 2021 Poster (2021)
Keywords: Neural Network Verification
Neural Network Bounding
Optimisation for Deep Learning

Abstract: Tight and efficient neural network bounding is of critical importance for the scaling of neural network verification systems. A number of efficient specialised dual solvers for neural network bounds have been presented recently, but they are often too loose to verify more challenging properties. This lack of tightness is linked to the weakness of the employed relaxation, which is usually a linear program of size linear in the number of neurons. While a tighter linear relaxation for piecewise linear activations exists, it comes at the cost of exponentially many constraints and thus currently lacks an efficient customised solver. We alleviate this deficiency via a novel dual algorithm that realises the full potential of the new relaxation by operating on a small active set of dual variables. Our method recovers the strengths of the new relaxation in the dual space: tightness and a linear separation oracle. At the same time, it shares the benefits of previous dual approaches for weaker relaxations: massive parallelism, GPU implementation, low cost per iteration and valid bounds at any time. As a consequence, we obtain better bounds than off-the-shelf solvers in only a fraction of their running time and recover the speed-accuracy trade-offs of looser dual solvers if the computational budget is small. We demonstrate that this results in significant formal verification speed-ups.

Published as a conference paper at ICLR 2021
SCALING THE CONVEX BARRIER WITH ACTIVE SETS
Alessandro De Palma,, Harkirat Singh Behl", Rudy Bunel, Philip H.S. Torr,
M. Pawan Kumar
University of Oxford {adepalma, harkirat, phst, pawann)robots..om ac. bunel. rudyyggadiaagaaa com com
ABSTRACT
Tight and efficient neural network bounding is of critical importance for the seal- ing of neural network verification systems. A number of efficient specialised dual solvers for neural network bounds have been presented recently, but they are often too loose to verify more challenging properties. This lack of tightness is linked to the weakness of the employed relaxation, which is usually linear program of size linear in the number of neurons. While tighter linear relaxation for piecewise linear activations exists, it comes at the cost of exponentially many constraints and thus currently lacks efficient customised solver. We alleviate this deficiency via a novel dual algorithm that realises the full potential of the new relaxation by op- erating small active set of dual variables. Our method recovers the strengths of the new relaxation in the dual space: tightness and a linear separation oracle. At the same time, it shares the benefits of previous dual approaches for weaker relax- ations: massive parallelism, GPU implementation, low cost per iteration and valid bounds at any time. As consequence, we obtain better bounds than aff-the-shelf solvers in only fraction of their running time and recover the speed-accuraee trade-offs of looser dual solvers if the computational budget is small. We demon- strate that this results in significant formal verification speed-ups.
INTRODUCTION
Verification requires formally proving disproving that given property of neural network holds over all inputs in a specified domain. We consider properties in their canonical form (Bunel et al. 2018), which requires us to either: (i) prove that input results in negative output (property is true); or (ii) identify a counter-example (property is false). The search for counter-examples is typically performed by efficient methods such as random sampling of the input domain (Webb al.) 2019), or projected gradient descent (Carlini & Wagner 2017). In contrast, establishing the veracity of property requires solving a suitable convex relaxation to obtain lower bound on the minimum output. If the lower bound is positive, the given property is true. If the bound is negative and no counter-example is found, either: (i) we make no conclusions regarding the property (incomplete verification); or (ii) we further refine the counter-example search and lower bound computation within a branch-and-boun framework until we reach concrete conclusion (complete verification).
The main bottleneck of branch and bound is the computation of the lower bound for each node of the cnumeration tree via convex optimization. While earlier works relied on off-the-shelf solvers (Ehlers 2017; Bunel et al. 20188) it was quickly established that such an approach does not scale-up elegantly with the size of the neural network. This has motivated researchers to design specialized dual solvers (Dvijotham et al. idence that verification can be realised in practice. However, the convex relaxation considered in the dual solvers is itself very weak (Ehlers) 2017), hitting what now commonly referred to as the "convex barrier" (Salman et al., 2019). In practice, this implies that either several properties remain undecided in incomplete verification, or take several hours to be verified exactly.
2019 Bunel et al. 2020a), thereby providing initial ev-
Multiple works have tried to overcome the convex barrier for piecewise linear activations (Raghu- nathan et al., 2018 Singh et al. 2019). Here, we focus on the single-neuron Linear Programming (LP) relaxation by Anderson et al. (2020). Unfortunately, its tightness comes at the price of expo- nentially many (in the number of variables) constraints. Therefore, existing dual solvers (Dvijotham) et al. 2018 Bunel et al. 2020a) are not easily applicable, limiting the scaling of the new relaxation.
"Equal contribution.

---

Published as a conference paper at ICLR 2021
We address this problem by presenting specialized dual solver for the relaxation by/Anderson et al. 2020), which realises its full potential by meeting the following desiderata:
By keeping active of dual variables, obtain sparse dual solver that recovers the strengths the original primal problem (Anderson et al., (020)) in the dual domain. In line with previous dual solvers, our approach yields valid bounds at anytime, leverages convolutional net- work structure and enjoys massive parallelism within a GPU implementation, resulting in better bounds in an order of magnitude less time than off- off-the-shelf solvers (Gurobi Optimization 2020
We present a unified dual treatment that includes both a linearly sized LP relaxation (Ehlers 2017) and the tighter formulation. As a consequence, our solver provides a wide range of speed- accuracy trade-offs: (i) it is competitive with dual approaches on the looser relaxation (Dvijotham et al., 2018; Bunel et al., 2020a); and (ii) it yields much tighter bounds if larger computa- tional budget is available. Owing to this flexibility, we show that our dual algorithm yields large complete verification gains compared to primal approaches (Anderson et al. 2020) and previous dual algorithms.
2 PRELIMINARIES: NEURAL NETWORK RELAXATIONS
We denote vectors by bold lower case letters (for example, x) and matrices by upper case letters (for example, W). We use @ for the Hadamard product, I-] for integer ranges, 1a for the indi- cator vector on condition a and brackets for intervals (Ik,, uz) and vector or matrix entries (x[i] Jii In addition, given W Rmxn and x E R" will employ W o x and W O as shorthands for respectively E, col, @ x and Σ, (W) x, where col, denotes the i-th column of matrix W.
Let C be the network input domain. Similar to Dvijotham et al. (2018) Bunel et al. (2020a), we assume that linear minimisation over C can be performed in closed-form. Our goal to compute bounds on the scalar output of a piecewise-lin linerr feedforward neural network. The tightest possible lower bound can be obtained by solving the following optimization problem:
min x,
In
S.t.
Xo E C,
(la)
RC+1 WK+1X8 + ba+1
k E [0, 11 k E [1, -1]-
(Ib)
X2 = a (8x)
(1c)
where the activation function a (8x) is seecewise-lineat &g E R. denote the outputs of the k- th linear layer ffully-coneccted convolutional) and activation function respectively, We and by denote its weight matrix and bias, ni is the number of activations at layer k. We will focus on the ReLU case (a (x) max (x, 0)), as common piecewise-linear functions can be expressed composition of ReLUs (Bunel et al. 20205).
Problem (I) is non-convex due to the activation function's non-linearity (Ic). As solving it is NP- hard (Katz et al., 2017), it is commonly approximated by a convex relaxation (see 44) The quality of the corresponding bounds, which is fundamental verification, depends on the tightness of the relaxation. Unfortunately, tight relaxations usually correspond to slower bounding procedures. We
first review a popular ReLU relaxation in 1. We then consider tighter one 2.2
2.1 PLANET RELAXATION
The so-called Planet relaxation (Ehlers) 2017) ) has enjoyed widespread use due to its amenability to efficient customised solvers (Dvijotham al. 2018 Bunel et al. (020a) and is the "relaxation of choice" for many works in the area (Bunel et al. 2020b Lu & Kumar (020).. Here, we describe it in its non-projected form Mr. the LP relaxation of the Big-M Mixed Integer Programming (MIP) formulation (Tjeng et al., 2019). Applying M& to problem [ results in:
min in s.t. x,2.2
Xo EC
* W Wx + b
k: [0,n - -1],
X > Ek. X4 S úg 3 Zki
(2)
X; S &a - 1 @1 zz),
(xt, Ec. 2x) E [l,, ux] x [e. an) x [0, 1] M₁ Ke [1, 1]
where 1. úg and lt Ug are intermediate bounds respectively on pre-activation variables & and post-activation variables Xt. These constants play important role in the structure of M& and,
2

---

Published as a conference paper at ICLR 2021
together with the relaxed binary constraints define box constraints the variables. We detail how to compute intermediate bounds in appendix E Projecting out auxiliary variables z results in the Planet relaxation (ef. appendix for details), which replaces (Tc) by its convex hull.
Problem (2). which is linearly-sized, can be easily solved via commercial black-box LP solvers (Bunel al. 2018). This does not scale-up well with the size of the neural network, mo- tivating the need for specialised solvers. Customised dual solvers have been designed by relaxing constraints (lb). (Ic) (Dvijotham al. (2018) or replacing ([c) by the Planet relaxation and employ-
ing Lagrangian Decomposition (Buncl ct al. 2020a). Both approaches result in bounds very close to optimality for problem (2) only a fraction of the runtime of off-the-shelf solvers.
2.2 A TIGHTER RELAXATION
A much tighter approximation of problem ①I than the Planet relaxation (: can be obtained by representing the convex hull of the composition (Tb) and (Icp rather than the convex hull of (Tc alone. A formulation of this type was recently introduced by Anderson et al. (2020).
Let us define Lx-1, UK-1 E R"=xx==- as: L [5,_] = 1-- -
wa[a.110 1 + [i==1[]] Walt,.j Additionally, let us introduce 2W. = (0, [0 xng 1, the set of all possible binary masks of weight matrix W,. and & 2 [0,1], which excludes the all-zero and all-one masks. The new representation
and Uk-1 = -u-1 w - 1 [w,,
results in the following primal problem:
min En s.t. x.8.x
Xo E C
Rkt1 WK+1IX + ba+1
k E [0, 1],
(xx, Rx. ZN) E ME
(We O In) x + 2 bb
= & 1]
X 5
WA O Te @ Ik-1) o 2x)
1Ia E &
WW& (1-11111) Ix) O ÚK-1)= o 2₁
Both Ma and A& yield valid MIP formulations for problem ①I when imposing integrality straints on z. However, the I.P relaxation of A& will yield tighter bounds. In the worst case, this tightness comes at the cost of exponentially many constraints: hand, given a set of primal assignments (x, z) that are not necessarily feasible for problem (3), one can efficiently compute the most violated constraint (if any) at that point. The mask associated to such constraint can be computed in linear-time (Anderson et al. 2020) as:
for each T & On the other
[i.j]] ((a 1i
We point out that As slightly differs from the original formulation of Anderson et al. 2020 which does not explicitly include pre-activation bounds i,, ûg (which we treat via M22). While this was implicitly addressed in practical applications (Botoeva et al.) (220)). not doing so has strong negative effect bound tightness, possibly the point of yielding looser bounds than problem (
In appendix [F we provide an example in which this is the case and extend the original derivation
by Anderson et al. (2000) recover A₁ as in problem 33)
Owing to the exponential number of constraints, problem cannot be solved as it is. As outlined
by Anderson et al.) (220), the availability of linear-time separation oracle (4) offers natural primal cutting plane algorithm, which can then be implemented in off-the-shelf solvers: solve the Big-M LP (2). then iteratively add the most violated constraints from Ar at the optimal solution. When applied to the verification of small neural networks via off-the-shelf MIP solvers, this leads
substantial gains with respect the looser Big-M relaxation (Anderson et al. 2020).
AN EFFICIENT DUAL SOLVER FOR THE TIGHTER RELAXATION
Inspired by the success of dual approaches looser relaxations (Bunel al. 2020a Dvijotham et al. (2019), we show that the formal verification gains by Anderson et al. (020)) (see 2.2) scale to larger networks if we solve the tighter relaxation the dual space. Due to the particular structure of the relaxation, a customised solver for problem (3) needs to meet a number of requirements.
Fact 1. In order to replicate the success of previous dual algorithms on looser relaxations, need a solver for problem 3 with the following properties: (i) sparsity: a memory cost linear in the number of network activations in spite of exponentially many constraints, (üi) tightness: the bounds
3

---

Published as a conference paper at ICLR 2021
should reflect the quality of those obtained in the primal space, (iti) anytime: low cost per iteration and valid bounds at each step.
The anytime requirement motivates dual solutions: any dual assignment yields a valid bound due to weak duality. Unfortunately, as shown in appendix [A] neither of the two dual derivations by [Bunel al. (2020a); Dvijotham al. (2018) readily satisty all desiderata at once. Therefore, we need
completely different approach. Let us introduce dual variables a, 3 and functions thereof:
fx) B) Okk W 0 -Et, BR, + Elar (WE+1 Ix++ BK
(k( EneS, (We O O (K-) BR, + BR,0 ④ ûr + BE, @ l
(5)
+ IIEE (Wa ② I o La-1) Bk. EieEs Ba. C bb:
where En, a shorthand for Σ ELE2.,. Starting from primal (3) we relax all constraints in AA, except box constraints (see 12.1). We obtain the following dual problem (derivation in appendix
where functions fk. 9k appear in inner products with primal variables Xt, Zk:
max d(c, B) (a,8)20
where:
d(a, 3) min 2, x
Ek= b a EKED (a(a,,)"xx Ek= (B)TT
C(x,2,0,3))
(6)
E Enee (W& @ La o 1x-1) Bx. + B3, (1x bx)) (bx))
s.t. Xg EC,
(k) [le. ux]= x 10, 1] E [1, 1].
This is again a challenging problem: the exponentially many constraints in the primal (3) are now associated to an exponential number of variables. Nevertheless, we show that the requirements of Fact[eea be met by operating on a restricted version of dual (6). To this end, we present Active Set, a specialised solver for the relaxation by Anderson et al. (2020) that is sparse, anytime and yields bounds reflecting the tightness of the new relaxation. Starting from the dual of problem (2), our solver iteratively adds variables to a small active set of dual variables 3s and solves the resulting reduced version of problem (6). We first describe our solver on a fixed BB and then outline how to iteratively modify the active set ($3.2). Pseudo-code can be found in appendix D
3.1 ACTIVE SET SOLVER
We want to solve a version of problem (6) for which the sums over the 4a masks of each layer k 0, we recover a novel dual solver for the Big-M relaxation @ (explicitly described in appendix B) which is employed as initialisation. Setting Br. 0, I & B& in (5) 6 and removing these from the formulation, we obtain:
are restricted B& E with B Ux By keeping
WT E+1 0 Σ1 E1.eB,u(o.)) BR.
f.. (o,
B (W @ B
EneB. (W& 6 (1 Ix) Ún-1) )0- Br... + 3k,0 @ úg + Br, O
(BB)
+ EI.B. (W& o I o Lx-1) 3x.1 + EneB, Br.l. @ by.
along with the reduced dual problem:
max dB(a,Bz)
where:
ds(a, BB) == min Lo(x, 2, x.2
(a.3e)20
Σ k= by Ex JB. (a,, Bz) xx Ex=1 95.6 (Bs) T2x
2,(x,2,a,, =
+ + EINE,, (Wn O I o La 1) D BR + BT,
(1. -ba)
s.t. EC
(Xx, 2x) [Rx, us] x [0, 1]
k k [1,
We can maximize (o, Bs), which is concave and non-smooth, via projected supergradient ascent or variants thereof, such as Adam (Kingma & Ba, 2015). In order to obtain a valid supergradient, need to perform the inner minimisation over the primals. Thanks to the structure of problem (8), the optimisation decomposes over the layers. For & E [I, n - 1]]. we can perform the minimisation in closed-form by driving the primals to their upper lower bounds depending the sign of their coefficients:
x 11. ûr + 1l (e. i.
z 11g. (l O 1.
(9)
As dual variables BR.1, indexed by Tg. B Up implicitly defines active set of variables Ag..
4

---

Published as a conference paper at ICLR 2021
The subproblem corresponding to X is different, as it involves linear minimization over Xy C:
X E argmin f Xo
s.t. Xg E C.
(10)
We assumed in 82 that (10) can be performed efficiently. We refer the reader to Bunel et al. (2020a)
for deseriptions of the minimisation when C' is a lo or &g ball, as common for adversarial examples.
Given (x", z²) above, the supergradient of dg(a,, (8) is a subset of the one for d( (), given by:
(a, B)- = W + bb: - x
Vo,,o((a,, (3) = x z o uR:
A d(a, (3) x (Wixk-1 + b&) + (1 - 2 3 1R.
(1I)
x (W& @ x + (We 4. Lx-1) o @ zit)
V Br d(a,
I BR
-2 by + (Wh o -1- (x) o Un-1) oz
for each ht E 10,n-1]] At each iteration, after taking step in the supergradient direction, the dual variables are projected the non-negative orthant by clipping negative values.
3.2 EXTENDING THE ACTIVE SET
We initialise the dual (6) with a tight bound on the Big-M relaxation by solving for (a, Bg) in (8) To satisfy the tightness requirement in Fact we then need to include constraints (via their Lagrangian multipliers) from the exponential family Ag into BN. Our goal is to tighten them as much as possible while keeping the active set small to save memory and compute.
The active set strategy is defined by a selection criterion for the 1 to be added to Bx. and the frequency of addition. In practice, we add the variables maximising the entries of supergradient V Vpr.,, d(ax, 3) after a fixed number of dual iterations. We now provide motivation for both choices.
Selection criterion ceed greedily and focus only on the immediate effect at the current iteration. Let us map a restricted set of dual variables Gg set of dual variables 3 for the full dual (6). We do by setting vari- ables not in the active set to 0: A8 = O, and 3 = Bs U Bg- Then, for each layer k, we add the set of variables Bk.rg maximising the corresponding entries of the supergradient of the full dual problem (6): Bx.rz E argmax a d(a, 1). Therefore, we use the subderivatives as a proxy for short-term improvement on the full dual objective d(a, 3)) Under primal interpretation, our selection criterion involves a call to the separation oracle (4) by Anderson et al. (2020).
The selection criterion needs to be computationally efficient. Thus, pro-
Proposition 1. BR.1z (as defined above) represents the Lagrangian multipliers associated to the most violated constraints from As at (x) E argminx, Lx(x, z, a, BB), the primal minimiser of the current restricted Lagrangian.
Proof. See appendix D.I
Frequency Finally, we need to decide the frequency at which to add variables to the active set.
Fact 2. Assume we obtained a dual solution (at, (at)))) argmax dg(a,, (a,, 3z) using Active Set on the current B. Then (x", .z") E argming,2 Lx(x, z, at (11) is not necessarily an optimal primal solution for the primal of the current restricted dual problem (Sherali & Choi 1996).
The primal of dg(a, (88) (restricted primal) is the problem obtained by setting &x <- B₁ in prob- lem (3). While the primal cutting plane algorithm by Anderson et al. (2020) calls the separation oracle at the optimal solution of the current restricted primal, Fact 2 shows that our selection criterion leads to different behaviour even at dual optimality for dæ (or, 33). Therefore, as we have no theoretical incentive to reach (approximate) subproblem convergence, we add variables after a fixed tunable number of supergradient iterations. Furthermore, we can add more than one variable "at once" by running the oracle (4) repeatedly for a number of iterations.
We conclude this section by pointing out that, while recovering primal optima is possible in prin- ciple (Sherali & Choil 1996). doing so would require dual convergence on each restricted dual problem (8). As the main advantage of dual approaches (Dvijotham et al. 2018 Bunel et al. (020a) is their ability to quickly achieve tight bounds (rather than formal optimality). adapting the selection criterion mirror the primal cutting plane algorithm would defeat the purpose of Active Set.
adding a single Ie mask B& extends 8g by ng variables: for each neuron layer
5

---

Published as a conference paper at ICLR 2021
3.3 IMPLEMENTATION DETAILS, TECHNICAL CHALLENGES
Analogously to previous dual algorithms Dvijotham et al. 2018 Bunel al. 2020a). our approach can leverage the massive parallelism offered by modern GPU architectures in three different ways. First, we execute in parallel the computations of lower and upper bounds relative to all the neurons of a given layer. Second, in complete verification, we can batch over the different Branch and Bound (BaB) subproblems. Third, as most of our solver relies standard linear algebra operations employed during the forward and backward passes of neural networks, we can exploit the highly optimised implementations commonly found in modern deep learning frameworks.
An exception are what we call "masked" forward/backwward passes: operations of the form (W& 3 11)xx or (W& O 1x)² XK+1, which are needed whenever dealing with constraints from Ag. In our solver, they appear B& & (see equations (8) ). Masked passes require a customised lower-level implementation for a proper treatment of convolutional layers, detailed in appendix G
4 RELATED WORK
In addition to those described in 2 many other relaxations have been proposed in the literature. In fact, all bounding methods are equivalent to solving some convex relaxation of neural network. This holds for conceptually different ideas such as bound propagation (Gowal et al. 2018). specific dual assignments (Wong& & Kolter 2018), dual formulations based on Lagrangian Relaxation (Dvi- jotham et al.) 2018) or Decomposition (Bunel et al., (020a). The degree of tightness varies greatly: from looser relaxations associated to dlsed-form methods (Gowal et al. 2018 Weng et al. 2018 Wong & Kolter 2018) to tighter formulations based on Semi-Defimit Programming (SDP) (Raghu-
nathan et al.1
The speed of closed-form approaches results from simplifying the triangle-shaped feasible region of the Planet relaxation (1. (Singh et al. 2018 Wang et al. 2018). On the other hand, tighter relaxations
expressive than the linearly- sized LP by (Ehlers 2017). The SDP formulation by Raghunathan et al. (2018) can represent interactions between activations in the same layer. Simi- larly, et al. (2019) tighten the Planet relaxation by considering the convex hull of the union of polyhedra relative to k ReLUs of a given layer at once. Alternatively, tighter LPs can be obtained by considering the ReLU together with the affine operator before it: standard MIP techniques (Jeroslow 1987) lead to formulation that is quadratic in the number of variables (see appendix (.2). The relax- ation by Anderson et al. (2020) detailed in (.. is a more convenient representation of the same set.
By projecting out the auxiliary z variables, Tjandraatmadja et al., 2020) recently introduced another formulation equivalent to the one by Anderson et al. (2020), with half as many variables and a linear factor more constraints compared to what described in (2.2 Therefore, the relationship between the two formulations mirrors the one between the Planet and Big-M relaxations (see appendix Our dual derivation and the Active Set algorithm
be adapted operate the projected relaxations.
Specialised dual solvers significantly improve in bounding efficiency with respect to off-the-shelf solvers for both I.P (Bunel al., 2020a) and SDPP formulations (Dvijotham et al.. (019).. Therefore, the design of similar solvers for other tight relaxations is an interesting line of future research. We contribute with a specialised dual solver for the relaxation by Anderson al. (2020) (33) In what follows, we demonstrate empirically that by seamlessly transitioning from the Planet relaxation to the tighter formulation, we can obtain large incomplete and complete verification improvements.
EXPERIMENTS
We empirically demonstrate the effectiveness of our method under two settings. On incomplete ver- ification (5..T). we assess the speed and quality of bounds compared to other bounding algorithms. On complete verification (55.)). we examine whether
apcud-acyuracy trade-offss correspond to faster exact verification. Our implementation is based on Pytorch (Paszke et al. (017) and is avail- able at httpt:///gttttttee htt.://t//aa/a//teac..ooooo//ooa
wwt htwpr//rooroooooooooooooooooooen
5.1 INCOMPLETE VERIFICATION
We evaluate incomplete verification performance by upper bounding the robustness margin (the dif- ference between the ground truth logit and the other logits) to adversarial perturbations (Szegedy ct al. 2014) on the CIFAR-10 test set (Krizhessky & Hinton 2009). the upper bound is negative, can certify the network's vulnerability to adversarial perturbations. We replicate the experimen-

---

Published as a conference paper at ICLR 2021
10
101
Method
Imprewentt inproemmett from oit
0.5
0.0
-0.5
-1.0
Gurohi
Gumbi
BDD
Bly-M
Active Set
Artive
Active steps
Artive
Active Set CPU
400
650
Figure 1: Upper plot: distribution of runtime in seconds. Lower plot: difference with the bounds obtained by Gurobi with a cut from As per neuron; higher is better. Results for the SGD-train trained network from Bunel] al. (2020a). The width given value represents the proportion of problems for which this is the result. Comparing Active Scts with 1650 steps to Gurobi Cut, tighter bounds are achieved with a smaller runtime.
tal setting from Bunel et al. (020a). The networks correspond to the small network architecture from WWong& & Kolter (2018). Here, present results for network trained via standard SGD and cross entropy loss, with no modification to the objective for robustness. Perturbations for this network lie in foo norm ball with radius ever 5/255 (which is hence lower than commonly em- ployed radii for robustly trained networks). In appendix [1] we provide additional CIFAR-10 results adversarially trained network using the method by Madrr et (2018).. and MNIST (LeCun iet al. 1998), for network adversarially trained with the algorithm by WWon & Kolter (2018).
We compare both against previous dual iterative methods and Gurobi (Gurobi Optimization 2020). the commercial black-box solver employed by Anderson et al. (2020). For GGrobi-basee baselines, Planet means solving the Planet Ehlers (2017) relaxation of the network. while Gurobi starts from the Big-M relaxation and adds constraints from A& in cutting-plane fashion, as the original primal algorithm by [Anderson et al. (2020). We run both on CPU threads. Amongst dual iterative methods, run an Nvidia Titan Xp GPU, we compare with BDD+, the recent proximal-based solver by Bunel et al. 2020a),, operating on Lagrangian Decomposition dual the Planet relaxation. As operate (a subset of) the data by Bunel et al. (2020a), based approach and the one by Dvijotham et al. (2018). as they both perform worse than BDD+ omit cheaper (and looser) methods, like interval propagation Gowal et al. (2018) and the one by WWng & Kolter (2018). Active Set denotes our Active Set reduces to Big-M, a solver for the non-projected Planet relaxation (appendix B which can be seen as Active Set's initialiser. In line with previous bounding algorithms (Bunel et al. 2020a), we employ Adam updates (Kingma) & Ba 2015) for supergradient-tynn type methods due to their faster empirical convergence. Finally, we complement the comparison with Gurobi-basee methods by rumning Active Set on 4 CPU threads (Active Set CPU). Further details, including hyper-paree paper-parammters, can be found in appendix
omit both their supergradient-
(Bunel al.] 2020a). For the same reason,
solver for problem (3) described in 43.1 By keeping B
Figure shows the distribution of runtime and the bound improvement with respect to Gurobi cut for the SGD-trained network. For Gurobi cut, we only add the single most violated cut from As per neuron, due to the cost of repeatedly solving the LP.
tuned BDD+ and Big-M, the dual methods operating the weaker relaxation ( to have the same average runtime. They obtain bounds comparable to Gurobi Planet in one order less time. Initialised from 500 Big-M iterations, at 600
Timigg th
Sae 15 Paret
Thiing In
Inprovement from: Panst
1.0-
3.5
Aodies stera
844
3.5
800+ slgss
stppt
BBN enss
0ig mass
(a) Comparison of runtime (left) and gap to Gurobi Planet bounds (right). For latter, lower better.
(b) Comparison of runtime (left) and difference with the Gurobi Planct bounds (right). For the latter, higher better.
Figure 2: Pointwise comparison for subset of the methods on the data presented in Figure Darker colour shades mean higher point density (on a logarithmic scale). The oblique dotted line corresponds the equality.

---

Published as a conference paper at ICLR 2021
iterations Active Set already achieves better bounds on average than Gurobi cut in around 1/20th of the time. With a computational budget twice as large (1050 iterations) or four times as large (1650 iterations). the bounds significantly improve over Gurobi cut in still a fraction of the time. As we empirically demonstrate in appendix [ the tightness of the Active Set bounds is strongly linked to our active set strategy (33.2). Remarkably, even if our method is specifically designed to take advantage of GPU acceleration, executing it on CPU proves to be strongly competitive with Gurobi cut, producing better bounds in less time for the benchmark of Figure |
Figure 2 shows pointwise comparisons for a subset of the methods of Figure I on the same data. Figure 2a shows the gap to the (Gurobi) Planet bound for BDD+ and our Big-M solver. Surprisingly, our Big-M solver is competitive with BDD+, achieving on average better bounds than BDD+, in the same time. Figure shows the improvement over Planet bounds for Big-M and Active Set. The latter achieves markedly better bounds than Big-M in the same time, demonstrating the benefit of operating (at least partly) on the tighter dual
5.2 COMPLETE VERIFICATION
We next evaluate the performance on complete verification, verifying the adversarial robustness of network to perturbations in fo norm on subset of the dataset by [u & Kumar (2020), replicating the experimental setting from Bunel et al. (2020a). The dataset associates a different perturbation ra- dius Everif to each CIFAR- 10 image, so to create challenging verification properties. Its difficulty makes the dataset an appropriate testing ground for tighter relaxations like the one by Anderron et al. (2020) 2.. Further details, including network architectures, can be found in appendix
Here, we aim to solve the non-convex problem ①T directly, rather than an approximation like 15.1 In order to do so,
BaSBR, branch and bound algorithm from [Bunel et al. (2020b). Branch and Bound works by dividing the problem domain into subproblems (branching) and bounding the local minimum over those domains. Any domain which cannot contain the global lower bound is pruned away, whereas the others are kept and branched over. In BaBSR, branching is carried out by splitting an unfixed ReLU into its passing and blocking phases. The ReLU which induces maximum change in the domain's 's lower bound, when made unambiguous, is selected for splitting.
A fundamental component of a BaB method is the bounding algorithm, which is, in general, the computational bottleneck (Lu & Kumar 2020). Therefore, we compare the effect on final verifica- tion time of using the different bounding methods in 35.T within BaBSR. In addition, we evaluate MIP Ag. which encodes problem (I) as a Big-M MIP (Tjeng et al. 2019) and solves it in Gurobi by adding cutting planes from Av. analogously to the original experiments by Anderson et al. (2020).
Finally, we also compare against ERAN (Singh al. 2020), a state-of-the-ar complete verification toolbox: results on the dataset by Lu & Kumar (2200) are taken from the recent VNN-COMP com- petition (VNN-COMP) 2020). We use 100 iterations for Active Set, 100 iterations for BDD+ and 180 iterations for Big-M. For dual iterative algorithms, we solve 300 subproblems at once for the base network and 200 for the deep and wide networks (see (3.3). Additionally, dual variables are initialised from their parent node's 's bounding computation. As in Bunel et al. (2020a), the time-limit kept at one hour. Due to the difference in computational cost between algorithms operating on the tighter relaxation by Anderson et al. (2020) and the other bounding algorithm we also experiment with a stratified version of the bounding within BaBSR. We devise a set of heuristics to determine
For Active Set, this is partly due to the masked forward/backward pass described in appendix G
Wide large
model
BCC- Bvtth
ACC- Buttl
ABC: Bvtth
baita s-reece to
Piaece mecue
15 of propertton vertfed
--n rianec cut owiee
68
of propertten
46-
46
Computation
Computation time
Computation
Figure 3: Cactus plots on properties from Lu & Kumar (2020), displaying the percentage of solved properties function of runtime. Baselines are represented by dotted lines.
8

---

Published as a conference paper at ICLR 2021
Bine
Wide
Deep
Method
sub-prohhlems
% limeout
sub-proolems
76mmmoout
sub-prohlemm
Timenut
BDD+ BABSR
883.55
82.699.00
22.00
568.25
43775188
13.00
281.47
10763 48
5,00
BOG-M BARSR
826.60
68 582.00
19.00
533,79
35 877.72224
12.00
253.37
9346.78
4,00
Set BABSS
+22.32
9471.99
7,00
169.73
1873.36
3,00
227.26
2302_16
2,00
BIG SE A. SEE BABSR
402.88
11 10/00..5 8/1
7.00
179.73
3712.62
3,00
197.99
3086 62
2,00
PLANET
BABSK
1191.48
2014 28
14,00
1272.99
1352,42
10.00
704.59
677. 74
3,00
.A.
3227.50
$2.00
100.93
3339.37
434.57
ERAN
805.89
5.00
632.12
9.00
545.72
0.00
Table I: We compare average solving time, average number of solved sub-problems and the percentage of timed out properties on data from Lu & Kumar (220)). The best dual iterative method is highlighted in bold.
whether given subproblem is easy (therefore looser bounds are sufficient) or whether we need to operate on the tighter relaxation. Instances of this approach are Big-M + Active Set and Gurobi
Planet + Gurobi cut. Further details are provided appendix H
Figure 3 and Table 1 show that Big-M performs competitively with BDD+. Active Set verifies
larger share of properties than the methods operating on the looser formulation (2), demonstrating the benefit of tighter bounds (55.1) in complete verification. On the other hand, the poor performance of MIP + As and Gurobi Planet + Gurobi cut, tied to scaling limitations of off-the-shelf solvers, shows that tighter bounds are effective only if they can be computed efficiently. Nevertheless, the difference in performance between the two Gurobi-based methods confirms that customised Branch and Bound solvers (BaBSR) preferable to generic MIP solvers, as observed by Bunell et al. (020b) on the looser Planet relaxation. Moreover, the stratified bounding system allows us to retain the speed of Big-M on easier properties, without excessively sacrificing Active Set's gains on the harder ones. Finally, while ERAN verifies 2% more properties than Active Set on two networks, BaBSR (with any dual bounding algorithm) is faster on most of the properties. BaBSR- -bsed results could be further improved by employing the learned branching strategy presented by Lu & Kumar 2220 in this work, we focused on the bounding component branch and bound.
DISCUSSION
The vast majority of neural network bounding algorithms focuses on (solving or loosening) a pop- ular triangle-shaped relaxation, referred to as the "convex barrier" for verification. Relaxations that are tighter than this convex harrier have been recently introduced, but their complexity hinders appli- cability. We have presented Active Set, sparse dual solver for one such relaxation, and empirically demonstrated that it yields significant formal verification speed-ups. Our results show that scalable tightness is key to the efficiency of neural network verification and instrumental in the definition of a more appropriate "convex barrier". We believe that new customised solvers for similarly tight relaxations are crucial avenue for future research in the area, possibly beyond piecewise-lines networks. Finally, as it is inevitable that tighter bounds will come at larger computational cost, future verification systems will be required to recognise a priori whether tight bounds are needed for given property. A possible solution to this problem could rely on learning algorithms.
ACKNOWLEDGMENTS
ADP was supported by the EPSRC Centre for Doctoral Training Autonomous Intelligent Ma- chines and Systems, grant EP//011987/,, and an IBM PhD fellowship. HSB was supported using Tencent studentship through the University of Oxford.
REFERENCES
Ross Anderson, Jocy Huchette, Christian Tjandraatmadja, and Juan Pablo Vielma. Strong mixed- integer programming formulations for trained neural networks. Conference Integer Program- ming and Combinatorial Optimization, 2019.
Ross Anderson, Joey Huchette, Will Ma, Christian Tjandraatmadja, and Juan Pablo Vielma. Strong mixed-integer programming formulations for trained neural networks. Mathematical Program- ming 2020.
Elena Botoeva, Panagiotis Kouvaros, Jan Kronqvist, Alessio Lomuscio, and Ruth Misener. Effi- cient verification of ReLU- neural networks via dependency analysis. AAAI Conference on Artificial Intelligence, 2020.

---

Published as a conference paper at ICLR 2021
Rudy Bunel, Ilker Turkaslan, Philip HS Torr, Pushmeet Kohli, and M Pawan Kumar. A unified view of piecewise linear neural network verification. Neural Information Processing Systems, 2018.
Rudy Bunel, Alessandro De Palma, Alban Desmaison, Krishnamurthy Dvijotham, Pushmeet Kohli, Philip HS Torr, and M Pawan Kumar. Lagrangian decomposition for neural network verification. Conference on Uncertainty in Artificial Intelligence, 2020a.
Rudy Bunel, Jingyue Lu, liker Turkaslan, P Kohli, P Torr, and M Pawan Kumar. Branch and bound for piecewise linear neural network verification. Journal of Machine Learning Research, 21(2020), 2020b.
Nicholas Carlini and David Wagner. Towards evaluating the robustness of neural networks. 2017 IEEE Symposium on Security and Privacy, 2017.
Krishnamurthy Dvijotham, Robert Stanforth, Sven Gowal, Timothy Mann, and Pushmeet Kohli. A dual approach scalable verification of deep networks. Uncertainty Artificial Intelligence, 2018.
Krishnamurthy Dvijotham, Robert Stanforth, Sven Gowal, Chongli Qin, Soham De, and Pushmeet Kohli. Efficient neural network verification with exactness characterization. Uncertainty in Arti- ficial Intelligence, 2019.
Ruediger Ehlers. Formal verification of piece-wise wise linear feed-forward neural networks. Automated Technology for Verification and Analysis, 2017.
Sven Gowal, Krishnamurthy Dvijotham, Robert Stanforth, Rudy Bunel, Chongli Qin, Jonathan Ue- sato, Timothy Mann, and Pushmeet Kohli. On the effectiveness of interval bound propagation for training verifiably robust models. SECML NeurlPS, 2018.
LLC Gurobi Optimization. Gurobi optimizer reference manual, 2020. URL http: //www.
gurobi com
R. G. Jeroslow. Representability in mixed integer programming, i: Characterization results. Discrete Applied Mathematics, 1987.
Guy Katz, Clark Barrett, David Dill, Kyle Julian, and Mykel Kochenderfer. Reluplex: An efficient SMT solver for verifying deep neural networks. International Conference Computer-Aided Verification, 2017.
Diederik P. Kingma and Jimmy Ba. Adam: A method for stochastic optmizzation. International Conference on Learning Representations, 2015.
A. Krizhevsky and G. Hinton. Learning multiple layers of features from tiny images. Master's thesis, Department of Computer Science, University of Toronto, 2009.
Yann LeCun, Léon Bottou, Yoshua Bengio, and Patrick Haffner. Gradient-basee learning applied to document recognition. IEEE, 1998.
Claude Lemaréchal. Lagrangian relaxation. In Computational combinatorial compiitttooacccooomiioocttto oatpiiiiaiiiiiooiiiiiiiiiiioo pp. 112- 156. Springer, 2001.
Jingyue Lu and M Pawan Kumar. Neural network branching for neural network verification. In International Conference Learning Representations, 2020.
Aleksander Madry, Aleksandar Makelov, Ludwig Schmidt, Dimitris Tsipras, and Adrian Vladu. To- wards deep learning models resistant to adversarial attacks. International Conference on Learning Representations, 2018.
Adam Paszke, Sam Gross, Soumith Chintala, Gregory Chanan, Edward Yang, Zachary DeVito, Zeming Lin, Alban Desmaison, Luca Antiga, and Adam Lerer. Automatic differentiation in pytorch. NIPS Autodiff Workshop, 2017.
Aditi Raghunathan, Jacob Steinhardt, and Percy S Liang. Semidefinite relaxations for certifying robustness to adversarial examples. Neural Information Processing Systems, 2018
10

---

Published as a conference paper at ICLR 2021
Hadi Salman, Greg Yang. Huan Zhang, Cho-Jui Hsieh, and Pengchuan Zhang. A convex relaxation barrier to tight robustness verification of neural networks. Neural Information Processing Systems, 2019.
Hanif D. Sherali and Gyunghyun Choi. Recovery of primal solutions when using subgradient op- timization methods to solve lagrangian duals of linear programs. Operations Research Letters, 1996.
Gagandeep Singh, Timon Gehr, Matthew Mirman, Markus Püschel, and Martin Vechev. Fast and effecttie robustness certification. Neural Information Processing Systems, 2018.
Gagandeep Singh, Rupanshu Ganvir, Markus Püschel, and Martin Vechev. Beyond the single neuron convex barrier for neural network certification. Neural Information Processing Systems, 2019.
Gagandeep Singh, Jonathan Maurer, Christoph Müller, Matthew Mirman, Timon Gehr, Adrian Hoff- mann, Petar Tsankov, Dana Drachsler Cohen, Markus Püschel, and Martin Vechev. ETH robust- ness analyzer for neural networks (ERAN), 2020. URL htts:::::: htt/ggttuubl com et sri/
eran
Christian Szegedy, Wojciech Zaremha, Ilya Sutskever, Joan Bruna, Dumitru Erhan, Ian Goodfellow, and Rob Fergus. Intriguing properties of neural networks. International Conference Learning Representations, 2014.
Christian Tjandraatmadja, Ross Anderson, Joey Huchette, Will Ma, Krunal Patel, and Juan Pablo Vielma. The convex relaxation barrier, revisited: Tightened single-neuron relaxations for neural network verification. arXiv preprint arXiv:2006. 14076, 2020.
Vincent Tjeng, Kai Xiao, and Russ Tedrake. Evaluating robustness of neural networks with mixed integer programming. International Conference on Learning Representations, 2019.
VNN-COMP. International verification of neural networks competition (VNN-COMP), Verification of Neural Networks workshop at the International Conference Computer-Aided Verification, 2020. utLpttps::// /ttsssttsssssse//ee///eeeeseeees eesiiess.ooogle....ooooo
google..o...ooooooooooooooo..... . com/ www./www//////..o....o........ee view/vvn20 / www.nn//v/n//npppp 'nnccopp
Shiqi Wang, Kexin Pei, Justin Whitehouse, Junfeng Yang, and Suman Jana. Efficient formal safety analysis of neural networks. Neural Information Processing Systems, 2018.
Stefan Webb, Tom Rainforth, Yee Whye Teh, and M Pawan Kumar. A statistical approach assess- ing neural network robustness. International Conference on Learning Representations, 2019.
Tsui- Weng, Huan Zhang, Hongge Chen, Zhao Song, Cho-Jui Hsich, Duane Boning, Inderjit Dhillon, and Luca Daniel. Towards fast computation of certified robustness for relu networks. International Conference Machine Learning, 2018.
Eric Wong and Zico Kolter. Provable defenses against adversarial examples via the convex outer adversarial polytope. International Conference on Machine Learning, 2018.

---

Published as a conference paper at ICLR 2021
A LIMITATIONS OF PREVIOUS DUAL APPROACHES
In this section, we show that previous dual derivations (Bunel et al. 2020a Dvijotham et al. 2018) violate Fact T Therefore, they not efficiently applicable problem motivating derivation in section
(011)), which relies on relaxing equality con- straints (Ib). Ic from the original non-convex problem ①I_ Dvijotham et al. (2018) prove that this relaxation corresponds to solving convex problem (2). which equivalent to the Planet relax- ation (Ehlers, 2017), to which the original proof refers. As we would like to solve tighter problem 33. the derivation is not directly applicable. Relying on intuition from convex analysis applied to
We start from the approach by Dvijotham
duality gaps (Lemaréchal 2001), conjecture that relaxing the composition (Ic) ([b) might tighten the primal problem equivalent to the relaxation, obtaining the following dual:
min
Wn bn
Σ k= k n
Mi T
(x+ max (WAXL-1 bR. 0})
1
(12)
s.t.
1 xx S Su
Xa ECC.
Unfortunately dual 12) requires an LP (the inner minimisation over x, which in this case does not decompose layers) to be solved exactly obtain supergradient and any time valid bound is needed. This is markedly different from the original dual by Dvijotham et al. (2018), which had an efficient closed-form for the inner problems.
The derivation by Bunel et al. (2020a). instead, operates by substituting (Ic) with its convex hull and solving its Lagrangian Decomposition dual. The Decomposition dual for the convex hull of (Ic) (1b (i.c., Ag) takes the following form:
max min A, DR
Pk (XB,1 XA,A)
P
X.2
(13)
s.t. Xo C, (XB... WnXAn-1 + b Zx) E Adec, E [1, 1].
where Adse,k corresponds to Ag with the following substitutions: Xx - XB.ks and & W.XAn-- b, It can be easily AAA:; XB,,, for each layer le >> 0) are an exponentially sized I.. Again, this differs from the original dual the Planet relaxation (Bunel et al. 2020a), which had efficient closed-form for the inner problems.
that the inner problems (the inner minimisation over
B DUAL INITIALISATION
Algorithm Big-M solver
1: function wwW.COOPEE..UUNDS(((((. bu. ls, ug lu, (=I..)
2: Initialise duals a ( using interval propagation bounds (Gowal al.
2018
3: for -1]
4:
x²,2" argming,x LM(x, a¹, 3 using [6
5: Compute supergradient using (18)
6: a!+1 B+++1 Adam's update rule (Kingma & Ba 2015)
7:
a (+1 3+++1 max(a 0), 1 0)
(dual projection)
8: end for
9:
return min LM a, aB,W%
10: end function
As shown in section 3 our Active Set solver yields dual solver for the Big-M relaxation 2 if the active set is kept empty throughout execution. As indeed B @ for the first Active Set iterations (see algorithm section D), the Big-M solver can be thought of as dual initialisation. Furthermore,
12

---

Published as a conference paper at ICLR 2021
demonstrate experimentally in S that, when used stand-alone solver, our Big-M solver is competitive with previous dual algorithms for problem (2). The goal of this section is to explicitly describe the Big-M solver, which is summarised in algorithm [ We point out that, in the notation of restricted variable sets from section 3.1 BM == fg.
We describe the equivalence between the Big-M and Planet relaxations, before presenting the solver in section B.3 and the dual it operates on section B.2
B.1 EQUIVALENCE TO PLANET
As previously shown (Bunel et al. 2018). the Big-M relaxation (Mg. when considering the k-th layer only) in problem 2 is equivalent to the Planet relaxation by [Ehlers (2017). Then, due to strong duality, our Big-M solver (section (3.2) and the solvers by Bunel al. (2020a); Dvijotham et al. 2018 will all converge to the bounds from the solution of problem (2). In fact, the Decomposition based method (Bunel et al. 2020a) directly operates on the Planet relaxation, while (vijothhi et al. (2018) prove that their dual equivalent to doing
On the k-th layer, the Planet relaxation takes the following form:
if l. S 0 and úg 2
X 2 0,
X; 2 &k:
X S ük @ (Rt i o (1/(0xx [a))
Pk =
ûg 0::
(14)
X = 0.
ig 2 0:
X₁ = &
It can be seen that Pk Projx,s (M&), where Projx.* denotes projection the & hyperplane. In fact, as Zς does not appear the objective of the primal formulation (2), but only in the constraints, this means assigning the value that allows the largest possible feasible region. This is trivial for passing or blocking ReLUs. For the ambiguous case, instead, Figure 4 (on a single ReLU) shows that zk fa=fa is the correct assignment.
B.2 BIG-M DUAL
As evident from problem 3 A S MM. If we relax all constraints in Mr (except, again, the box constraints). we are going to obtain a dual with striet subset of the variables in problem (6). The Big-M dual is specific instance of the Active Set dual ( where B 0, and it takes the following
27
34 1 8
0. ④
C,(āg. x))
26
Figure 4: Mr plotted on the (zu, x1) plane, under the assumption that i. < and ú 2
13

---

Published as a conference paper at ICLR 2021
form:
max (a BM) (a,9)20
where:
dm(a, BM) = min (x, 2, a, BM), x,2
k=0 (akk W 1++10++11000 (3k.0 + +B,x11 W -1x+1(++111)) " X
LM(x, 2, Cz BM)
+ k= bf Σk=1 B3 U + Bk,1 2 lx)
Z
M (he b&)² 1x))3Ba. Sx.1
s.t.
Xo E C,,
(xk 2k) E [1s uu x 00, 1]
(15)
B.3 BIG-M SOLVER
We initialise dual variables from interval propagation bounds (Gowal al. 2018): this can be easily done by setting all dual variables except On, to 0. Then, we can maximize dyt (ox, (3 via projected supergradient ascent, exactly described in section 3.1 on a generic active set B. All the
computations in the solver follow from keeping B = 0 in 3.1 We explicitly report them here for the reader's convenience.
Let us define the following shorthand for the primal coefficients:
fM,,(a, BM) = (aus- W k+1 0++1 (Bx.0 + Ar.1
W (1,,))
g.,.((3m) 3k.0 C ú; + Br, ll.
The minimisation of the Lagrangian L. (x, z, a over the primals for [[. n 1] is as follows:
x (a. ûg + (M.+(. l z
ll. @
(16)
For 0, instead (assuming, as [3.1 that this can be done efficiently):
xij
argmin Xg
fM, BM)
EC
(17)
The supergradient over the Big-M dual variables a, Bk,u, Br, is computed exactly as $3.1 and is again a subset of the supergradient of the full dual problem (6). We report it for completeness. For each A E [0, n - Il:
V Va, a(a 3) = WaxF-1 + ba -x-
V V 3) Xx Zh ? ús:
(18)
V Xx (WAXA-1 + b&) (1 2xk) ④ i
C DUAL DERIVATIONS
We now derive problem (6), the dual of the full relaxation by Anderson (2020) described
in equation (3). The Active Set (equation (D) and Big-M duals (equation (15) can be obtained
by removing BR.I. T E & B& and Be.I. Ti E &k. respectively. We employ the following Lagrangian multipliers:
X₁ 2 8 = Olk.
X < U 25 7 BL
X S 84 - 1 o (1 - 2₁)) = B
(W₂ & 46) X xx- + bj
Xb <
WA @ Is 2 L-11
Zz)
= 1k.
We @ (1 Ix) 1 2k
14

---

Published as a conference paper at ICLR 2021
and obtain, Lagrangian (using & WAXx-1 bb):
nリー k= ag (WaXk-1 + by -x6) -k=l 0% (x: - 2 2 ūs)
(Wh @ 11 2 (--1) o 2z) (Wn 2 12) Xk-1
[(x, 2, a,
Dk=1 EinEE, Sh..e
-by 4 22 (Wx @ (1 Ix) o
24 Xx
Ex=l BK, Xv: -WAXK- + b&) + (1 2₁) @ 1s) Wn + b,.
Let use Σr, shorthand for EIεEEu(OO.1) If collect the terms with respect to the primal variables and employ dummy variables On = 0, Bo 0, On 1, Bn = 0, obtain:
W W+100++1 1111-En, i Bu.te
Ek=0
+ + Jet (Wa+1 @ (k+11 ( Bb+1,
Xx
T
Eneer Bkl ① bb.- + Bk,1 ① 1u+to 3k,0 ① us
L(x, a)
>
Eneen (We 3 1x ② Lx-1) BR
km
Zx
EinEE, (W& (1 1x) о ÚK-1) BE,
+ + az + ERSl EieEs, (W& @ L @ Lx-1) O Br. + 18, - by))
which corresponds to the form shown in problem (6).
D IMPLEMENTATION DETAILS FOR ACTIVE SET METHOD
From the high-level perspective, our Active Set solver proceeds by repeatedly solving modified instances of problem 6 where the exponential set Eg is replaced by a fixed (small) set active variables B. The full solver procedure is summarised in algorithm 2
Algorithm 2 Active Set solver
1: function ACTIVESET COMPUTE BOUNDS([W,, by Le, ug lR uR 1==
2: Initialise duals a P B using Algorithm
3: Set BE,,, 0, T E &
4: ==
5: for ibbadditions do for / E 1
6:
7:
x",z" argminx,x Lx(x, a',3%) using 9 10)
8:
if S nb vars to add then
9:
For each layer k, add output of 44 called (x", 2²) to Be
10:
end
Compute supergradient using II
12:
a ++1 3 +1 <- Adam's 's update rule Kingma & Ba 2015) (dual projection)
13:
a ++1 B < m 0). max (
end for end for
15:
16: return minx Lx(x, z, aT, .
17:
end function
We conclude this section by proving the primal interpretation of the selection criterion for adding new set of variables to B.
D. ACTIVE SET SELECTION CRITERION
We map a restricted set of dual variables By to a set of dual variables /3 for the full dual (6) by setting variables not in the active set to Agg = 0, and 3 = Gs U Ag.
Proposition. Let BR.Iz be the of dual variables maximising the corresponding entries of the supergradient of the full dual problem ( Bk.J E argmax (V d(a,, 13) GR. represents the Lagrangian multipliers associated to the most violated constraints from A at (x" (2") ) E argminxx, Ls(x, 2, at, Bs), the primal minimiser of the current restricted Lagrangian.
15

---

Published as a conference paper at ICLR 2021
Proof. In the following, (x*,, z=) denotes the points introduced in the statement. Recall the definition of V Br,,, d(a, 3) in equation which applies beyond the current active set:
X (W& o 1x) - + (Wa @ la o LB-1) o (i1) 1 @ by + (W4 の Tx) Uk-1) o
z ) In E &k.
V Ba. d(
We want to compute [ E argmax 12 (V [V , that is:
- の 1 x - (Wx I 1 - 1Z-1)
)
T 1.
I argmaax IxEE2
-2% by + ( o (1 I&) Uk-1) o ozk
By removing the terms that do not depend on Lg. we obtain:
((W (W& xx)x (Wa O H Ia @ 1 2 ) 1.
max E& (Wa O 1 1 (
Let us denote the i-th row of We and Is by wis and i,a, respectively, and define Ex]i] 2 \\ {0, 10. The optimisation decomposes over each such row: we thus focus on the optimi- sation problem for the supergradient's i-th entry. Collecting the mask, we get:
((1-zi)1)) @ 17,1l + A 11-11 - x xi- 11) (ill Ia Aali..l. 5l.
Σ
max E
As the solution to the problem above is obtained by setting [i if its coefficient is positive otherwise,
and 1. / Anderson et al. (2020) on (x',z) z²). Hence, in addition to being the mask associated to (k.rg= the variable set maximising the supergradient, I corresponds to the most violated constraint from A at (x*, 2.).
see that the optimal L corresponds calling oracle 4 by
E INTERMEDIATE BOUNDS
A crucial quantity in both ReLU relaxations (M₃ and A,)) are intermediate pre-activation bounds le, ûr In practice, they are computed by solving relaxation Cs. (which might be Ma. A., or something looser) of ( over subsets of the network (Bunel al. 2020a). For 1, this means solving the following problem (separately, for each entry 1[7)):
min [ x.8.2
Xo C
(19)
& = W + bb-1,
E [0, 1] A E [1, 1]
(Xx, 82, E Cs
As (19) needs to be solved twice for each neuron (lower and upper bounds, changing the sign of the last layer's weights) rather than once in 3 depending
the computational budget, Cx might be looser than the relaxation employed for the last layer bounds (in our case, AA).. In all our experiments, we compute intermediate bounds the tightest bounds between the method by Wong & Kolter (2018) and Interval Propagation (Gowal et al. 2018).
Once pre-activation bounds are available, post-activation bounds can be simply computed as 1 max (1,, 0), u, = max(ū;, 0).
F PRE-ACTIVATION BOUNDS IN AA
We now highlight the importance of explicit treatment of pre-activation bounds in the context of the relaxation by/ Anndrrso et al. (2020). In F.I we will show through an example that, without a separate pre-activation bounds treatment, A could be looser than the less computationally expen- sive M2 relaxation. We then ($F2) justify our specific pre-activation bounds treatment by extending the original proof by Anderson al. 2020
16

---

Published as a conference paper at ICLR 2021
The original formulation by Anderson al. 2020) is the following:
Xk WAXK-1 + be
(Wa & x + Zμ: bb
We o T o 2₁.) + ( @
41 E 2
A'c.
(20)
Ix) Ux-1) o
(Xk. &t. 2R) E [k,, us)
x [R. 0g] x 0, 1]
The difference with respect to A& defined in equation (3) exclusively lies in the treatment of pre-activation bounds. While As explicitly employs generic lx: üg in the constraint set via MA. A implicitly sets 4. üg to the value dictated by interval propagation bounds (Gowal et al. 2018 via the constraints in = 0 and & 1 from the exponential family. In fact, setting 4 0 and Lx 1, we obtain the following two constraints:
X S &x ME @ (1-22) 24)
Xx S M Zk
where:
Mg
min
W Wi + by W4 @ Lx-1 + b
(21)
11-1 uk-1l
M =
max
W7 X + by = Wx @ Ük-1 + b
which correspond to the upper bounding ReLU constraints in Mr if we set i. + M ûN + Mr. While ic, úgg (potentially) computed solving an optimisation problem the entire network (problem 19h. the optimisation for M2 M4;"" involves only the layer before the current. Therefore, the constraunts in (2I) might be much looser than those in Ma.
In practice, the effect of [ u on the resulting set so significant that Mr might yield better bounds than A/e even on very small networks. We now provide simple example.
F.1 MOTIVATING EXAMPLE
Input-
Hidden 1
Hidden 2
Output -
(1)
( 1)
Xy [0]
x1 [0]
x2 [0]
(1)
-1
2)
-2
X3 [0)
(-1)
(2)
(-1)
+1
(1)
X [1]
[1
x2
Figure 5: Example network architecture in which Ma Ale, with pro-activation bounds computed with Ca Mg. For the bold nodes (the two hidden layers) ReLU activation follows the linear function. The numbers between parentheses indicate multiplicative weights, the others additive biases (if any).
Figure 5 illustrates the network architecture. The size of the network is the minimal required to reproduce the phenomenon. M₁ and A₁ coincide for single-neuron layers (Anderson et al. 2020), and i. = M û, = M on the first hidden layer (hence, second layer is needed). Let us write the example network as (not yet relaxed, as in problem (T) optimization problem for the lower bound on the output node X3.
1a arg min
[2
x2
(22a)
s.t.
Xy
[11, 1]2
(22b)
17

---

Published as a conference paper at ICLR 2021
81
Xo +
X = max(0, &1)
(22c)
-1
- 2
-2
82
X1 +
X2 max(0,
(22d)
- 1
X3 = 12 - 1 X2
(22e)
Let us compute pre-activation bounds with C₁ Ma (see problem (9)). For this network, the final output lower bound tighter if the employed relaxation is MM rather than A& (hence, in this case, Mr C AQ)). Specifically: I3.A4 =1 11.7777, 13.M. 1.2273. In fact:
In order to compute 1 and u1 the post-activation bounds of the first-layer, it suffices to solve a box-constrained linear program for 1 and ú, which at this layer coincide with interval propagation bounds, and to clip them to be non-negative. This yields = =0 oj",
u1 = [1 3/=
Computing M. [1] ma m.] 1-2 x|x = we assuming that x1 = [0]
and x₁ u1 [1]. These two assignments are in practice conflicting, they imply differ-
ent values for Xg. Specifically, x1 [ = uj [1] requires X = [u][0] [o[1]] = -1]. but this would also imply x, [0] u, [0]. yielding {2[1] 1 3.
Therefore, explicitly solving LP relaxation of the network for the value of ú2 [1] will tighten the bound. Using MR. the LP for this intermediate pre-activation bound is:
üg [1]
arg min
[-2 1]x1
(23a)
x,8,
s.t.
Xu
[-1,
1]
21
1], R2
(23b)
1 -1
-
8 =
Xo +
(23c)
-1
2 &1
(23d)
[
X1 ú ④ 21 =
② 21
(23e)
3
-33
XI < 8 - @1 - 21) = =
-1
④ (1 (())
21)
(23f)
Yielding ú2[1] = 2.25 < M2 [1]. An analogous reasoning holds for M₂ [1] and i [1].
In Mr. we therefore added the following two constraints:
x < &2 [1] i [1](1 (1)
(24)
x < [1]]]
that in A'e correspond to the weaker:
x [ M2 [1](1 22(1))
(25)
x [1] M₁ [1]
As the last layer weight corresponding to x2 x²[1] negative (W3(0, 1] -1), these con- straints are going to influence the computation of 1g.
In fact, the constraints in 24 both active when optimizing for 1s, whereas their counterparts for 1s.A4 in 25 are not. The only active upper constraint at neuron x2[1] for the Anderson relaxation [2[1] < [1[1], corresponding to the constraint from Al, with [1 j0 1]. Evidently, its effect is not sufficient
cantte-ballace the effect of the tighter constraints (24) for -2[1,-] = [1 1] and Ig 11 = [o Oj, yielding a weaker lower
bound for the network output.
18

---

Published as a conference paper at ICLR 2021
F.2 DERIVATION OF Ae
Having motivated an explicit pre-activation bounds treatment for the relaxation by Anderson et al.
2020 now extend the original proof for A's (equation (20)) to obtain formulation A& (as defined in equation 3 For simplicity, we will operate on single neuron xx
A self-contained way derive Ale is by applying Fourier-Motzkin elimination standard MIP formulation referred to as the multiple choice formulation (Anderson et al. 2019), which is defined
follows:
(Xg- 1, Xx [x (x
[xx-1xx + -x,,,,
x 0 w
[b] (1xx
x
w x
bx
Sk
(26)
1-11
Zk x
su-11 [][[ (x)
< x u 1 [7
[x[0] E [0, 1]
Where wi,k denotes the i-th row of Wk. and x and x are copies of the previous layer vari- ables. Applying (26 the entire neural network results quadratic number of variables (relative
to the number of neurons). The formulation can be obtained from well-known techniques from the MIP literature (Jeroslow 1987) (it the union of the two polyhedra for passing and blocking ReLU, operating in the space of x)- Anderson et al. (2019) show that A = - 2z (S)).
If pre-activation bounds le, ûR (computed as described in section are available, we can naturally add them 266) as follows:
(Xc (x xx 2x E Sk
}
1 [1 2x(1)) I / x k-1 + by [1] [+) S úg (xl1
Sk.
(27)
1s O zx [i] I / + by [ilzz ûR
[+]
We now prove that this formulation yields Ag when projecting out the copies of the activations.
Proposition. Sets Sy. from equation 27 and Av. from problem are equivalent, in the sense that A = Proj P (Sx).
Proof. In order to prove the equivalence, will rely on Fourier-Motzkin elimination as in the original Anderson relaxation proof (Anderson et al. 2019). Going along the lines of the original proof, we start from (26) and eliminate x x and x| exploiting the equalities. We then re-write all the inequalities as upper or lower bounds on [0] in order to eliminate this variable. As Anderson et al. (2019), we assume w; > 0. The proof generalizes by using L and Ú for w. [0] < 0, whereas the coefficient is the variable easily eliminated. We get the following system:
x² x2_110
W w
Σ 121
Wisk Libl L1] bulilzn x
(28a)
wi,k [0]
x2 10] <
Σ
Win Winnbill--l10
be [+] 2n 2x/(i))
(28b)
wis [0]
J>1
10_ <
w W
Σ
Wi,k 1 by buli|zz
(28c)
wi.k [0]
3>1
4 [0](1 zu [?])
< x < u [0 (x(1))
(28d)
[0] < xx-1(0 [0] - 1 2z
(28c)
[0 x [0] - 10
(28f)
19

---

Published as a conference paper at ICLR 2021
x [0] S
Wi
Σ
Wi.h
(ba( 1x
(28g)
Wi,k [0]
1>1
x0-1
w Wi
Σ
wi,h 111x1=1x1100 [bxli] u. [3]))=
(28h)
wi.x [0]
j21
1
(0)
(1x|i) bx(i))(1 =x(1l) Σ WA. [
(281)
wi.k [0]
- 10]
0
1
(ú bb(1)(( zx(il)
Σ s biix -bl
(28))
w [0]
where only inequalities (28g to are not present in the original proof. We therefore focus on the part of the Fourier-Motzkin elimination that deals with them, and invite the reader to refer to Anderson et al. (2019) for the others. The combination of these inequalities yields trivial
constraints. For instance:
281
28g
1 S 7 + bu(d] [
(29)
which holds by the definition of pre-activation bounds.
Let us recall that x 2 0 and x 2 [t. the latter constraint resulting from (8a
28b
Then, it can be easily verified that the only combinations of interest (i.e., those that do not result in constraints that are obvious by definition or are implied by other constraints) are those containing the equality (28a). In particular, combining inequalities (28g) to (28) with inequalities (28d) to (28)) generates constraints that are (after algebraic manipulations) superfluous with respect to those (30). We are now ready to show the system resulting from the elimination:
X [ Z 0
(30a)
[
2 & [i]
(30b)
x
S Wi,A Wi; 00x=-100
zx(i)) Σ
wrbbbx 151 + bxllzz]
(30c)
xx(i] < u (0) +
wi w [5] + by [
(30d)
3>1
[ 2 (0x
wis (0) ug [0](1 (x) + Σ wo,+[>x2 151 + by [1
(30e)
X+|i] 2 wi WH [0] +
w.bbixbk- 16] by [i]
(30f)
1>1
1=- [0] S [i] [
(30g)
x 2 [ [i]
(30h)
x < á [ [i]
(30i)
xx]i] < [+] - [x] - (
(30j)
[ 2 [ [1] zz ii)
(30k)
Constraints from 30a to (30g) are those resulting from the original derivation of A'. (see (Anderson
al. 2019)). The others result from the inclusion of pre-activation bounds in (27). Of these, (30h)
implied by (30a) if 1, [i]
and by the definition of pre-activation bounds (together with (30b)
i > O. Analogously.
(3K is implied by (30b) ú
[i] 2 0 and by (30a) otherwise.
By noting that no auxiliary variable is left in (30i) and in (30j). we can conclude that these will not be affected by the remaining part of the elimination procedure. Therefore, the rest of the proof
20

---

Published as a conference paper at ICLR 2021
(the elimination of x2-111 [1].11111-111
12),
proceeds as in (Anderson et al. 2019), leading Ak,i-
Repeating the proof for each neuron at layer k. we get A
Pro) Na. (Sk).
G MASKED FORWARD AND BACKWARD PASSES
Crucial to the practical efficiency of our solvers is to represent the various operations standard forward/backward passes over a neural network. This way, we can leverage the engineering efforts behind popular deep learning frameworks such as PyTorch (Paszke al.) 2017). While this can be trivially done for the Big-M solver (appendix [). the Active Set method (3..)) requires a spe- cialised operator that we call "masked" forward/backward pass. Here, we provide the details to our implementation.
The masked forward and backward passes respectively take the following forms: (W& O I&) Xt.
(W& @ T&)" Xk+1 and they are needed when dealing with the exponential family of constraints from the relaxation by Anderson al.) (2020). A crucial property of the operator is that the Tx mask may take on a different value for each input/output combination. While this is straightforward implement for fully connected layers, we need to be more careful when handling convolutional operators, which rely on re-applying the same weights (kernel) to many different parts of the image. naive solution is to convert convolutions into equivalent linear operators, but this has a high cost in terms of performance, as it involves much redundancy.
convolutional operator can be represented via matrix matrix multiplication if the input is un- folded and the filter is appropriately reshaped. The multiplication output can then be reshaped to the correct convolutional output shape. Given a filter w E Rexk xka an input x E Rú Xia xis and the convolutional output conv (x) = ceexxxxxx we need the following definitions:
17, I
{ R2 Xd - TR. x.xxx 1xd511 x
(31)
(nfoldy()) R Xigxii >> TR ky
foldw()) Tp Tgki2xxop00
R Xig
where the brackets simply reshape the vector from shape I to O while the braces sum the j-th dimension. unfold decomposes the input image into the (possibly overlapping) 0203 blocks the sliding kernel operates on, taking padding and striding into account. fold brings the output of unfold to the original input space.
Let us define the following reshaped versions of the filter and the convolutional output:
Wi nt
[w] geck x 82
YR y x xes, Rexeges
The standard forward/backward convolution (neglecting the convolutional bias, which can added at the end of the forward pass) then be written
conv (x) [Rn unfold unf (x)]Rexxza pex pex Xxy
(32)
back_conv,,()) fold (WT yR).
We need to mask the convolution with a different mask for each input-output pair. This means employing mask I E 12x Therefore, assuming vectors are broadcast to the correct output ". we can write the masked forward and backward passes as:
conv = [[WR @ I O unfold unfold((x))22eer) (xlop bwck.convw...... (y) ([W,, @ T yn))).
(33)
Owing the avoided redundancy with respect to the equivalent linear operation (e.g., copying of the kernel matrix, zero-padding in the linear weight matrix), this implementation of the masked for-
'if we want to perform an element- wise product @ b between a E xdxxxx and b E dddxds, the operation is implicitly performed as a @ b', where b' E Rdixd2xds dy is an extended version of b obtained by copying along the missing dimension.
21

---

Published as a conference paper at ICLR 2021
ward/backward pass reduces both the memory footprint and the number of floating point operations (FLOPs) associated to the passes computations by a factor (f12243) ( k2). In practice, this ratio
might be significant: the incomplete verification networks [5.1 1 it ranges from 16 to 64.
H STRATIFIED BOUNDING FOR BRANCH AND BOUND
As seen in the complete verification results in Figure 3 (.2), the use of a tighter bounding algorithm results in the verification of larger number of properties. In general, tighter bounds come at larger computational cost, which might negatively affect performance where a small number of domain splits with loose bounds suffices to verify the property (hence, the tightness is not needed). As general complete verification system needs to be efficient priori, both easy and hard properties, we employ a stratified bounding system for use within Branch and Bound-like complete verification algorithms.
"easy" verification properties,
We design a set of heuristics to determine whether the children (i.e., the subproblems arising after splitting it) of given subproblem will require tight bounds not. The graph describing parent- child relations between sub- sub-problems is referred to as the BaB search tree. Given a problem we individually mark it as difficult (with all its future offspring) if it meets all the conditions in the following set:
The lower bound on the minimum over the subproblem, PLB. should be relatively "close" to the decision threshold (0, in the standard form by |unell et al. (2018)). That is, PLB S CLB- We argue that tighter hounds are worth computing only if they are likely to result in a crossing of the decision threshold, thus cutting away a part of the BaB search tree.
The depth Pdepth in the BaB search tree should not be below certain threshold Cdepth. This is a way to ensure a certain number of splits performed before adopting tighter bounds.
The running average of the lower bound improvement from parent child Pimpr should fall below threshold Compr- This reflects the idea that if splitting seems to be effective on a given section the BaB tree, is perhaps wiser to invest computational budget more splits (with cheaper bounding computations, more splits can be performed in given time) than tighter bounds. Empirically, this is the criterion with the largest effect on performance.
Additionally, we not employ tighter bounds unless the verification problem itself (in addition to the individual subproblems) is marked as "hard". We do so when the size of the un-pruned domains reaches certain threshold Chord-
The set of criteria above needs to address the following decision trade-off: one should postpone the use of tighter bounds until sure the difficulty of the verification task requires it. At the same time, if the switch to tighter bounds is excessively delayed, one might end up computing tight bounds on a large number of subproblems (as sections of the BaB tree were not pruned in time), hurting performance for harder tasks. In practice,
the criteria as follows for all experiments: Cdepth 0, CL,B = 0.5, Cimpr = 10 Chard = 200. As seen in Figure 3 this resulted a reasonable trade-off between losing efficiency the harder properties (wide model) and gaining on the easier ones (base and deep models).
22

---

Published as a conference paper at ICLR 2021
I EXPERIMENTAL APPENDIX
We conclude the appendix by presenting supplementary experiments and additional details with respect to the presentation in the main paper.
1.1 EXPERIMENTAL SETTING, HYPPER-PARMEER
All the experiments and bounding computations (including intermediate bounds) were run on single Nvidia Titan Xp GPU, except Gurobi-based methods and "Active Set CPU". These were run on 17-6850K CPUs, utillising 4 cores for the incomplete verification experiments, and 6 cores for the more demanding complete verification experiments. The experiments were run under Ubuntu 16.04.2 LTS. Complete verification results for ERAN, the method by /Singh et al. (2020), from the recent VNN-COMP competition (VNN-COMP 2020). These were executed by [Singh] al. (2020) on a 2.6 GHz Intel Xcon CPU E5-2690 with 512 GB of main memory, utilising 14 cores.
taken
Gused methods make use of LP incrementalism (warm-starting) when possible. In the ex- periments of [.. where each image involves the computation of 9 different output upper bounds, we warm-start each LP from the LP the previous neuron. For "Gurobi 1 cut", which involves two LPs per neuron, we first solve all Big-M LPs, then proceed with the LPs containing a single cut.
Hyper-parameter tuning for incomplete verification was done on a small subset of the CIFAR-10 test set, on the SGDD-trained network from Figures I 2 BDD+ is run with the hyper-parameters found by Bunel et al. (2020a) on the same datasets, for both incomplete and complete verification. For all supergradienttbasse methods (Big-M, Active Set), employed the Adam update rule (Kingma & 2015), which showed stronger empirical convergence. For Big-M, replicating the findings by Bunel et al. (2020a) their supergradient method, we linearly decrease the step size from 10-3
Ba
to 10 Active Set is initialized with 500 Big-M iterations, after which the step size is reset and linearly scaled from 10-3 3 to 10-6, We found the addition of variables the active set to be effective
before convergence: we add variables every 450 iterations, without re-scaling the step size again. Every addition consists of new variables (see pseudo-code in appendix D). which was found to be good compromise between fast bound improvement and computational cost. On complete verification, we re-employed the same hyper-parameter for both Big-M and Active Set, except the number of iterations. For Big-M, this was tuned to employ roughly the same time per bounding computation as BDD+.
As the complete verification problem is formulated as a minimisation (as in problem (1), in Branch and Bound we need a lower and an upper bound to the minimum of each sub-problem. The lower bound is computed by running the bounding algorithm, while the upper bound on the minimum is obtained by evaluating the network at the last primal output by the bounding algorithm in the lower bound computation (rumning the bounding algorithm to get upper bound would result in a much looser bound, as it would imply having an upper bound on a version of problem ( with maximisation instead of minimisation).
1.2 DATASET DETAILS
We now provide further details on the employed datasets. For incomplete verification. the architec- ture was introduced by WWong& & Kolter (2018) and re-employed by Buneel et al. (2020a). It corre- sponds the "Wide" complete verification architecture, found in Table 2 Due to the additional computational cost of bounds obtained via the tighter relaxation (3. we restricted the experiments to the first 3450 CIFAR-IO test set images for the experiments on the SGD-trained network (Figures 2 and to the first 4513 images for the network trained via the method by Madry et al. (2018)
(F (Figures 6
For complete verification, we employed a subset of the adversarial robustness dataset presented by Lu & Kumar (2020) and used by Bunel et al.) (2020a). where the set of properties per network has been restricted to 100. The dataset provides, for a subset of the CIFAR-10 test set, a verification radius Ever defining the small region over which to look for adversarial examples (input points for which the output of the network is not the correct class) and (randomly sampled) non-correct class to verify against. The verification problem is formulated as the search for an adversarial example, carried out by minimizing the difference between the ground truth logit and the target logit. If the minimum is positive, we have not succeeded in finding a counter-exant exampe,, and the network is ro- bust. The ever radius was tuned to meet certain "problem difficulty" via binary search, employing
23

---

Published as a conference paper at ICLR 2021
Network Name
No. of Properties
Network Architecture
Conv2d(3,8,4) stride=2, padding=1)
Conv2d(8, stride=2, padding=1)
BASE Model
100
linear layer of 100 hidden units
linear layer of 10 hidden units
(Total ReLU activation units: 3172)
Conv2d(3, 16,4, stride=2, padding=1) Conv2d 16,22,4 stride=2, paddlig==1)
WIDE
100
linear layer of 100 hidden units
hidden units (Total ReLU activation units: 6244)
linear layer
Conv2d(3,88 stride=2, padding=1) Conv2d stride=1, packing==)
Conv2d(8,8, stride=1, padding=1)
DEEP
Conv2d(8,8,2 stride=2, padling=1)
linear layer of 100 hidden units linear layer 10 hidden units (Total ReLU activation units: 6756)
Table each complete verification experiment, network architecture used and number of verifica- tion properties tested, a subset of the dataset by Lu & Kumar (2020). Each layer but the last followed by ReLUU activation function.
Gurobi-based bounding algorithm (Lu & Kumar 2020). The networks robust all the prop- erties we employed. Three different network architectures of different sizes are used. A "Base" network with 3172 ReLU activations, and two networks with roughly twice as many activations: "Deep" network, and a "Wide" network. Details can be found in Table 2
1.3 ADVERSARIALLY- TRAINED INCOMPLETE VERIFICATION
Figure 6: Upper plot: distribution runtime in seconds. Lower plot: difference with the bounds obtained by Gurobi with from A& per neuron; higher is better. Results for the network adversarially trained with the method by Madry al. (218)), from Bunel et al. 2020aa). The width at given value represents the proportion of problems for which this is the result.
10
Lo1
Inproverrent eereeeentti fromm cut
-0.5
-1.0
Gurobi Panet
Gursbi
Active B0 sheps
1 cut
400 steps
850 steps
600 steps
1050 steps
1650 steps
600 steps
In addition to the tecture, trained with adversarial training method by Madry et al. (2018) for robustness to per- turbations of étrain 8/255. Each adversarial sample for the training was obtained using 50 steps of projected gradient descent. For this network, we measure the robustness margin to perturbations with Ever 12/255.
network in $.11 we now present results relative to the same archi-
Figures 6 7 confirm most of the observations carried out for the SGD-trained network $5.1 with fewer variability around the bounds returned by Gurobi cut. Big-M is competitive with BDD+, and switching to Active Set after 500 iterations results much better bounds the same time. Increasing the computational budget for Active Set still results in better bounds than Gurobi cut in
24

---

Published as a conference paper at ICLR 2021
Tiiing
Gap
Timtts 2m
epprovmmett
22-1
Sut 30000
0.4
03
B
15
00
0.0
64
depo
(a) Comparison runtime (left) gap to Gurobi Planct bounds. (b) Compariscn rontime (left) improvement from the Gurobi For the latter. lower better.
Planet bounds. For the latter, higher better.
Pointwise comparison for subset the methods the data presented in Figure Darker colour shades mean higher point density (on a logarithmic scale). The oblique dotted line corresponds the equality.
Figure 7:
fraction of its running time, even though the performance gap is on average smaller than on the SGD-trained network.
1.4 SENSITIVITY TO SELECTION CRITERION AND FREQUENCY
In section 3.2. we deseribe how to iteratively modify B, the active set of dual variables on which
our Active Set solver operates. In short, Active Set adds the variables corresponding to the output
Figure 8: Upper plot: distribution of runtime in seconds. Lower plot: difference with the bounds obtained by Gurobi with a cut from As per neuron; higher better. Results for the SGD-trained network from Bunel ct al. (2020a). The width at given value represents the proportion of problems for which this is the result. Sensitivity of Active Set to selection criterion (see 3.2).
3x10
Tima
107
Method
0.5
frem
0.2
Ran, Selection
Selection 1050
Active Sut,
Pianet
1050
1650
Figure 9: Upper plot: distribution of runtime in seconds. Lower plot: difference with the bounds obtained by Gurobi with a cut from As per neuron; higher is better. Results for the SGD-trained network from Bunel al. (2020a). The width given value represents the proportion of problems for which this is the result. Sensitivity of Active Set variable addition frequency w, with the selection criterion presented in 3.2
3x0
2x000
10a
Method
improvement from
0.0
.02
Active Set
A5, o- 300 steps
A5. 600 steps
Active Set
AS um 300
A5 uu 600
Active Set
AS. - 300
A6, 600 stops
steps
25

---

Published as a conference paper at ICLR 2021
of oracle (4) invoked at the primal minimiser of LB(x, z, a, (B), at fixed frequeney We now investigate the empirical sensitivity of Active Set to both the selection criterion and the frequency of addition.
We test against Ran. Selection, a version of Active Set for which the variables to add are selected random by uniformly sampling from the binary Ja masks. As expected, Figure g shows that good selection criterion key to the efficiency of Active Set. In fact, random variable selection only marginally improves upon the Planet relaxation bounds, whereas the improvement becomes significant with our selection criterion from :3.2
In addition, we investigate the sensitivity of Active Set (AS) to variable addition frequency w. In order to do so, we cap the maximum number of cuts at 7 for all runs, and vary w while keeping the time budget fixed (we test on three different time budgets). Figure 9 compares the results for = 450 (Active Set), which were presented in 15.1 with the bounds obtained by setting 600). Our solver proves to be relatively
300
and w 600 (respectively AS w = 300 and AS robust across all the three budgets, with the difference in obtained bounds decreasing with the number of iterations. Moreover, early cut addition tends to yield better bounds in the same time, suggesting that our selection criterion is effective before subproblem convergence.
1.5 MNIST INCOMPLETE VERIFICATION
We conclude the experimental appendix by presenting incomplete verification results (the experi- mental setup mirrors the one employed in section 5.1 and appendix I.3) on the MNIST dataset (Le-
Cun et al. 1998).
We report results on the "wide" MNIST network from Lu Kumar (2020), whose architecture identical to the "wide" network in Table 2 except for the first layer, which has only one input channel to reflect the MNIST specification (the total number of ReLU activation units 4804). As those employed for the complete verification experiments (55.) and differently from the incomplete
10a
Time
Method
Improvement from cut
-0.5
-1.0
Gurobi
Gurobi
800 steps
M
Active Sut steps
Active Set
Active Set
Actie Set
Active Set steps
Figure 10: Upper plot: distribution of runtime in seconds. Lower plot: difference with the bounds obtained by Gurobi with a cut from An neuron; higher is better. MNIST results for a network adversarially trained with the method by [Wong Kolter (208)), from Lu & Kumar 2020). The width a given value represents the proportion of problems for which this is result.
Tinkng 1h
Instevament from Hlant
143
Trring [mna
Instevetment from Plant
:2
3.3
S6.
02
60
6.3
02
CUrob L 001
Gurcal out
Guxbb 1
Figure 11: Pointwise comparison for a subset the methods on the data presented in Figure 10 Comparison of runtime (left) and improvement from the Gurobi Planet bounds. For the latter. higher is better. Darker colour shades mean higher point density (on a logarithmic scale). The oblique dotted line corresponds to the equality.
26

---

Published as a conference paper at ICLR 2021
verification experiments in section 5.1 and appendix 1.3 the network was adversarially trained with the method by Wong & Kolter (2)18). We compute the robustness margin to Ever 0.15 on the first 2960 images of the MNIST test set. All hyper-parameters are kept to the values employed for the CIFAR-10 10 networks, except the Big-M step size, which was linearly decreased from 10-1 to 10-3, and the weight of the proximal terms for BDD+, which was linearly increased from to 50.
As seen on the CIFAR-10 10 networks, Figures TO T1 show that Active Set yields comparable or better bounds than Gurobi cut in less average runtime. However, more iterations are required to reach the same relative bound improvement over Gurobi 1 cut (2500 as opposed to 600 in Figures [ O Furthermore, the smaller average gap between the bounds of Gurobi Planet and Gurobi cut (especially with respect to Figure I suggests that the relaxation by Anderson et al. (2020) is less effective on this MNIST benchmark.
27


Reviews:
[{"review_id": "t9PpC4eAqMt", "reviewer": "AnonReviewer3", "summary": "", "questions": "", "limitations": "", "rating": 7, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "The paper proposes a dual approach for solving the exponentially-sized linear programs that arise from the relaxation of the single ReLUs by Anderson et al. (2020). The algorithm is demonstrably faster than existing methods in experiments resembling those from Bunel et al. (2020). The idea of applying a dual active set method is arguably quite straightforward from an optimization perspective. What makes the contribution in this paper strong is \n\n(1) A good implementation in multiple ways. The performance of active set methods highly depends on the how the sets are maintained (section 3.2). Also, making the algorithm run efficiently in practice requires taking advantage of GPUs and the neural network structure (section 3.3, G). I believe that without doing these well, it is likely that such an approach would be significantly slower than existing methods. I look forward to seeing the authors code released.\n\n(2) Providing the proper context for the work -- The paper builds upon prior contributions from Anderson et al. and other works that study the dual aspect of the problem. The authors do a good job explaining how their work fits in and how it compares to other approaches. This both provides a nice theoretical grounding for this work and is also useful pedagogically.\n\n(3) Detailed experiments on reasonably sized neural networks and providing hyperparameters.\n\nTwo things I would like the authors to address:\n\n(1) Experimental results - CPU-only\n\nHow does an all CPU version of the algorithm compare to Gurobi? This would be a better apples-to-apples comparison since commercial solvers do not make use GPUs. I would understand if a CPU-only implementation is slower since Gurobi is highly-tuned.\n\n(2) Comparison with Tjandraatmadja et al. (2020)\n\nI would like the authors to contrast their approach more against the the one by Tjandraatmadja et al. The formulations used there are different (and the focus of that paper is different), they do describe a cutting plane approach, which in a sense can be viewed as an approach that incrementally increases the number of dual variables. The authors can either do so here in the response, or if they think it would benefit the paper add it in.\n\n-------------------------\nUpdate after author response:\nThanks to the authors for addressing my questions!", "extra_scores": {}, "extra_text": {}}, {"review_id": "cB6EuUc3gS", "reviewer": "AnonReviewer1", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "1) Summary\nThe manuscript proposes an optimisation scheme to compute bounds on the output of piecewise linear feedforward networks. Improving upon the planet relaxation, the authors provide a tighter relaxation scheme with exponentially man constraints that they tackle by a dual solver using active sets. Some experiments are conducted in order to benchmark the relaxation with others.\n\n2) Strengths\n+ The paper is mostly well written.\n+ The evaluation seems pretty exhaustive as it comes to competing optimisation schemes.\n\n3) Concerns\n- The practical utility for deriving a proper bound for a relevant architecture is still limited.\n- It does not seem that the paper makes code and data available to the public.\n\n4) Remarks/Questions\n  a) Section 2: The optimization problem in Equation (1) needs to be sharpened. On the one hand, $x^hat_n$ needs to be a scalar quantity in order to have a meaningful objective and on the other hand, the only free variable is $x_0$ as all the other variables depend on it. The first part also applies to Equations (2) and (3). Also, the solution is not necessarily unique in the ReLu case if many neurons are silent.\n  b) References: capitalization not correct e.g. \"smt\"\n  c) Please provide a more concise statement what you require from the oracle. Do you require that $min_{x in C} a'*x$ can be computed in $O(1)$?\n  d) What restricts the approach to piecewise linear functions?", "extra_scores": {}, "extra_text": {}}, {"review_id": "Ue-TsCMawlK", "reviewer": "AnonReviewer2", "summary": "", "questions": "", "limitations": "", "rating": 7, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "The authors present a custom solver for verifying properties of neural networks (such as robustness properties). Prior work for neural network verification relies on generating bounds by solving convex relaxations (\"convex barrier\"). The authors describe a sparse dual solver for a new relaxation which is tighter (but has higher computational complexity). The solver is represented (for the most part) as standard operations built into pytorch, and so it can be easily run on GPUs (they do require a specialized operator to support masked forward/backward passes, and they describe how this is done efficiently for convolutional networks). The solver involves repeatedly solving modified instances of a problem, where only a small active set of dual variables (instead of exponentially many) is considered at each step.\n\nExperimental results are promising in that it outperforms generic solvers in terms of both the bounds achieved and the time taken to do so. This does seem to be a promising approach.\n\n", "extra_scores": {}, "extra_text": {}}, {"review_id": "P6Dpu81SbrS", "reviewer": "AnonReviewer6", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This paper describes a dual solver for neural network verification. In particular, the authors consider a linear relaxation of neural networks with relu activations and develop an active set based method. Numerical comparisons show that the proposed method provides speed-ups in verifying deep networks.\n\nAlthough the authors present promising empirical results, a precise and rigorous analysis and of the active set strategy is lacking.  It would be great if the authors can clarify several points raised below.\n\n\n1. Does the active set approach provide any guarantees on the tightness of the solution? It looks like exponentially many optimization variables in eq 6 are initialized at zero, which provides a lower-bound on eq 2. However, it's not clear if the produced lower-bounds are effective.\n\n2. Does the greedy active set extension strategy converge to an optimal solution of the relaxation? Is the algorithm sensitive to the selection criterion and frequency?\n\n3. The authors employ projected gradient with Adam to maximize the dual function d(\\alpha,\\beta). Does this approach provably converge to the solution of eq 9? It would be nice to describe the performance of other optimizers, e.g., plain SGD.\n\n4. In eq (1a), is \\hat x_n a scalar? It would be better to specify that n_n=1.\n\n5. It looks like the integer constraints z\\in \\{0,1\\} are missing in eq 2? \n\n6. Last paragraph on page 2, by an optimal solution of the problem 2 are you referring to the linear relaxation?", "extra_scores": {}, "extra_text": {}}, {"review_id": "E7HW2HnN12_", "reviewer": "AnonReviewer5", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This paper presents an algorithm for verification of neural networks with ReLU activations. Essentially, it takes the tightened ReLU relaxation of (Anderson et al. 2020), builds its Lagrangian dual, and then applies a column generation scheme to accommodate the explosion of decision variables. The authors present computational results showing that, due to the amenability of the methods to GPU accelaration, they can produce stronger verification bounds than comparable methods working with weaker relaxations in a modest amount more time. Moreover, they show that the method can be successfully embedded in a branch-and-bound-like framework for exact verification.\n\nI like the paper and think it makes a good contribution to the literature. The paper builds heavily on (Anderson 2020), but the algorithm approach is very different and the authors clearly had to do some work (e.g. rederiving the dual to ensure efficient inner problem solves, making convolutions+masking work on the GPU) to make everything work out. The computational results also seem compelling, though I have some potential concerns about the comparison being made.\n\nMy main concern is the use of \"Gurobi 1 cut\" as the baseline for comparison. Given that there is a one-to-one mapping between cuts in the primal formulation (Gurobi 1 cut) and variables in the dual formulation (the new approach), I am curious why the authors did not choose symmetric generation schemes for the two. Would the solve times be significantly lower if only one cut per layer is added (as in ActiveSet), instead of one per neuron? If so, what benefit do multiple iterations of cut generation offer? Is LP incrementalism or warm-starting used, or is the second LP solved from scratch? Even with these changes, I would imagine that the ActiveSet method still runs (much?) faster than the primal approach, but it's quite possible that the bound improvement would shrink.\n\nMinor comments:\n* p1: The phrasing \"The main bottleneck of the above approach\" is ambigious (which approach?). If the approach includes (ii), then the bottleneck will be the enumeration tree, not the convex subproblems.\n* p2: The \"which is linearly-sized\" reads like it applies to the optimal solution of (2), not the formulation (2) itself.\n* p18: I think there are some extra primes in the text of Appendix F.2.", "extra_scores": {}, "extra_text": {}}, {"review_id": "18XxitswLpL", "reviewer": "AnonReviewer4", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 1, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "The paper focuses on the notion of dual solvers, which are solvers that have the ability to prove as well as disprove that a specified property of neural networks holds over all inputs in a specified domain. Unsurprisingly, the bottle neck for this approach is the computation of lower bounds, which can be translated into a non-linear (and NP-hard) optimization problem. Dual solvers typically instead attempt to solve a linear formulation of the problem, and recently Anderson et al. introduced a new MIP formulation for the problem which is more accurate than previously used formulations at the cost of potentially requiring a significantly larger number of constraints.\n\nIn this work, the authors introduce a dual solver called \"Active Set\" that is designed to solve Anderson's tighter formulation, with the aim of achieving increased accuracy while keeping the computational costs low. The solver seems correct, even though I couldn't verify the details; generally speaking, it operates by starting with the easier but less accurate formulation that is classically used, and making it more precise by the gradual introduction of (an \"active set\" of) variables from a modified version of Anderson et al's more accurate formulation. Experimental evaluations show that the new solver can be used both for incomplete verification and complete verification; in both cases, it is possible to use Active Set to obtain an implementation that can either achieve higher accuracy (for incomplete verification) or avoid a timeout for a greater percentage of instances (for complete verification). \n\nThe paper is well written. The preliminaries and writeup in general are clearly aimed at experts with some background knowledge in the given area (which prevented me from verifying the technical details), but I believe that it fits well within the scope of ICLR. The main contribution is the design of a new solver and its experimental evaluation; the theoretical contribution is negligible. It is also easy to imagine that there could be many alternate ways one could use Anderson et al.'s recent formulation to improve the state of the art - in that sense, the results achieved by the authors are not surprising. But that does not mean that there were no challenges left to overcome, and the experimental evaluations seem to be reasonable. So I think that the paper's contribution is sufficient to warrant a presentation at ICLR.\n\nQuestion:\n-Are there other \"standard\" datasets that could be considered instead of CIFAR-10, and would it be difficult to also use these for experimental evaluations?\n\nMinor remarks:\n-page 2: \"if more compute budget is...\"  should be \"if a larger computational budget is...\"\n\nPost-rebuttal comment:\nI acknowledge having read the authors' response and I have also glanced over the updated version of the paper.", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: Dear authors,

as you have noticed this paper was not easy to review. I have hence invited 2 additional reviewers which I strongly respect and are very knowledgeable. After carefully reading the paper myself, I have to agree with one of the reviewers who said "... it [your paper] makes a good contribution to the literature ....". To be honest, we were working in my group on a very similar approach but did not manage to finish it (and I know how hard it is).

To conclude, when preparing to the final version, please try to go over the reviews, I am sure they can make your paper even stronger :) 

Decision: Accept (Poster)