--- title: "Reachability Analysis of Deep Neural Networks with Provable Guarantees" authors: - "Wenjie Ruan" - "Xiaowei Huang" - "Marta Kwiatkowska" doi: "10.24963/ijcai.2018/368" url: "https://ijcai.org/proceedings/2018/368" pdf: "https://ijcai.org/proceedings/2018/0368.pdf" venue: "Proceedings of the Twenty-Seventh International Joint Conference on Artificial Intelligence" year: 2018 pages: "2651-2659" track: "Main track" keywords: - "Machine Learning: Neural Networks" - "Machine Learning: Deep Learning" - "Computer Vision: Computer Vision" - "Agent-based and Multi-agent Systems: Formal Verification, Validation and Synthesis" --- # Reachability Analysis of Deep Neural Networks with Provable Guarantees Wenjie Ruan, Xiaowei Huang, Marta Kwiatkowska *Proceedings of the Twenty-Seventh International Joint Conference on Artificial Intelligence*, pages 2651-2659. DOI: [10.24963/ijcai.2018/368](https://doi.org/10.24963/ijcai.2018/368) ## Abstract Verifying correctness for deep neural networks (DNNs) is challenging. We study a generic reachability problem for feed-forward DNNs which, for a given set of inputs to the network and a Lipschitz-continuous function over its outputs computes the lower and upper bound on the function values. Because the network and the function are Lipschitz continuous, all values in the interval between the lower and upper bound are reachable. We show how to obtain the safety verification problem, the output range analysis problem and a robustness measure by instantiating the reachability problem. We present a novel algorithm based on adaptive nested optimisation to solve the reachability problem. The technique has been implemented and evaluated on a range of DNNs, demonstrating its efficiency, scalability and ability to handle a broader class of networks than state-of-the-art verification approaches. ## Full text ## 1 Introduction Concerns have been raised about the suitability of deep neural networks (DNNs), or systems with DNN components, for deployment in safety-critical applications, see e.g., [Amodei et al., 2016; Sun et al., 2018]. To ease this concern and gain users' trust, DNNs need to be certified similarly to systems such as airplanes and automobiles. In this paper, we propose to study a generic reachability problem which, for a given DNN, an input subspace and a function over the outputs of the network, computes the upper and lower bounds over the values of the function. The function is generic, with the only requirement that it is Lipschitz continuous. We argue that this problem is fundamental for certification of DNNs, as it can be instantiated into several key correctness problems, including adversarial example generation [Szegedy et al., 2013; Goodfellow et al., 2014], safety verification [Huang et al., 2017; Katz et al., 2017; Ruan et al., 2018b], output range analysis [Lomuscio and Maganti, 2017; Dutta et al., 2017], and robustness comparison. To certify a system, a certification approach needs to provide not only a result but also a guarantee over the result, such as the error bounds. Existing approaches for analysing DNNs with a guarantee work by either reducing the problem to a constraint satisfaction problem that can be solved by MILP [Lomuscio and Maganti, 2017; Cheng et al., 2017; Bunel et al., 2017; Xiang et al., 2017], SAT [Narodytska et al., 2017] or SMT [Katz et al., 2017; Bunel et al., 2017] techniques, or applying search algorithms over discretised vector spaces [Huang et al., 2017; Wicker et al., 2018]. Even though they are able to achieve guarantees, they suffer from two major weaknesses. Firstly, their subjects of study are restricted. More specifically, they can only work with layers conducting linear transformations (such as convolutional and fully-connected layers) and simple non-linear transformations (such as ReLU), and cannot work with other important layers, such as the sigmoid, max pooling and softmax layers that are widely used in state-of-the-art networks. Secondly, the scalability of the constraint-based approaches is significantly limited by both the capability of the solvers and the size of the network, and they can only work with networks with a few hundreds of hidden neurons. However, state-of-the-art networks usually have millions, or even billions, of hidden neurons. This paper proposes a novel approach to tackle the generic reachability problem, which does not suffer from the above weaknesses and provides provable guarantees in terms of the upper and lower bounds over the errors. The approach is inspired by recent advances made in the area of global optimisation [Gergel et al., 2016; Grishagin et al., 2018]. For the input subspace defined over a set of input dimensions, an adaptive nested optimisation algorithm is developed. The performance of our algorithm is not dependent on the size of the network and it can therefore scale to work with large networks. Our algorithm assumes certain knowledge about the DNN. However, instead of directly translating the activation functions and their parameters (i.e., weights and bias) into linear constraints, it needs a Lipschitz constant of the network. For this, we show that several layers that cannot be directly translated into linear constraints are actually Lipschitz continuous, and we are able to compute a tight Lipschitz constant by analysing the activation functions and their parameters. We develop a software tool DeepGO 1 and evaluate its performance by comparing with existing constraint-based approaches, namely, SHERLOCK [Dutta et al., 2017] and Reluplex [Katz et al., 2017]. We also demonstrate our tool on DNNs that are beyond the capability of existing tools. 1 Available on https://github.com/trustAI/DeepGO. ## 2 RelatedWorks We discuss several threads of work concerning problems that can be obtained by instantiating our generic reachability problem. Their instantiations are explained in the paper. Due to space limitations, this review is by no means complete. Safety Verification There are two ways of achieving safety verification for DNNs. The first is to reduce the problem into a constraint solving problem. Notable works include, e.g., [Pulina and Tacchella, 2010; Katz et al., 2017]. However, they can only work with small networks with hundreds of hidden neurons. The second is to discretise the vector spaces of the input or hidden layers and then apply exhaustive search algorithms or Monte Carlo tree search algorithm on the discretised spaces. The guarantees are achieved by establishing local assumptions such as minimality of manipulations in [Huang et al., 2017] and minimum confidence gap for Lipschitz networks in [Wicker et al., 2018]. Adversarial Example Generation Most existing works, e.g., [Szegedy et al., 2013; Goodfellow et al., 2014; Nguyen et al., 2014; Moosavi-Dezfooli et al., 2016; Carlini and Wagner, 2016], apply various heuristic algorithms, generally using search algorithms based on gradient descent or evolutionary techniques. [Papernot et al., 2015] construct a saliency map of the importance of the pixels based on gradient descent and then modify the pixels. In contrast with our approach based on global optimisation and works on safety verification, these methods may be able to find adversarial examples efficiently, but are not able to conclude the nonexistence of adversarial examples when the algorithm fails to find one. Output Range Analysis The safety verification approach can be adapted to work on this problem. Moreover, [Lomuscio and Maganti, 2017] consider determining whether an output value of a DNN is reachable from a given input subspace, and propose an MILP solution. [Dutta et al., 2017] study the range of output values from a given input subspace. Their method interleaves local search (based on gradient descent) with global search (based on reduction to MILP). Both approaches can only work with small networks. ## 3 Lipschitz Continuity of DNNs This section shows that feed-forward DNNs are Lipschitz continuous. Let f : Rn → Rm be a N-layer network such that, for a given input x ∈ Rn , f(x) = {c1, c2, ..., cm} ∈ Rm represents the confidence values for m classification labels. Specifically, we have f(x) = fN(fN-1(...f1(x;W1, b1);W2, b2); ...);WN, bN) where Wi and bi for i = 1, 2, ...,N are learnable parameters and fi(zi-1;Wi-1, bi-1) is the function mapping from the output of layer i - 1 to the output of layer i such that zi-1 is the output of layer i - 1. Without loss of generality, we normalise the input to lie x ∈ [0, 1]n . The output f(x) is usually normalised to be in [0, 1]m with a softmax layer. Definition 1 (Lipschitz Continuity) Given two metric spaces (X, dX) and (Y, dY ), where dX and dY are the metrics on the sets X and Y respectively, a function f : X → Y is called Lipschitz continuous if there exists a real constant K ≥ 0 such that, for all x1, x2 ∈ X: K is called the Lipschitz constant for the function f. The smallest K is called the Best Lipschitz constant, denoted as Kbest . [Szegedy et al., 2013] show that deep neural networks with half-rectified layers (i.e., convolutional or fully connected layers with ReLU activation functions), max pooling and contrast-normalization layers are Lipschitz continuous. They prove that the upper bound of the Lipschitz constant can be estimated via the operator norm of learned parameters W. Next, we show that the softmax layer, sigmoid and Hyperbolic tangent activation functions also satisfy Lipschitz continuity. First we need the following lemma [Sohrab, 2003]. Lemma 1 Let f : Rn → Rm, if ||∂f(x)/∂x|| ≤ K for all x ∈ [a, b]n , then f is Lipschitz continuous on [a, b]n and K is its Lipschitz constant, where ||∗|| represents a norm operator. Based on this lemma, we have the following theorem. Theorem 1 Convolutional or fully connected layers with the sigmoid activation function s(Wx + b), Hyperbolic tangent activation function t(Wx + b), and softmax function p(x)j are Lipschitz continuous and their Lipschitz constants are 12 ‖W‖,‖W‖, and supi,j(‖xi‖+∥∥xixj ∥∥), respectively. Proof 1 First of all, we show that the norm operators of their Jacobian matrices are bounded. (1) Layer with sigmoid activation s(q) = 1/(1+e-q ) with q = Wx + b: (2) Layer with Hyperbolic tangent activation function t(q) = 2/(1 + e-2q ) - 1 with q = Wx + b: (3) Layer with softmax function p(x)j = exj /(∑nk=1 exk ) for j = 1, ..., m and n = m (dimensions of input and output of softmax are the same): Since the softmax layer is the last layer of a deep neural network, we can estimate its supremum based on Lipschitz constants of previous layers and box constraints of DNN's input. The final conclusion follows by Lemma 1 and the fact that all the layer functions are bounded on their Jacobian matrix. ## 4 ProblemFormulation Let o : [0, 1]m → R be a Lipschitz continuous function statistically evaluating the outputs of the network. Our problem is to find its upper and lower bounds given the set X′ of inputs to the network. Because both the network f and the function o are Lipschitz continuous, all values between the upper and lower bounds have a corresponding input, i.e., are reachable. Definition 2 (Reachability of Neural Network) Let X′ ⊆ [0, 1]n be an input subspace and f : Rn → Rm a network. The reachability of f over the function o under an error tolerance ϵ ≥ 0 is a set R(o, X′ , ϵ) = [l, u] such that We write u(o,X′ , ϵ) = u and l(o,X′ , ϵ) = l for the upper and lower bound, respectively. Then the reachability diameter is Assuming these notations, we may write D(o,X′ , ϵ; f) if we need to explicitly refer to the network f. In the following, we instantiate o with a few concrete functions, and show that several key verification problems for DNNs can be reduced to our reachability problem. Definition 3 (Output Range Analysis) Given a class label j ∈ [1, ..,m], we let o = Πj such that Πj((c1, ..., cm)) = cj . We write cj(x) = Πj(f(x)) for the network's confidence in classifying x as label j. Intuitively, output range [Dutta et al., 2017] quantifies how a certain output of a deep neural network (i.e., classification probability of a certain label j) varies in response to a set of DNN inputs with an error tolerance ϵ. Output range analysis can be easily generalised to logit 2 range analysis. We show that the safety verification problem [Huang et al., 2017] can be reduced to solving the reachability problem. Definition 4 (Safety) A network f is safe with respect to an input x and an input subspace X′ ⊆ [0, 1]n with x ∈ X′ , written as S(f, x,X′ ), if We have the following reduction theorem. Theorem 2 A network f is safe with respect to x and X′ s.t. x ∈ X′ if and only if u(⊕, X′ , ϵ) ≤ 0, where ⊕(c1, ..., cm) = maxi∈{1..m}(Πi(c1, ..., cm) - Πj(c1, ..., cm)) and j = argmaxj cj(x). The error bound of the safety decision problem by this reduction is 2ϵ. It is not hard to see that the adversarial example generation [Szegedy et al., 2013], which is to find an input x′ ∈ X′ such that argmaxj cj(x′ ) ̸= argmaxj cj(x), is the dual problem of the safety problem. The following two problems define the robustness comparisons between the networks and/or the inputs. 2 Logit output is the output of the layer before the softmax layer. The study of logit outputs is conducted in, e.g., [Papernot et al., 2015; Dutta et al., 2017]. Definition 5 (Robustness) Given two homogeneous 3 networks f and g, we say that f is strictly more robust than g with respect to a function o, an input subspace X′ and an error bound ϵ, written as Ro,X′ ,ϵ(f, g), if D(o,X′ , ϵ; f) < D(o, X′ , ϵ; g). Definition 6 Given two input subspaces X′ and X′′ and a network f, we say that f is more robust on X′ than on X′′ with respect to a statistical function o and an error bound ϵ, written as Rf,o,ϵ(X′ , X′′ ), if D(o, X′ , ϵ) < D(o, X′′ , ϵ). Thus, by instantiating the function o, we can quantify the output/logit range of a network, evaluate whether a network is safe, and compare the robustness of two homogeneous networks or two input subspaces for a given network. ## 5 Confidence Reachability with Guarantees Section 3 shows that a trained deep neural network is Lipschitz continuous regardless of its layer depth, activation functions and number of neurons. Now, to solve the reachability problem, we need to find the global minimum and maximum values given an input subspace, assuming that we have a Lipschitz constant K for the function o·f. In the following, we let w = o·f be the concatenated function. Without loss of generality, we assume the input space X′ is a box-constraint, which is clearly feasible since images are usually normalized into [0, 1]n before being fed into a neural network. The computation of the minimum value is reduced to solving the following optimization problem with guaranteed convergence to the global minimum (the maximization problem can be transferred into a minimization problem): However, the above problem is very difficult since w(x) is a highly non-convex function which cannot be guaranteed to reach the global minimum by regular optimization schemes based on gradient descent. Inspired by an idea from optimisation, see e.g., [Piyavskii, 1972; Torn and Zilinskas, 1989], we design another continuous function h(x, y), which serves as a lower bound of the original function w(x). Specifically, we need Furthermore, for i ≥ 0, we let Yi = {y0, y1, ..., yi} be a finite set containing i + 1 points from the input space [a, b]n , and let Yi ⊆ Yk when k > i, then we can define a function H(x; Yi) = maxy∈Yi h(x, y) which satisfies the following relation: We use li = infx∈[a,b]n H(x; Yi) to denote the minimum value of H(x; Yi) for x ∈ [a, b]n . Then we have 3 Here, two networks are homogeneous if they are applied on the same classification task but may have different network architectures (layer numbers, layer types, etc) and/or parameters. Figure 1: A lower-bound function designed via Lipschitz constant ![Image](https://ijcai.org/proceedings/2018/fulltext/368/image_000000_0d651fcf619df1009362577e89a934cb3fc6395f3747d26217b5fde55883b117.png) Similarly, we need a sequence of upper bounds ui to have By Expression (12), we can have the following: Therefore, we can asymptotically approach the global minimum. Practically, we execute a finite number of iterations by using an error tolerance ϵ to control the termination. In next sections, we present our approach, which constructs a sequence of lower and upper bounds, and show that it can converge with an error bound. To handle the high-dimensionality of DNNs, our approach is inspired by the idea of adaptive nested optimisation in [Gergel et al., 2016], with significant differences in the detailed algorithm and convergence proof. ## 5.1 One-dimensional Case We first introduce an algorithm which works over one dimension of the input, and therefore is able to handle the case of x ∈ [a, b] in Eqn. (8). The multi-dimensional optimisation algorithm will be discussed in Section 5.2 by utilising the onedimensional algorithm. We define the following lower-bound function. where K > Kbest is a Lipschitz constant of w and H(x; Yi) intuitively represents the lower-bound sawtooth function shown as Figure 1. The set of points Yi is constructed recursively. Assuming that, after (i - 1)-th iteration, we have Yi-1 = {y0, y1, .., yi-1}, whose elements are in ascending order, and sets The elements in sets w(Yi-1), Li-1 and Ui-1 have been defined earlier. The set Zi-1 records the smallest values zk computed in an interval [yk-1, yk]. In i-th iteration, we do the following sequentially: - Compute yi = arg infx∈[a,b] H(x; Yi-1) as follows. Let z∗ = minZi-1 and k be the index of the interval [yk-1, yk] where z∗ is computed. Then we let and have that yi ∈ (yk-1, yk). - Let Yi = Yi-1 ∪ {yi}, then reorder Yi in ascending order, and update w(Yi) = w(Yi-1) ∪ {w(yi)}. - Calculate and update Zi = (Zi-1 \ {z∗ }) ∪ {zi-1, zi}. - Calculate the new lower bound li = infx∈[a,b] H(x; Yi) by letting li = minZi, and updating Li = Li-1 ∪ {li}. - Calculate the new upper bound ui = miny∈Yi w(y) by letting ui = min{ui-1, w(yi)}. We terminate the iteration whenever |ui - li| ≤ ϵ, and let the global minimum value be y∗ = minx∈[a,b] H(x; Yi) and the minimum objective function be w∗ = w(y∗ ). Intuitively, as shown in Fig. 1, we iteratively generate lower bounds (by selecting in each iteration the lowest point in the saw-tooth function in the figure) by continuously refining a piecewise-linear lower bound function, which is guaranteed to below the original function due to Lipschitz continuity. The upper bound is the lowest evaluation value of the original function so far. ## Convergence Analysis In the following, we show the convergence of this algorithm to the global minimum by proving the following conditions. - Convergence Condition 1: lim i→∞ li = min x∈[a,b] w(x) - Convergence Condition 2: limi→∞(ui - li) = 0 ## Proof 2 (Monotonicity of Lower/Upper Bound Sequences) First, we prove that the lower bound sequence Li is strictly monotonic. Because and li-1 = minZi. To show that li > li-1, we need to prove zi-1 > z∗ and zi > z∗ . By the algorithm, z∗ is computed from interval [yk-1, yk], so we have We then have Since yi < yk and K > Kbest, by Lipschitz continuity we have zi-1 > z∗ . Similarly, we can prove zi > z∗ . Thus li > li-1 is guaranteed. Second, the monotonicity of upper bounds ui can be seen from the algorithm, since ui is updated to min{ui, w(yi)} in every iteration. ## Proof 3 (Convergence Condition 1) Since Yi-1 ⊆ Yi, we have H(x; Yi-1) ≤ H(x; Yi). Based on Proof 2, we also have li-1 < li. Then since the lower bound sequence {l0, l1, ..., li} is strictly monotonically increasing and bounded from above by minx∈[a,b] w(x). Thus limi→∞ li = minx∈[a,b] w(x) holds. ## Proof 4 (Convergence Condition 2) Since limi→∞ li = minx∈[a,b] w(x), we show limi→∞(ui - li) = 0 by showing that limi→∞ ui = minx∈[a,b] w(x). Since Yi = Yi-1∪{yi} and yi ∈ X = [a, b], we have limi→∞ Yi = X. Then we have limi→∞ui = limi→∞ infy∈Yi w(y) = inf w(X). Since X = [a, b] is a closed interval, we can prove limi→∞ ui = inf w(X) = minx∈[a,b] w(x). ## Dynamically Improving the Lipschitz Constant A Lipschitz constant closer to Kbest can greatly improve the speed of convergence of the algorithm. We design a practical approach to dynamically update the current Lipschitz constant according to the information obtained from the previous iteration: where η > 1. We emphasise that, because this dynamic update does not compromise the convergence. ## 5.2 Multi-dimensional Case The basic idea is to decompose a multi-dimensional optimization problem into a sequence of nested one-dimensional subproblems. Then the minima of those one-dimensional minimization subproblems are back-propagated into the original dimension and the final global minimum is obtained. We first introduce the definition of k-th level subproblem. Definition 7 The k-th level optimization subproblem, written as φk(x1, ..., xk), is defined as follows: for 1 ≤ k ≤ n - 1, and for k = n, Combining Expression (23) and Definition 7, we have that which is actually a one-dimensional optimization problem and therefore can be solved by the method in Section 5.1. However, when evaluating the objective function φ1(x1) at x1 = a1, we need to project a1 into the next one-dimensional subproblem We recursively perform the projection until we reach the n-th level one-dimensional subproblem, Once solved, we back-propagate objective function values to the first-level φ1(a1) and continue searching from this level until the error bound is reached. ## Convergence Analysis We use mathematical induction to prove convergence for the multi-dimension case. - Base case: for all x ∈ R, limi→∞ li = infx∈[a,b] w(x) and limi→∞(ui - li) = 0 hold. - Inductive step: if, for all x ∈ Rk , limi→∞ li = infx∈[a,b]k w(x) and limi→∞(ui - li) = 0 are satisfied, then, for all x ∈ Rk+1 , limi→∞ li = infx∈[a,b]k+1 w(x) and limi→∞(ui - li) = 0 hold. The base case (i.e., one-dimensional case) is already proved in Section 5.1. Now we prove the inductive step. Proof 5 By the nested optimization scheme, we have Since miny∈[ai,bi]k w(x, y) is bounded by an interval error ϵy, assuming Φ∗ (x) is the accurate global minimum, then we have So the k + 1-dimensional problem is reduced to the onedimensional problem minx∈[a,b] Φ(x). The difference from the real one-dimensional case is that evaluation of Φ(x) is not accurate but bounded by |Φ(x) - Φ∗ (x)| ≤ ϵy, ∀x ∈ [a, b], where Φ∗ (x) is the accurate function evaluation. Assuming that the minimal value obtained from our method is Φ∗min = minx∈[a,b] Φ∗ (x) under accurate function evaluation, then the corresponding lower and upper bound sequences are {l∗ 0 , ..., l∗ i } and {u∗0 , ..., u∗i }, respectively. For the inaccurate evaluation case, we assume Φmin = minx∈[a,b] Φ(x), and its lower and bound sequences are, respectively, {l0, ..., li} and {u0, ..., ui}. The termination criteria for both cases are |u∗i - l∗ i | ≤ ϵx and |ui - li| ≤ ϵx, and φ∗ represents the ideal global minimum. Then we have φ∗ - ϵx ≤ li. Assuming that l∗ i ∈ [xk, xk+1] and xk, xk+1 are adjacent evaluation points, then due to the fact that l∗ i = infx∈[a,b] H(x; Yi) we have Since |Φ(xi) - Φ∗ (xi)| ≤ ϵy, ∀i = k, k + 1, we thus have Based on the search scheme, we know that and thus we have φ∗ - li ≤ ϵy + ϵx. Similarly, we can get so ui-φ∗ ≤ ϵx+ϵy. By φ∗ -li ≤ ϵy+ϵx and the termination criteria ui-li ≤ ϵx, we have li-ϵy ≤ φ∗ ≤ ui+ϵy, i.e., the accurate global minimum is also bounded. The proof indicates that the overall error bound of the nested scheme only increases linearly w.r.t. the bounds in the one-dimensional case. Moreover, an adaptive approach can be applied to optimise its performance without compromising convergence. The key observation is to relax the strict subordination inherent in the nested scheme and simultaneously consider all the univariate subproblems arising in the course of multidimensional optimization. For all the generated subproblems that are active, a numerical measure is applied. Then an iteration of the multidimensional optimization consists in choosing the subproblem with maximal measurement and carrying out a new trial within this subproblem. The measure is defined to be the maximal interval characteristics generated by the one-dimensional optimisation algorithm. ## 5.3 Proof of NP-completeness We prove NP-completeness of our method. For space reasons we only describe the proof idea; for the full proof see [Ruan et al., 2018a]. For the upper bound, we first show that finding the optimal value for the one-dimensional case can be done in polynomial time with respect to the error bound ϵ. Then, for the multi-dimensional case, we have a non-deterministic algorithm to first guess a subset of dimensions and then conduct the one-dimensional optimisation one by one. The entire procedure can be done in polynomial time with a nondeterministic automaton, i.e., in NP. For the lower bound, we show a reduction from the 3-SAT problem. For any instance ϕ of 3-SAT, we can construct a network f and an evaluation function o, such that the satisfiability of ϕ is equivalent to non-reachability of value 0 for the function w = o · f . ## 6 Experiments ## 6.1 Comparison with State-of-the-art Methods Two methods are chosen as baseline methods in this paper: - Reluplex [Katz et al., 2017]: an SMT-based method for solving queries on DNNs with ReLU activations; we apply a bisection scheme to compute an interval until an error is reached - SHERLOCK [Dutta et al., 2017]: aMILP-basedmethod dedicated to output range analysis on DNNs with ReLU activations. Our software is implemented in Matlab 2018a, running on a notebook computer with i7-7700HQ CPU and 16GB RAM. Figure 2: Comparison with SHERLOCK and Reluplex Since Reluplex and SHERLOCK (not open-sourced) are designed on different software platforms, we take their experimental results from [Dutta et al., 2017], whose experimental environment is a Linux workstation with 63GB RAM and 23Cores CPU (more powerful than ours) and ϵ = 0.01. Following the experimental setup in [Dutta et al., 2017], we use their data (2-input and 1-output functions) to train six neural networks with various numbers and types of layers and neurons. The input subspace is X′ = [0, 10]2 . The comparison results are given in Fig. 2. They show that, while the performance of both Reluplex and SHERLOCK is considerably affected by the increase in the number of neurons and layers, our method is not. For the six benchmark neural networks, our average computation time is around 5s, 36 fold improvement over SHERLOCK and nearly 100 fold improvement over Reluplex (excluding timeouts). We note that our method is running on a notebook PC, which is significantly less powerful than the 23-core CPU stations used for SHERLOCK and Reluplex. ## 6.2 Safety and Robustness Verification by Reachability Analysis We use our tool to conduct logit and output range analysis. Seven convolutional neural networks, represented as DNN1,...,DNN-7, were trained on the MNIST dataset. Images are resized into 14×14 to enforce that a DNN with deeper layers tends to over-fit. The networks have different layer types, including ReLu, dropout and normalization, and the number of layers ranges from 5 to 19. Testing accuracies range from 95% to 99%, and ϵ = 0.05 is used in our experiments. We randomly choose 20 images (2 images per label) and manually choose 4 features such that each feature contains 8 pixels, i.e., X′ = [0, 1]8 . Fig. 3 (a) illustrates the four features and the architecture of two DNNs with the shallowest and deepest layers, i.e., DNN-1 and DNN-7. Safety Verification Fig. 4 (a) shows an example: for DNN1, Feature-4 is guaranteed to be safe with respect to the image x and the input subspace X′ . Specifically, the reachability interval is R(Π0,X′ , ϵ) = [74.36%, 99.98%], which means that l(Π0,X′ , ϵ) = 74.36%. By this, we have u(⊕-0,X′ , ϵ) ≤ (1 - 0.7436) < 0.7436 = l(Π0,X′ , ϵ). Then, by Theorem 2, we have S(DNN-1, x,X′ ). Intuitively, no matter how we manipulate this feature, the worst case is to reduce the confidence of output being '0' from 99.95% (its original confidence probability) to 74.36%. Statistical Comparison of Safety Fig. 4 (b) compares the ratios of safe images for different DNNs and features. It shows that: i) no DNN is 100% safe on those features: DNN-6 is the safest one and DNN-1, DNN-2 and DNN-3 are less safe, which means a DNN with well chosen layers are safer than those DNNs with very shallow or deeper layers; and ii) the safety performance of different DNNs is consistent for the same feature, which suggests that the feature matters - some features are easily perturbed to yield adversarial examples, e.g., Feature-1 and Feature-2. Figure 3: (a) The four features and the architecture of DNN-1 and DNN-7. (b) Left: boxplots of confidence reachability diameters for 7 DNNs, based on 4 × 20 analyses of each DNN. Right: boxplot of confidence reachability diameters for 4 features, based on 7 × 20 analyses of each feature. The red line represents the median value: a lower value indicates a more robust model or feature. ![Image](https://ijcai.org/proceedings/2018/fulltext/368/image_000001_10b772c155b73eac6160338c78ce19c0291a7b4ab645c63e179023ab78092d7f.png) Statistical Comparison of Robustness Fig. 3 (b) compares the robustness of networks and features with two boxplots over the reachability diameters, where the function o is Πj for a suitable j. We can see that DNN-6 and DNN-5 are the two most robust, while DNN-1, DNN-2 and DNN-3 are less robust. Moreover, Feature-1 and Feature-2 are less robust than Feature-3 and Feature-4. We have thus demonstrated that reachability analysis with our tool can be used to quantify the safety and robustness of deep learning models. In the following, we perform a comparison of networks over a fixed feature. Safety Comparison of Networks By Fig. 4 (c), DNN-4 and DNN-6 are guaranteed to be safe w.r.t. the subspace defined by Feature-3. Moreover, the output range of DNN-7 is [1.8%, 100.0%], which means that we can generate adversarial images by only perturbing this feature, among which the worst one is as shown in the figure with a confidence 1.8%. Thus, reachability analysis not only enables qualitative safety verification (i.e., safe or not safe), but also allows benchmarking of safety of different deep learning models in a principled, quantitive manner (i.e., how safe) by quantifying the 'worst' adversarial example. Moreover, compared to retraining the model with 'regular' adversarial images, these 'worst' adversarial images are more effective in improving the robustness of DNNs [Kolter and Wong, 2017]. Robustness Comparison of Networks The bar chart in Fig. 4 (c) shows the reachability diameters of the networks over Feature-3, where the function o is Πj. DNN-4 is the most robust one, and its output range is [94.2%, 100%]. ## 6.3 AComprehensive Comparison with the State-of-the-arts This section presents a comprehensive, high-level comparison of our method with several existing approaches that have been used for either range analysis or verification of DNNs, including SHERLOCK [Dutta et al., 2017], Reluplex [Katz et al., 2017], Planet [Ehlers, 2017], MIP [Cheng et al., 2017; Lomuscio and Maganti, 2017] and BaB [Bunel et al., 2017], as shown in Fig. 5. Core Techniques Most existing approaches (SHERLOCK, Reluplex, Planet, MIP) are based on reduction to constraint solving, except for BaB which mixes constraint solving with local search. On the other hand, our method is based on global optimization and assumes Lipschitz continuity of the networks. As indicated in Section 3, all known layers used in classification tasks are Lipschitz continuous. Workable Layer Types While we are able to work with all known layers used in classification tasks because they are Lipschitz continuous (proved in Section 3 of the paper), Planet, MIP and BaB can only work with Relu and Maxpooling, and SHERLOCK and Reluplex can only work with Relu. Running Time on ACAS-Xu Network We collect running time data from [Bunel et al., 2017] on the ACAS-Xu network, and find that our approach has similar performance to BaB, and better than the others. No experiments for SHERLOCK are available. We reiterate that, compared to their experimental platform (Desktop PC with i7-5930K CPU, 32GB RAM), ours is less powerful (Laptop PC with i7-7700HQ CPU, 16GB RAM). We emphasise that, although our approach performs well on this network, the actual strength of our approach is not the running time on small networks such as ACAS-Xu, but the ability to work with large-scale networks (such as those shown in Section 6.2). Computational Complexity While all the mentioned approaches are in the same complexity class, NP, the complexity of our method is with respect to the number of input dimensions to be changed, as opposed to the number of hidden neurons. It is known that the number of hidden neurons is much larger than the number of input dimensions, e.g., there are nearly 6.5 × 106 neurons in AlexNet. Applicable to State-of-the-art Networks We are able to work with state-of-the-art networks with millions of neurons. However, the other tools (Reluplex, Planet, MIP, BaB) can only work with hundreds of neurons. SHERLOCK can work with thousands of neurons thanks to its interleaving of MILP Figure 4: (a) Left: an original image (logit is 11.806, confidence of output being '0' is 99.95%), where area marked by dashed line is the feature. Middle: an image on the confidence lower bound. Right: an image on the confidence upper bound; for the output label '0', the feature's output range is [74.36%, 99.98%], and logit reachability is [7.007, 13.403]. (b) Ratios of safe images for 7 DNNs and 4 features. (c) A detailed example comparing the safety and robustness of DNNs for image '9' and Feature-3: the top number in the caption of each figure is logit and the bottom one is confidence; the unsafe cases are all misclassified as '8'; the last bar chart shows their confidence reachability diameters. ![Image](https://ijcai.org/proceedings/2018/fulltext/368/image_000002_095add6eaf10e67563c79fc50315cca0ef52eed9c3aa8ef20f78b9945186ff0d.png) with local search. Maximum Number of Layers in Tested DNNs We have validated our method on networks with 19 layers, whereas the other approaches are validated on up to 6 layers. In summary, the key advantages of our approach are as follows: i) the ability to work with large-scale state-of-theart networks; ii) lower computational complexity, i.e., NPcompleteness with respect to the input dimensions to be changed, instead of the number of hidden neurons; and iii) the wide range of types of layers that can be handled. ## 7 Conclusion We propose, design and implement a reachability analysis tool for deep neural networks, which has provable guarantees and can be applied to neural networks with deep layers and nonlinear activation functions. The experiments demonstrate that our tool can be utilized to verify the safety of deep neural networks and quantitatively compare their robustness. We envision that this work marks an important step towards a practical, guaranteed safety verification for DNNs. Future work includes parallelizing this method in GPUs to improve its scalability on large-scale models trained on ImageNet, and a generalisation to other deep learning models such as RNNs and deep reinforcement learning. ## Acknowledgements WR and MK are supported by the EPSRC Programme Grant on Mobile Autonomy (EP/M019918/1). XH acknowledges NVIDIA Corporation for its support with the donation of the Titan Xp GPU, and is partially supported by NSFC (no. 61772232). ## References 1. [Amodei et al., 2016] Dario Amodei, Chris Olah, Jacob Steinhardt, Paul Christiano, John Schulman, and Dan Mané. Concrete problems in ai safety. arXiv preprint arXiv:1606.06565, 2016. 2. [Bunel et al., 2017] Rudy Bunel, Ilker Turkaslan, Philip HS Torr, Pushmeet Kohli, and M Pawan Kumar. Piecewise linear neural network verification: A comparative study. arXiv preprint arXiv:1711.00455, 2017. 3. [Carlini and Wagner, 2016] Nicholas Carlini and David A. Wagner. Towards evaluating the robustness of neural net- works. CoRR, abs/1608.04644, 2016. 4. [Cheng et al., 2017] Chih-Hong Cheng, Georg Nührenberg, and Harald Ruess. Maximum resilience of artificial neu- ral networks. In Deepak D'Souza and K. Narayan Kumar, editors, Automated Technology for Verification and Anal- ysis, pages 251-268, Cham, 2017. Springer International Publishing. 5. [Dutta et al., 2017] Souradeep Dutta, Susmit Jha, Sriram Sanakaranarayanan, and Ashish Tiwari. Output range analysis for deep neural networks. arXiv preprint arXiv:1709.09130, 2017. 6. [Ehlers, 2017] Ruediger Ehlers. Formal verification of piece- wise linear feed-forward neural networks. In International Symposium on Automated Technology for Verification and Analysis, pages 269-286. Springer, 2017. Figure 7. A high-level comparison with state-of-the-art methods: SHERLOCK [Dutta et al., 2017], Reluplex [Katz et al., 2017], Planet [ 8. Ehlers, 2017], MIP [ 9. Cheng et al., 2017; 10. [Lomuscio and Maganti, 2017] and BaB [Bunel et al., 2017]. 11. [Gergel et al., 2016] Victor Gergel, Vladimir Grishagin, and Alexander Gergel. Adaptive nested optimization scheme for multidimensional global search. Journal of Global Op- timization, 66(1):35-51, 2016. 12. [Goodfellow et al., 2014] I. J. Goodfellow, J. Shlens, and C. Szegedy. Explaining and Harnessing Adversarial Ex- amples. ArXiv e-prints, December 2014. 13. [Grishagin et al., 2018] Vladimir Grishagin, Ruslan Is- rafilov, and Yaroslav Sergeyev. Convergence conditions and numerical comparison of global optimization methods based on dimensionality reduction schemes. Applied Mathematics and Computation, 318:270-280, 2018. 14. [Huang et al., 2017] Xiaowei Huang, Marta Kwiatkowska, Sen Wang, and Min Wu. Safety verification of deep neu- ral networks. In Computer Aided Verification, pages 3-29. Springer Berlin Heidelberg, 2017. 15. [Katz et al., 2017] Guy Katz, Clark Barrett, David Dill, Kyle Julian, and Mykel Kochenderfer. Reluplex: An efficient smt solver for verifying deep neural networks. arXiv preprint arXiv:1702.01135, 2017. 16. [Kolter and Wong, 2017] J Zico Kolter and Eric Wong. Provable defenses against adversarial examples via the convex outer adversarial polytope. arXiv preprint arXiv:1711.00851, 2017. 17. [Lomuscio and Maganti, 2017] Alessio Lomuscio and Lalit Maganti. An approach to reachability analysis for feed- forward relu neural networks. CoRR, abs/1706.07351, 2017. 18. [Moosavi-Dezfooli et al., 2016] Seyed-Mohsen Moosavi- Dezfooli, Alhussein Fawzi, Omar Fawzi, and Pascal Frossard. Universal adversarial perturbations. CoRR, abs/1610.08401, 2016. 19. [Narodytska et al., 2017] Nina Narodytska, Shiva Prasad Kasiviswanathan, Leonid Ryzhyk, Mooly Sagiv, and Toby Walsh. Verifying properties of binarized deep neural net- works. CoRR, abs/1709.06662, 2017. 20. [Nguyen et al., 2014] Anh Mai Nguyen, Jason Yosinski, and Jeff Clune. Deep neural networks are easily fooled: High confidence predictions for unrecognizable images. CoRR, abs/1412.1897, 2014. 21. [Papernot et al., 2015] Nicolas Papernot, Patrick D. Mc- Daniel, Somesh Jha, Matt Fredrikson, Z. Berkay Celik, and Ananthram Swami. The limitations of deep learning in adversarial settings. CoRR, abs/1511.07528, 2015. 22. [Piyavskii, 1972] SA Piyavskii. An algorithm for finding the absolute extremum of a function. USSR Computa- tional Mathematics and Mathematical Physics, 12(4):57- 67, 1972. 23. [Pulina and Tacchella, 2010] Luca Pulina and Armando Tac- chella. An abstraction-refinement approach to verification of artificial neural networks. In Computer Aided Verifica- tion, pages 243-257. Springer Berlin Heidelberg, 2010. 24. [Ruan et al., 2018a] Wenjie Ruan, Xiaowei Huang, and Marta Kwiatkowska. Reachability analysis of deep neu- ral networks with provable guarantees. arXiv preprint arXiv:1805.02242, 2018. 25. [Ruan et al., 2018b] Wenjie Ruan, Min Wu, Youcheng Sun, Xiaowei Huang, Daniel Kroening, and Marta Kwiatkowska. Global robustness evaluation of deep neu- ral networks with provable guarantees for L0 norm. arXiv preprint arXiv:1804.05805, 2018. 26. [Sohrab, 2003] Houshang H Sohrab. Basic real analysis, volume 231. Springer, 2003. 27. [Sun et al., 2018] Youcheng Sun, Min Wu, Wenjie Ruan, Xi- aowei Huang, Marta Kwiatkowska, and Daniel Kroening. Concolic testing for deep neural networks. arXiv preprint arXiv:1805.00089v1, 2018. 28. [Szegedy et al., 2013] Christian Szegedy, Wojciech Zaremba, Ilya Sutskever, Joan Bruna, Dumitru Erhan, Ian Goodfellow, and Rob Fergus. Intriguing properties of neural networks. arXiv:1312.6199v4, 2013. 29. [Torn and Zilinskas, 1989] Aimo Torn and Antanas Zilin- skas. Global Optimization. Springer-Verlag New York, Inc., New York, NY, USA, 1989. 30. [Wicker et al., 2018] Matthew Wicker, Xiaowei Huang, and Marta Kwiatkowska. Feature-guided black-box safety testing of deep neural networks. In Proc. 24th Interna- tional Conference on Tools and Algorithms for the Con- struction and Analysis of Systems (TACAS'18), pages 408- 426, 2018. 31. [Xiang et al., 2017] Weiming Xiang, Hoang-Dung Tran, and Taylor T Johnson. Output reachable set estimation and verification for multi-layer neural networks. arXiv preprint arXiv:1708.03322, 2017.