DEEPSPLIT: An Efficient Splitting Method for Neural Network Verification via Indirect Effect Analysis

1 Introduction

AI-based, safety-critical applications, including autonomous driving, require reliable systems that are amenable to analysis leading to their certification. At present even the most complex neural classifiers suffer from misclassification rates that render the technology potentially unsafe. It therefore remains of importance to develop methods to formally verify the correctness of neural networks and identify unwanted features.

The area of formal verification of neural networks [Liu et al., 2019] is concerned with the development of methods to ascertain whether a neural model is correct with respect to a given specification. For example, in image classification, we might need to establish whether for a given image a classifier is robust with respect to noise [Szegedy et al., 2014] or geometric perturbations [Kouvaros and Lomuscio, 2018; Balunovic et al., 2019]. In closed-loop systems with neural controllers, it is of interest to establish whether the system may violate some safety constraints expressed as temporal specifications [Akintunde et al., 2018; Akintunde et al., 2019; Tran et al., 2020b; Akintunde et al., 2020].

Given their noteworthy application domains, particular attention has been devoted to verifying high-dimensional, ReLU-based, feed-forward neural networks; this is also the aim of this paper. While a number of increasingly scalable methods have been proposed recently, the present state-ofthe-art (SoA) cannot yet address the models used in practice. The key contribution in this paper is to define a procedure, within a symbolic interval propagation framework, for deciding which nodes to split based on the direct and indirect effect that the splitting has on the relaxations in subsequent nodes in the network. As we show, this leads to 1-2 orders of magnitude speedup over current SoA verifiers.

Related Work

Verification methods for neural networks can be divided into complete and incomplete methods. Complete methods are theoretically guaranteed to provide an answer to any verification query; incomplete methods may not reach a conclusion on a query but may be more scalable when a solution can be found. Both streams of work are considered important as they can complement each other.

Complete approaches can be further divided into Mixed Integer Linear Programming (MILP), reachable set and relaxation based approaches. MILP based approaches encode the verification problem into MILP solvers [Lomuscio and Maganti, 2017; Akintunde et al., 2018; Anderson et al., 2020; Botoeva et al., 2020; Singh et al., 2019b; Tjeng et al., 2019]. Reachable-set based approaches [Tran et al., 2020a; Tran et al., 2020b; Bak et al., 2020] verify properties by propagating exact symbolic representations of the reachable states through the network and reasoning over the resulting output states. Relaxation based methods [Henriksen and Lomuscio, 2020; Katz et al., 2019; Rubies-Royo et al., 2019; Wang et al., 2018a; Bunel et al., 2020; Singh et al., 2019a] compute a linear relaxation of the network and iteratively refine it by splitting ReLU nodes and branching; an efficient relaxation can be obtained with Symbolic Interval Propagation (SIP) [Wang et al., 2018b; Wang et al., 2018a; Singh et al., 2019c; Singh et al., 2018; Henriksen and Lomuscio, 2020]. The work presented here also adopts the SIPbased approach, but differs from the SoA in that it (i) employs a novel procedure for splitting which is based on dependencies between ReLU relaxations, (ii) utilises a new succinct LP-encoding of the split-constraints and (iii) combines two different SIP algorithms to improve the linear relaxation. This results in significant speedups, as we report below.

[Lu and Kumar, 2020; Bak, 2020] also consider procedures for branching decisions. [Lu and Kumar, 2020] uses a graph-neural network for deciding on branching nodes; this approach was only shown to be able to verify relatively small networks with < 7000 ReLU nodes. In contrast we use intermediate results from the SIP phase, leading to a more scalable approach. The procedure in [Bak, 2020] starts from a nonrelaxed network and gradually increases the nodes to relax; in contrast, we start with a relaxed network and gradually decrease the abstraction through branching. In our experiments this results in a more scalable approach.

Node-dependency analysis is also present in [Botoeva et al., 2020; Rubies-Royo et al., 2019]. However, [Botoeva et al., 2020] determines binary relations between nodes, while we estimate fractional dependencies between the nodes' relaxations. Dependencies are considered in [Rubies-Royo et al., 2019] via the solution of a linear program encoding all the network's ReLUs; in common with other SIP-based approaches, we avoid encoding individual ReLUs for efficiency.

Incomplete verification algorithms [Singh et al., 2019c; Singh et al., 2018; Wang et al., 2018b] differ from complete approaches in that they reason over a relaxed, or abstract network and do not implement a refinement strategy which guarantees that all overestimation is removed in a finite number of steps and are therefore not directly comparable to complete approaches. However, we show later that SoA complete approaches significantly outperformed a SoA incomplete approach on several benchmarks.

2 Preliminaries

In this section we introduce feed-forward neural networks, symbolic interval propagation and some related concepts.

We assume a Feed-Forward Neural Network (FFNN) to be composed of an input layer, one or more hidden layers and an output layer [Goodfellow et al., 2016]. Each layer consists of multiple nodes and hidden layers are governed by non-linear activation functions σ : R → R. We here consider networks with the ReLU(z) = max(0, z) activation function only. A layer's output is calculated by applying the layer's activation function element-wise to its input vector where the input vector is an affine transformation of the output from previous layers. So, a network with n input nodes and m output nodes is associated with a function f : Rn → Rm where f is the composition of the layer-by-layer transformations.

In this paper we consider verification problems defined by box constraints on the network's input and linear constraints on the output as defined in Definition 1. While much of the current body of work on neural network verification does not explicitly define the verification problem, the problem presented here is in line with the main focus of the field.

Definition 1. A verification problem is a tuple 〈f, ψin, ψout〉 where f : Rn → Rm is a FFNN with ReLU activation functions for all hidden nodes, ψin = ∧mk=0 lk ≤ xk ≤ uk is the input constraint, ψout = ∧i (∨j aᵀi,j f(x) + bi,j < 0) is the output constraint, lk, uk, bi,j ∈ R and ai,j ∈ Rm. If ψout is satisfied for all x that satisfy ψin, then we say that the f is safe for the given constraints; otherwise, f is unsafe. An x that satisfies ψin and violates ψout is a counterexample.

Verification problems are challenging to solve due to nonlinear activation functions. A scalable approach to solve verification problems involves considering a linear relaxation of the FFNN [Wang et al., 2018b; Wang et al., 2018a]. The relaxed FFNN can be obtained by combining Symbolic Interval Propagation (SIP) with two-constraint linear relaxations of the activation functions. In the rest of this paper we follow the presentation of two-constraint linear relaxations from [Henriksen and Lomuscio, 2020] summarised below.

Definition 2 (Two-constraint relaxation). Let zl, zu ∈ R be concrete lower and upper bounds on the input of an activation function σ : [zl, zu] → R, where [zl, zu] ⊆ R denotes the closed real interval between zl and zu. A two-constraint linear relaxation is a tuple 〈rl, ru〉, where rl : [zl, zu] →R is a lower bounding linear function such that rl(z) ≤ σ(z) and ru : [zl, zu] → R is an upper bounding linear function such that ru(z) ≥ σ(z) for all z ∈ [zl, zu].

For ReLU functions where the input is bounded by lower and upper bounds zl, zu ∈ R, the upper relaxation is normally defined by ru(z) = zu zu-zl (z -zl) [Wang et al., 2018a; Singh et al., 2019c]. The formulation of the lower relaxation depends on the SIP algorithm used. Some approaches use rl(z) = zu zu-zl z [Wang et al., 2018a]; others minimise the overestimation area by considering rl(z) = 0 iff zu < |zl| and rl(z) = z otherwise [Singh et al., 2019c].

Relaxations are used in conjunction with SIP to calculate bounds on the network's output. Relevant to this work are two distinct SIP algorithms: Error-based Symbolic Interval Propagation (ESIP) [Wang et al., 2018a; Henriksen and Lomuscio, 2020] and the algorithm here referred to as Reversed Symbolic Interval Propagation (RSIP) [Singh et al., 2019c].

The ESIP algorithm propagates bounding equations and concrete errors layer-by-layer through the network. Activation functions are handled by propagating the terms through the lower relaxation and adding a new concrete error initialised as ϵ = maxz ∈[zl,zu](ru(z) - rl(z)) for each node; these errors are then propagated through the network with the equations. In the following ϵl1,l2 n1,n2 denotes the error from a node in layer l1 and position n1 at a node in layer l2 and position n2 where l1 < l2; however, we sometimes use the simplified notation ϵn to mean the error from a node n when this does not cause confusion. For a node in layer l2 and position n2 with equation ql2 n2 (x), the concrete lower and upper bounds zl2 n2,low , zl2 n2,up (denoted as yn2,low and yn2,up at the output layer) are calculated as zl2 n2,low = minx ql2 n2 (x) + ∑l2-1 l1=0 ∑n1|ϵl1,l2 n1,n2 <0 ϵl1,l2 n1,n2 and zl2 n2,up = maxx ql2 n2 (x) + ∑l2-1 l1=0 ∑n1|ϵl1,l2 n1,n2 >0 ϵl1,l2 n1,n2 .

Different from ESIP, the RSIP algorithm calculates bounding equations by substituting variables layer-by-layer backwards through the network. For non-linear activation functions, the upper symbolic bounds are back-propagated by applying the upper relaxation on positive coefficients and lower relaxation for negative coefficients. Correspondingly, lower bounds are calculated by applying the lower relaxation to positive and upper relaxations to negative coefficients. RSIP often produces more succinct bounds than ESIP; however, in the next section we show that intermediate calculations from ESIP are particularly well suited to analyse the network.

The bounds from the symbolic interval propagation can be used to calculate bounds for the linear equations in the verification problem's output constraint ψout [Singh et al., 2018;

Singh et al., 2019c]. Complete verification approaches [Wang et al., 2018b; Wang et al., 2018a; Henriksen and Lomuscio, 2020] combine SIP with a branch and bound phase and linear programming to constrain the ReLU nodes' input to > 0 and < 0 in separate branches to achieve completeness. Since the bounds produced by both SIP algorithms are exact whenever all ReLU nodes operate in one of their linear domains, the verification problem can be solved by exploring at most 2N branches where N is the number of ReLU nodes in the network. In Section 4 we show how ESIP and RSIP can be combined into an efficient verification algorithm.

Due to the branch and bounds phases' exponential complexity, it is infeasible to branch on all ReLU nodes even for small networks. So, splitting the correct nodes is essential for an efficient implementation. In the next section we propose a novel method for identifying good splitting candidates thereby resulting in a an efficient verification procedure.

3 Splitting Score Functions

A key aspect differentiating all symbolic interval propagation methods lies in the choice of nodes to be split. Technically, splitting a ReLU node is realised by creating two separate branches in the resolution tree: in the first the node's input is constrained to negative values, in the second to positive ones. In the ESIP algorithm this operation has three distinct effects.

• Direct effect. The error cascading from the split node becomes 0 in both branches, thus improving the bounds of successor nodes.
• Indirect effect. Since the bounds for all nodes succeeding the split node improve, the corresponding linear relaxations are also improved which, in turn, tightens the bounds in all layers after those nodes.
• Propagation effect. After the split, the equations and errors are propagated through the ReLU function instead of the lower relaxation at the split node.

While the categorisation above is novel to this paper, the direct effect has previously been used to determine splitting candidates [Henriksen and Lomuscio, 2020] and the gradientbased methods from [Wang et al., 2018b; Wang et al., 2018a] are closely related to the direct and propagation effects. However, neither of the methods account for the indirect effect.

In the following we propose a novel method that estimates the effect a split has on the SIP-bounds relevant to the verification problem. The method considers both the direct and indirect effects, but not propagation effects as empirical evidence obtained suggests that they do not have a significant impact in ReLU networks.

The effect of splitting a node is estimated with a score function s : N → R where N is the set of all hidden nodes in the FFNN. This score function is composed of the direct score sdir : N → R and the indirect score sindir : N → R, so s = sdir + sindir. These functions are formalised below, but first we introduce some necessary notation and concepts.

In the rest of this section let 〈f, ψin, ψout〉 be the verification problem as defined in Definition 1 where f : Rn → Rm is an FFNN with l layers. Moreover, let ϵl1,l2 n1,n2 denote the ESIP error from a node in layer l1 and position n1 at a node in layer l2 and position n2. Finally, let yn,low, yn,up be the concrete lower and upper bounds for output node n as calculated by ESIP, so yn,low ≤ f(x)n ≤ yn,up for all x satisfying ψin where f(x)n ∈ R is the n-th output of the neural network.

With this in mind, we now introduce a upper bound for the linear equations in ψout and motivate the score functions. Let aᵀ f(x) + b ≤ 0 be a linear constraint in ψout and notice that for all x satisfying ψin the equation gup = ∑n|(an>0) yn,upan + ∑n|(an<0) yn,lowan + b is an upper bound for aᵀ f(x)+b. Thus gup ≤ 0 is a sufficient condition for aᵀ f(x) + b ≤ 0 for all x satisfying ψin. The score function here proposed aims to identify splitting candidates that minimise the bound gup by minimising yn,up if an > 0 and maximising yn,low if an < 0 for all n.

To estimate the direct effect, recall from Section 2 that negative errors are added to yn,low and positive errors to yn,up and that the error from a node hl1 n1 in layer l1 and position n1 is removed after splitting the node. So, a good estimate of the reduction of gup due to the direct effect from splitting hl1 n1 is given by ˆ sdir (hl1 n1 ) = ∑n|(sign(ϵl1,l n1,n)=sign(an)) ϵl1,l n1,nan. where l is the network's output layer.

The score function ˆ sdir (hl1 n1 ) therefore estimates how much closer we are to proving that a single linear constraint aᵀ f(x) + b ≤ 0 is satisfied for all x satisfying ψin after the split. However, the output constraint ψout as defined in Definition 1, is generally a conjunction of disjunctions of linear constraints. To estimate the cumulative effect, we sum the coefficients in each individual linear constraint aᵀi,j f(x) + bi,j ≤ 0 in ψout, so ˆa = ∑i,j ai,j and:

Section 4 covers the details of which constraints are used at different stages of the algorithm. Notice that the score function above computes the same value as that in [Henriksen and Lomuscio, 2020] in the case of local robustness for classification networks, but generalises it and covers all verification problems conforming to Definition 1. An example of such a problem can be found in the commonly used ACAS Xu benchmarks, where the task is to verify that an output node is never minimal; this is not supported in [Henriksen and Lomuscio, 2020], but does indeed conform to Definition 1.

In order to estimate the indirect effect, first notice that nodes in the last hidden layer l - 1 do not have an indirect effect as there are no intermediate layers between layer l -1 and the output layer, so we take the score function to be s(hl-1 n2 ) = sdir(hl-1 n2 ). However, indirect effects are present for layers preceding l - 1.

Recall that s(hl-1 n2 ) estimates howmuch gup is reduced by removing the error from node hl-1 n2 . Thus, if splitting node hl-2 n1 removes an αl-2,l-1 n1,n2 ∈ [0, 1] fraction of the error from hl-1 n2 , then sindir(hl-2 n1 ) = ∑n2 αl-2,l-1 n1,n2 s(hl-1 n2 ) is a good estimate of the indirect effects from splitting hl-2 n1 .

We now estimate αl-2,l-1 n1,n2 . Since the node hl-1 n2 is governed by the ReLU activation function, the error is completely removed if the lower bound zl-1 n2,low is larger than 0, or if the upper bound zl-1 n2,up smaller than 0. Similar to the direct score function, we estimate that splitting node hl-2 n1 reduces the upper bound of node hl-1 n2 by ϵl-2,l-1 n1,n2 if ϵl-2,l-1 n1,n2 ≥ 0 (respectively, increases lower bound if ϵl-2,l-1 n1,n2 ≤ 0). This leads us to the estimation ˆαl-2,l-1 n1,n2 = ϵl-2,l-1 n1,n2 /zl-1 n2,up if ϵl-2,l-1 n1,n2 ≥ 0 otherwise ˆαl-2,l-1 n1,n2 = ϵl-2,l-1 n1,n2 /zl-1 n2,low . Generalising this reasoning to arbitrary layers l1 and l2 with l1 < l2, we have ˆαl1,l2 n1,n2 = ϵl1,l2 n1,n2 /zl2 n2,up if ϵl1,l2 n1,n2 ≥ 0 and ϵl1,l2 n1,n2 /zl2 n2,low otherwise; which results in the following score function:

The score function as presented above can also be used to estimate the effects of splitting input nodes. Input nodes do not have an error, thus sdir(h0n1 ) = 0; however, they do have an indirect effect that can be estimated with a small modification to αl1,l2 n1,n2 . Calculating the lower and upper bounds at node hl2 n2 requires minimising and maximising the linear equation ql2 n2 (x), respectively (See Section 2). So, if cl2 n1,n2 is the coefficient for xn1 in ql2 n2 (x), then it is reasonable to assume that by splitting xn1 at (max(xn1 ) +min(xn1 ))/2 the input interval of node hl2 n2 is reduced by cl2 n1,n2 (max(xn1 )-min(xn1 ))/2 in each branch. Thus, we use the estimate ˆα0,l2 n1,n2 = cl2 n1,n2 (max(xn1 ) - min(xn1 ))/2(zl2 n2,up - zl2 n2,low ).

In the next section we present an algorithm that uses the functions here described to improve the scalability of the verification problem.

4 Verification Algorithm

We now present a verification algorithm that extends previous symbolic interval propagation-based approaches [Wang et al., 2018a; Henriksen and Lomuscio, 2020] by introducing three novel contributions. Firstly, the method considers both direct and indirect effects of a split in order to locate the most influential split candidates. Secondly, when compared to the current SoA encoding for the adaptive splitting paradigm [Henriksen and Lomuscio, 2020], splits are encoded in a way that reduces the search space for counterexamples. Thirdly, we combine the relaxation from RSIP and ESIP to produce a precise relaxation that is also suitable for the novel splitting score-function and split-encoding.

The procedure is outlined in Algorithm 1. It takes as input an FFNN f, input constraints ψin and output constraints ψout; the output is the result of the verification problem. The main loop consists of three phases.

The symbolic interval propagation phase (lines 5 to 6) of the algorithm uses RSIP and ESIP (see Section 2) to calculate a linear relaxation for each term ˆ fi,j(x) = aᵀi,j f(x)+bi,j in the output constraint ψout = ∧i (∨j aᵀi,j f(x) + bi,j < 0). Notice that ˆ fi,j : Rn → R is an affine transformation of f(x) so ˆ fi,j is an FFNN in itself, and that the output con-

Algorithm 1 Verification Algorithm

1: f, ψin, ψout ←FFNN, inputConstraint, outputConstraint 2: queue← [Branch(splitConstraints=None, bounds=None)] 3: while (not queue.isEmpty) do 4: foundCandidateCex, branch ← False, queue.pop() 5: rsip ← RSIP(f, ψin, branch) 6: esip←ESIP(f, ψin, branch, rsip.bounds) 7: for each clause ψjout in ψout do 8: candidateCex ← LPSolver(esip, rsip, ψin, ψjout , ψsplit) 9: if (candidateCex is not None) then 10: foundCandidateCex = T rue 11: cex ← LocalSearch(candidateCex, f, ψin) 12: if (cex is not None) then return (UNSAFE, cex) 13: if (foundCandidateCex) then 14: splitNode←argmax(splitScore(f, ψin, ψout, esip))) 15: queue.add(Branch(branch, splitNode ≥ 0, esip.bounds)) 16: queue.add(Branch(branch, splitNode ≤ 0, esip.bounds)) 17: return SAFE

straint of the verification problem can be written as ψout = ∧i (∨j ˆ fi,j (x) < 0).

For each FFNN ˆ fi,j, the algorithm first runs RSIP, and then ESIP augmented with the bounds from RSIP. During ESIP, calculated bounds are compared with the bounds from RSIP on the fly and the tightest bounds are used to relax the ReLU nodes, resulting in a better relaxation than what can be achieved by ESIP alone.

Recall that RSIP produces linear bounds yrsip low (x)i,j and yrsip up (x)i,j such that yrsip low (x)i,j ≤ ˆ fi,j(x) ≤ yrsip up (x)i,j for all x satisfying ψout. In contrast, ESIP produces a linear equation q(x)i,j and concrete errors ϵni,j at the output node, such that q(x)i,j + ∑n|(ϵni,j <0) ϵni,j ≤ ˆ fi,j(x) ≤ q(x)i,j + ∑n|(ϵni,j >0) ϵni,j for all x satisfying ψout. Note that we here take ϵni,j to mean the errors at the output node of network ˆ fi,j

where n is an enumeration of all ReLU nodes in the network.

The search phase (lines 7-12) uses an LP-satisfiability call and a gradient descent-based local search in an attempt to locate a counterexample as first introduced in [Henriksen and Lomuscio, 2020]; however, our work differs in that we introduce a generalised version supporting any output constraint ψout, not just local robustness for classification problems.

A counterexample is an input x satisfying ψin that violates at least one of the clauses ψiout from ψout or, equivalently, satisfies ¬ψiout = (∧j ˆ fi,j(x) ≥ 0). Since ˆ fi,j(x) is not linear, the relaxations from the previous phase are used to encode necessary conditions for ¬ψiout .

The satisfiability call is performed for each clause ψiout with the constraints {∧j yrsip up (x)i,j ≥ 0, ∧j yesip (x,σ)i,j ≥ 0, ψin, ψsplit, (σk ∈ [0, 1])∀k}. Here, yrsip up (x)i,j is the upper bound produced by RSIP, ψsplit are the split constraints covered later in this section, and yesip (x,σ)i,j = q(x) + ∑n σ n ϵni,j is a novel encoding of the ESIP bound introducing new auxiliary variables σ n ∈ [0, 1]. Indeed, maxσ yesip i,j (x,σ) produces the same upper bound as in [Henriksen and Lomuscio, 2020]; however, we clarify below that σ can be further constrained during splitting thus resulting in tighter bounds.

If the calls are unsatisfiable for all clauses, then the network is safe by Theorem 1. However, if an assignment xi is found for at least one clause ψiout , then this assignment is a potential counterexample. The candidate xi is tested by evaluating f(xi ); if the result does not satisfy ψout, then xi is a valid counterexample and the network is unsafe; otherwise, the algorithm proceeds to the local search.

The local search phase (lines 11 to 12) performs a gradient descent with the loss function L(x) = - ∑j aᵀi,j f(x) for each potential counterexample xi and corresponding clause ψiout = ∨j aᵀi,j f(x) + bi,j ≤ 0. The gradient descent initialises at x0 = xi and after each step the resulting xk is clipped to the input box ψx and checked to see whether it is a valid counterexample. If a counterexample is found, then the network is unsafe under the given constraints; otherwise, the branch and bound phase is initiated.

The branch and bound phase (lines 13 to 16) starts by determining the next split with the method described in Section 3. A ReLU-split is performed by constraining the split node's input to ≤ 0 and ≥ 0 in a lower and upper branch, respectively. During ESIP and RSIP the ReLU constraints are handled trivially by only considering the corresponding linear parts of the ReLU; however, enforcing constraints in the LP-solver is somewhat more involved.

In the LP-solver, splitting a ReLU node with ESIP-errors ϵn and equation q(x) is handled by adding the constraint (z(x, σ) = q(x)+ ∑n σ n ϵn ) ≤ 0 in the lower branch (upper branch ≥ 0, respectively) where σ n ∈ [0, 1] are the auxiliary variables introduced in the search-phase. This is in contrast to previous algorithms that either split hierarchically from the first layer [Wang et al., 2018a], so split nodes have no errors, or use a worst-case bound (z(x) = q(x) + ∑n|ϵn>0 ϵn ) ≤ 0 [Henriksen and Lomuscio, 2020]. The encoding here presented is always at least as succinct as the one in the latter reference, but it is usually smaller as the σ variables encapsulate dependencies between the constraints.

Splitting an input node xk is performed by adjusting the split node's lower and upper bounds to (min(xk) + max(xk))/2 in ψin in the upper and lower branch, respectively. To ensure that the algorithm eventually splits all ReLUs, and thus is complete from Theorem 1, we only consider input nodes for splitting at branching-depths up to d ∈ N where d is an adjustable hyper-parameter.

After adding the new branch constraints, the main loop restarts with a few optimisations. Any clause ψiout that was proven not to have a counterexample in the previous branch cannot admit a counterexample in the new, stronger constrained branch and is disregarded. Moreover, the bounds of ESIP and RSIP only change for layers after the layer containing the split node and are thus only recalculated for those layers. Finally, we perform a low-cost relaxation utilising only ESIP in the first branch; RSIP is applied in later branches.

We complete this section by showing that Algorithm 1, as presented, is sound and complete.

Theorem 1 (Soundness and Completeness). Let 〈f, ψin, ψout〉 be a verification problem as defined in Definition 1. Algorithm 1 returns "unsafe" iff there exists an x satisfying ψin that violates ψout. It returns "safe" iff ψout is satisfied for all x satisfying ψin Proof sketch. This result follows from the fact that (i) the linear program is sound in the sense that all counterexamples are satisfiable solutions to the program and (ii) the algorithm is guaranteed to terminate after splitting all ReLU nodes and exploring the corresponding branches.

5 Implementation and Experimental Results

We implemented the algorithm proposed in Section 4 in a Python toolkit DEEPSPLIT [Henriksen and Lomuscio, 2021].

For benchmarking we used three MNIST fully-connected, two CIFAR10 convolutional and two MNIST convolutional networks from the Verification of Neural Networks Competition (VNN-COMP) [VNN20, 2020]. The convolutional nets were originally published in [Balunovic and Vechev, 2020]), while the FC nets were introduced in the VNN-COMP.

The performance of DEEPSPLIT was evaluated against all CPU-based toolkits that finished in the top 3 in at least one of the categories (VERINET [Henriksen and Lomuscio, 2020], VENUS [Botoeva et al., 2020], NNENUM [BAK, 2020], ERAN [Singh et al., 2018; Singh et al., 2019c]), as well as the MARABOU toolkit [Katz et al., 2019]. The remaining SIPbased toolkits were not used as [Henriksen and Lomuscio, 2020] outperformed [Wang et al., 2018a] and [Wang et al., 2018a] outperformed [Wang et al., 2018b] in experiments.

For all toolkits, we used the latest release of the implementations and ran them with the same settings as in the VNNCOMP for all toolkits with the exception of MARABOU and ERAN. MARABOU did not compete in VNN-COMP and we were unable to find the actual codebase that was used for ERAN in VNN-COMP; thus, we used the default settings for both. All toolkits are complete, except ERAN, which was used both in its complete and incomplete setting. However, note that incomplete and complete algorithms are not directly comparable, as discussed in the introduction. For the largest CIFAR10 network we reduced the ERAN incomplete LP and MILP timeouts per call to 1 second and the k-ReLU nodes to 2 due to runtimes of more than one hour per query.

While the results here are mostly in line with the VNNCOMP, note that each toolkit was run under different hardware in VNN-COMP, often with different results; in contrast, we here use the same experimental setup for all toolkits. Moreover, according to [VNN20, 2020], ERAN was augmented for the VNN-COMP with a PGD attack; this was not considered here for the reason above. Also, we limited the benchmarks to CPU-based toolkits. Because of this, we did not benchmark OVAL as it makes considerable use of GPUs.

The verification problem in all benchmarks concerns establishing whether the classification remains correct for all input perturbations within some l∞-ball. We used the same inputs and perturbation radii as in the VNN-COMP (0.02 and 0.05 for the MNIST FC networks and 0.1, 0.3, 2/255 and 8/255 for the two MNIST and CIFAR10 convolutional networks, respectively). For the FC networks we also used an additional intermediate perturbation radius of 0.03, which had a more balanced number of safe and unsafe cases. All benchmarks used a 1800 second timeout and were performed on a workstation with an Intel Core i9-10920X 3.5 GHz 12-core CPU, 128 GB Ram and Fedora 31 with Linux kernel 5.4.20.

Speedups reported below are calculated as follows: For toolkits A and B with nA and nB successfully verified cases, the speedup of toolkit A is tAnA /tBnA if nA < nB and tAnB /tBnB otherwise; here tAnA , tBnA is the total amount of time toolkit A and B used to verify their nA fastest verified cases. Thus we compare the cases where each toolkit performed the best.

In the following, we present the results of the experiments.

Fully-connected networks (Figure 1). Compared to the two most performing toolkits, DEEPSPLIT achieved speedups of 90 and 157 times (VERINET, NNENUM) and had 21% and 41% fewer timeouts (NNENUM, ERAN-COMPLETE). DEEPSPLIT was also the most performing considering safe and unsafe cases in isolation; DEEPSPLIT verified 109 safe and 68 unsafe cases, while the second most performing toolkits verified 107 safe (NNENUM) and 66 unsafe cases (VERINET).

Convolutional networks (Figure 2). Compared to the two most performing toolkits, DEEPSPLIT achieved speedups of 7.4 and 12.3 times (NNENUM, VERINET) and had 34% and 35% fewer timeouts (ERAN-COMPLETE, NNENUM). ERANCOMPLETE verified 251 safe cases while allowing long runtimes, whereas DEEPSPLIT only found 240 cases; however, ERAN-COMPLETE was among the slowest toolkits in general and only identified 3 unsafe cases, while DEEPSPLIT found 42. VERINET performed well on unsafe cases, also verifying 42; however, it only identified 188 safe-cases.

Ablation tests (Table 1). To determine the impact of the individual novel techniques proposed in this paper, we performed ablation experiments on the two-layer MNIST network. The results are reported in Table 1 for the following settings from top to bottom row: (i) all techniques enabled (ii) indirect score disabled (iii) RSIP disabled, only ESIP was used, (iv) worst-case encoding of split-constraints were used in the LP-Solver as discussed in Section 4. Each contribution resulted in significant overall gains. This is due to the fact that each method contributed to significantly reducing the number of branches that had to be explored during verification.

6 Conclusions

SoA methods for the verification of neural networks at present do not scale to the network sizes that are commonly used in applications such as computer vision.

In this paper we proposed a complete symbolic interval propagation-based verification algorithm which extends previous approaches with several contributions. We introduced a novel technique to determine promising split candidates, succinct LP-encodings of the split-constraints and a relaxation technique combining two methods for symbolic interval propagation. The experiments demonstrate that these improvements lead to a speedup of 1-2 orders of magnitude over the current SoA for most networks.

In future work we plan to investigate how similar techniques can be used to improve the performance for networks with Sigmoid and Tanh activation functions.

Acknowledgements

This work is partly funded by the UKRI Centre for Doctoral Training in Safe and Trusted Artificial Intelligence (grant number EP/S023356/1). Alessio Lomuscio is supported by a Royal Academy of Engineering Chair in Emerging Technologies.