NEURAL NETWORK BRANCHING FOR NEURAL NETWORK VERIFICATION
Year: 2020

Abstract:
Formal verification of neural networks is essential for their deployment in safetycritical areas. Many available formal verification methods have been shown to be instances of a unified Branch and Bound (BaB) formulation. We propose a novel framework for designing an effective branching strategy for BaB. Specifically, we learn a graph neural network (GNN) to imitate the strong branching heuristic behaviour. Our framework differs from previous methods for learning to branch in two main aspects. Firstly, our framework directly treats the neural network we want to verify as a graph input for the GNN. Secondly, we develop an intuitive forward and backward embedding update schedule. Empirically, our framework achieves roughly 50% reduction in both the number of branches and the time required for verification on various convolutional networks when compared to the best available hand-designed branching strategy. In addition, we show that our GNN model enjoys both horizontal and vertical transferability. Horizontally, the model trained on easy properties performs well on properties of increased difficulty levels. Vertically, the model trained on small neural networks achieves similar performance on large neural networks. Code for all experiments is available at https://github.com/oval-group/GNN_branching.

heading: 1 INTRODUCTION
text: Despite their outstanding performances on various tasks, neural networks are found to be vulnerable to adversarial examples (Goodfellow et al., 2015; Szegedy et al., 2013). The brittleness of neural networks can have costly consequences in areas such as autonomous driving, finance and healthcare. When one requires robustness to adversarial examples, traditional model evaluation approaches, which test the trained model on a hold-out set, do not suffice. Instead, formal verification of properties such as adversarial robustness becomes necessary. For instance, to ensure self-driving cars make consistent correct decisions even when the input image is slightly perturbed, the required property to verify is that the underlying neural network outputs the same correct prediction for all points within a norm ball whose radius is determined by the maximum perturbation allowed.
Several methods have been proposed for verifying properties on neural networks (NN). Bunel et al. (2018) showed that many of the available methods can be viewed as instances of a unified BaB framework. A BaB algorithm consists of two key components: branching strategies and bounding methods. Branching strategies decide how the search space is recursively split into smaller spaces. Bounding methods compute bounds of each subspace to tighten the bounds of the final objective function over the whole search space. In this work, we focus on improving the branching strategies. By directly working with a general framework, our identified algorithmic improvements can be combined with any bounding method, leading to potential performance improvement for BaB based verification algorithms.
Branching strategies have significant impacts on the overall problem-solving process, as it directly decides the total number of steps, and consequently the total time, required to solve the problem at hand. The quality of a branching strategy is even more important when NN verification problems are considered, which generally have a very large search space. Each input dimension or each activation unit can be a potential branching option and neural networks of interest often have high dimensional inputs and thousands of hidden activation units. With such a large search space, an effective branching strategy could mean a large reduction of the total number of branches required, and consequently of the time required to solve a problem. Developing an effective strategy is thus of significant importance to the success of BaB based NN verification.
So far, to the best of our knowledge, branching rules adopted by BaB based verification methods are either random selection (Katz et al., 2017; Ehlers, 2017) or hand-designed heuristics (Wang et al.,
2018b; Bunel et al., 2018; Royo et al., 2019; Bunel et al., 2019). Random selection is generally inefficient as the distribution of the best branching decision is rarely uniform. In practice, this strategy often results in exhaustive search to make a verification decision. On the other hand, hand designed heuristics often involve a trade-off between effectiveness and computational cost. For instance, strong branching is generally one of the best performing heuristics for BaB methods in terms of the number of branches, but it is computationally prohibitive as each branching decision requires an expensive exhaustive search over all possible options. The heuristics that are currently used in practice are either inspired by the corresponding dual problem when verification is formulated as an optimization problem (Bunel et al., 2018; Royo et al., 2019) or incorporating the gradient information of the neural network (Wang et al., 2018b). These heuristics normally have better computational efficiency. However, given the complex nature of the search space, it is unlikely that any hand-designed heuristic is able to fully exploit the structure of the problem and the data distribution encountered in practice. As mentioned earlier, for large size NN verification problems, a slight reduction in the quality of the branching strategy could lead to substantial increase in the total number of branches required to solve the problem. A computationally cheap but high quality branching strategy is thus much needed.
In order to exploit the inherent structure of the problem and the data, we propose a novel machine learning framework for designing a branching strategy. Our framework is both computationally efficient and effective, giving branching decisions that are of a similar quality to that of strong branching. Specifically, we make the following contributions:
• We use a graph neural network (GNN) to exploit the structure of the neural network we want to verify. The embedding vectors of the GNN are updated by a novel schedule, which is both computationally cheap and memory efficient. In detail, we mimic the forward and backward passes of the neural network to update the embedding vectors. In addition, the proposed GNN allows a customised schedule to update embedding vectors via shared parameters. That means, once training is done, the trained GNN model is applicable to various verification properties on different neural network structures.
• We train GNNs via supervised learning. We provide ways to generate training data cheaply but inclusive enough to represent branching problems at different stages of a BaB process for various verification properties. With the ability to exploit the neural network structure and a comprehensive training data set, our GNN is easy to train and converges fast.
• Our learned GNN also enjoys transferability both horizontally and vertically. Horizontally, although trained with easy properties, the learned GNN gives similar performance on medium and difficult level properties. More importantly, vertically, given that all other parts of BaB algorithms remain the same, the GNN trained on small networks performs well on large networks. Since the network size determines the total cost for generating training data and is positively correlated with the difficulty of learning, this vertical transferability allows our framework to be readily applicable to large scale problems.
• We further enhance our framework via online learning. For a learned branching strategy, it is expected that the strategy can fail to output satisfactory branching decisions from time to time. To deal with this issue, we provide an online scheme for fine-tuning the GNN along the BaB process in order to best accommodate the verification property at hand.
• Finally, we supply a dataset on convolutional NN verification problems, covering problems at different difficulty levels over neural networks of different sizes. We hope that by providing a large problem dataset it could allow easy comparisons among existing methods and additionally encourage the development of better methods.
Since most verification methods available work on ReLU-based deep neural networks, we focus on neural networks with ReLU activation units in this paper. However, we point out that our framework is applicable to any neural network architecture.

heading: 2 RELATED WORKS
text: Learning has already been used in solving combinatorial optimization problems (Bello et al., 2016; Dai et al., 2017) and mixed integer linear programs (MILP) (Khalil et al., 2016; Alvarez et al., 2017; Hansknecht et al., 2018; Gasse et al., 2019). In these areas, instances of the same underlying structure are solved multiple times with different data values, which opens the door for learning. Among them, Khalil et al. (2016), Alvarez et al. (2017), Hansknecht et al. (2018), and Gasse et al. (2019) proposed learned branching strategies for solving MILP with BaB algorithms. These meth-
ods imitate the strong branching strategy. Specifically, Khalil et al. (2016) and Hansknecht et al. (2018) learn a ranking function to rank potential branching decisions while Alvarez et al. (2017) uses regression to assign a branching score to each potential branching choice. Apart from imitation, Anderson et al. (2019) proposed utilizing Bayesian optimization to learn verification policies. There are two main issues with these methods. Firstly, they rely heavily on hand-designed features or priors and secondly, they use a generic learning structure which is unable to exploit the neural network architecture.
The approach most relevant to ours is the concurrent work by Gasse et al. (2019). They managed to reduce feature reliance by exploiting the bipartite structure of an MILP through a GNN. The bipartite graph is capable of capturing the network architecture, but cannot exploit it effectively. Specifically, it treats all the constraints the same and updates them simultaneously using the same set of parameters. This limited expressiveness can result in a difficulty in learning and hence in a high generalization error for NN verification problems. Our proposed framework is specifically designed for NN verification problems. By exploiting the neural network structure, and designing a customized schedule for embedding updates, our framework is able to scale elegantly both in terms of computation and memory. Finally, we mention that the recently proposed verification methods (Katz et al., 2019; Singh et al., 2018; Anderson et al., 2019) are not explicitly formulated as BaBs. Since our focus is on branching, we mainly use the methods in Bunel et al. (2019) for comparison.

heading: 3 BACKGROUND
text: Formal verification of neural networks refers to the problem of proving or disproving a property over a bounded input domain. Properties are functions of neural network outputs. When a property can be expressed as a Boolean expression over linear forms, we can modify the neural network in a suitable way so that the property can be simplified to checking the sign of the neural network output (Bunel et al., 2018). Note that all the properties studied in previous works satisfy this form, thereby allowing us to use the aforementioned simplification. Mathematically, given the modified neural network f , a bounded input domain C, formal verification examines the truthfulness of the following statement: ∀x ∈ C, f(x) ≥ 0. (1) If the above statement is true, the property holds. Otherwise, the property does not hold.

heading: 3.1 BRANCH AND BOUND
text: Verification tasks are often treated as a global optimization problem. We want to find the minimum of f(x) over C in order to compare it with the threshold 0. Specifically, we consider an L layer feedforward neural network, f : R|x| → R, with ReLU activation units such that for any x0 ∈ C ⊂ R|x|, f(x0) = x̂L ∈ R where
x̂i+1 =W i+1xi + b i+1, for i = 0, . . . , L− 1, (2a) xi = max(x̂i, 0), for i = 1, . . . , L− 1. (2b)
The termsW i and bi refer to the weights and biases of the i-th layer of the neural network f. Domain C can be an `p norm ball with radius . Finding the minimum of f is a challenging task, as the optimization problem is generally NP hard (Katz et al., 2017). To deal with the inherent difficulty of the optimization problem itself, BaB (Bunel et al., 2018) is generally adopted. In detail, BaB based methods divide C into sub-domains, each of which defines a new sub-problem (branching). They then compute a relaxed lower bound of the minimum on each sub-problem (bounding). The minimum of the lower bounds of all the generated sub-domains constitutes a valid global lower bound of the global minimum over C. As a recursive process, BaB keeps partitioning the subdomains to tighten the global lower bound. The process terminates when the computed global lower bound is above zero (property is true) or when an input with a negative output is found (property is false). A detailed description of the BaB is provided in the appendices. In what follows, we provide a brief description of the two components, bounding methods and branching strategies, that is necessary for the understanding of our novel learning framework.

heading: 3.2 BOUNDING
text: For NN verification problems, bounding consists of finding upper and lower bounds for the final output, the minimum of f(x) over C. An effective technique to compute a lower bound is to transform the original optimization problem into a linear program (LP) by introducing convex relaxations over ReLU activation units. As we can see in Eq. (2b), ReLU activation units do not define a convex feasible set, and hence, relaxations are needed. Denote the j-th element of the vector xi as xi[j].
Possible convex relaxations for a hidden node xi[j] that have been introduced so far are shown in Figure 1. We replace ReLU with the shaded green area. The tighter the convex relaxation introduced, the more computational expensive it is to compute a bound but the tighter the bound is going to be. From Figure 1, we note that in order to introduce a convex relaxation, we need intermediate bounds li[j] and ui[j]. Thus intermediate bounds are required for building the LP for the final output lower bound. Given their purpose and the large number of intermediate bound computations, rough estimations are mainly used. On the other hand, the final output lower bound is directly used in making the pruning decision and hence a tighter lower bound is preferred as it avoids further unnecessary splits on the sub-problem.

heading: 3.3 BRANCHING
text: Branching is of equal importance as bounding in the BaB framework. Especially for large scale networks f , each branching step has a large number of putative choices. In these cases, the effectiveness of a branching strategy directly determines the possibility of verifying properties over these networks within a given time limit. On neural networks, two types of branching decisions are used: input domain split and hidden activation unit split.
Assume we want to split a parent domain D. Input domain split selects an input dimension and then makes a cut on the selected dimension while the rest of the dimensions remain the same. The common choice is to cut the selected dimension in half and the dimension to cut is decided by the branching strategy used. Available input domain split strategies are Bunel et al. (2018) and Royo et al. (2019). Royo et al. (2019)’s is based on sensitivity test of the LP onD while Bunel et al. (2018) use the formula provided in Wong & Kolter (2018) to estimate final output bounds for sub-domains after splitting on each input dimension and selects the dimension that results in the highest output lower bound estimates.
In our setting, we refer to a ReLU activation unit xi[j] = max(x̂i[j], 0) as ambiguous over D if the upper bound ui[j] and the lower bound li[j] for x̂i[j] have different signs. Activation unit split chooses among ambiguous activation units and then divides the original problem into cases of different activation phase of the chosen activation unit. If a branching decision is made on xi[j], we divide the ambiguous case into two determinable cases: {xi[j] = 0, li[j] ≤ x̂i[j] ≤ 0} and {xi[j] = x̂i[j], 0 ≤ x̂i[j] ≤ ui[j]}. After the split, the originally introduced convex relaxation is removed, since the above sets are themselves convex. We expect large improvements on the output lower bounds of the newly generated sub-problems if a good branching decision is made. Apart from random selection, employed in Ehlers (2017) and Katz et al. (2017), available ReLU split heuristics are Wang et al. (2018a) and Bunel et al. (2019). Wang et al. (2018a) compute scores based on gradient information to prioritise ambiguous ReLU nodes. Similarly, Bunel et al. (2019) use scores to rank ReLU nodes but scores are computed with a formula developed on the estimation equations in Wong & Kolter (2018). We note that for both branching strategies, after the split, intermediate bounds are updated accordingly on each new sub-problem. For NN verification problems, either domain split or ReLU split can be used at each branching step. When compared with each other, ReLU split is a more effective choice for large scale networks, as shown in Bunel et al. (2019).
All the aforementioned existing branching strategies use hand-designed heuristics. In contrast, we propose a new framework for branching strategies by utilizing a GNN to learn to imitate strong branching heuristics. This allows us to harness the effectiveness of strong branching strategies while retaining the efficiency of GPU computing power.

heading: 4 GNN FRAMEWORK
text: Overview We begin with a brief overview of our overall framework, followed by a detailed description of each of its components. A graph neural network G is represented by two components: a set of nodes V and a set of edges E, such that G = (V,E). Each node and each edge has its set of features. A GNN uses the graph structure and node and edge features to learn a representation vector (embedding vector) for each node v ∈ V . The GNN is a key component of our framework, in which we treat the neural network f as a graph Gf . A GNN takes Gf as an input and initializes an embedding vector for each node in V . The GNN updates each node’s embedding vector by aggregating its own node features and all its neighbours’ embedding vectors. After several rounds of updates, we obtain a learned representation (an updated embedding vector) of each node. To make a branching decision, we treat the updated embedding vectors as inputs to a score function, which assigns a score for each node that constitutes a potential branching option. A branching decision is made based on the scores of potential branching decision nodes. Our framework is visualised in Figure 2. We now describe each component in detail.

heading: Neural Network Graph Neural Network
text: Nodes Given a neural network f , V consists of all input nodes v0[j], all hidden activation nodes vi[j] and an output node vL. In our framework, we combine every pre-activation variable and its associated post-activation variable and treat them as a single node. Pre- and post-activation nodes together contain the information about the amount of convex relaxation introduced at this particular activation unit, so dealing with the combined node simplifies the learning process. In terms of the Eq. (2), let x′i[j] denote the combined node of x̂i[j] and xi[j]. The nodes v0[j], vi[j] and vL are thus in one-to-one correspondence with x0[j], x′i[j] and xL. We note that V is larger than the set of all potential branching decisions as it includes unambiguous activation nodes and output nodes.
Node features Different types of nodes have different sets of features. In particular, input node features contain the corresponding domain lower and upper bounds and the primal solution. For activation nodes, the node features consist of associated intermediate lower and upper bounds, the layer bias, primal and dual solutions and new terms computed using previous features. Finally, the output node has features including the associated output lower and upper bounds, the layer bias and the primal solution. Other types of features could be used and some features could be excluded if they are expensive to compute. We denote input node features as z0[j], activation node features as zi[j] and output node features as zL. Our framework uses simple node features and does not rely on extensive feature engineering. Nonetheless, by relying on the powerful GNN framework, it provides highly accurate branching decisions.
Edges E consists of all edges connecting nodes in V , which are exactly the connecting edges in f . Edges are characterized by the weight matrices that define the parameters of the network f such that for an edge eijk connecting x ′ i[j] and x ′ i+1[k], we assign e i jk =W i jk.
Embeddings We associate a p-dimensional embedding vector µv for each node v ∈ V . All embedding vectors are initialised as zero vectors.
Forward and Backward embedding updates In general, a graph neural network learns signals from a graph by acting as a function of two inputs: a feature matrix X ∈ R|V |×p, where each row is the embedding vector µv for a node v ∈ V , and an adjacency matrix A representing the graph structure. Under this formulation, all node embedding vectors are updated at the same time and there is no particular order between nodes. In this work, instead, we propose an update scheme where only the nodes corresponding to the same layer of the network f are updated at the same time, so embedding vector updates are carried out in a layer-by-layer forward-backward way.
We argue that the forward-backward updating scheme is a natural fit for our problem. In more detail, for a given problem D, each branching decision (an input node or an ambiguous activation node) will generate two sub-problems s1 and s2, with each sub-domain having an output lower bound lLs1 and l L s2 respectively, equal to or higher than l L D the lower bound that of D. Strong branching heuristic uses a predetermined function to measure the combined improvement of lLs1 and l L s2 over lLD and makes the final branching decision by selecting the node that gives the largest improvement. Thus, to maximise the performance of a graph neural network, we want a node embedding vector to maximally capture all information related to the computation of lLs1 and l L s2 . For estimating l L s1 , l L s2 of splitting on a potential branching decision node v, we note that these values are closely related to two factors. The first factor is the amount of convex relaxations introduced at a branching decision node v, when v corresponds to an ambiguous activation node. The second factor considers that the impact that splitting node v will have on the convex relaxations introduced to nodes on layers after that of v. Recall that, if there are no ambiguous activation nodes, the neural network f is simply a linear operator, whose minimum value can be easily obtained. When ambiguous activation nodes are present, the total amount of relaxation introduced determines the tightness of the lower bound to f . We thus treat embedding vectors as a measure of local convex relaxation and its contribution to other nodes’ convex relaxation.
As shown in Figure 1, at each ambiguous activation node x′i[j], the area of convex relaxation introduced is determined by the lower and upper bounds of the pre-activate node x̂i[j]. We observe that intermediate lower and upper bounds of a node x̂i[j] are significantly affected by the layers prior to it and have to be computed in a layer-by-layer fashion. Based on the observation, we utilise a forward layer-by-layer update on node embedding vectors. This should allow these embedding vectors to capture the local relaxation information. In terms of the impact of local relaxation change to that of other nodes, we note that by splitting an ambiguous node into two fixed cases, all intermediate bounds of nodes on later layers will be affected, leading to relaxation changes at those nodes. We thus employ a backward layer-by-layer update to account for the impact the local change has over other nodes. Theoretically, by fixing an ambiguous ReLU node, intermediate bounds of nodes at previous layers and on the same layer might change as well. For a naturally trained neural network, the changes for these nodes should be relatively small compared to nodes on the later layers. To account for these changes, we rely on multiple rounds of forward-and-backward updates.
In summary, during the forward update, for i = 1, . . . , L− 1, we have, for all possible j,
µ0[j] ←− Finp(z0[j];θ0), if µ0[j] = 0, (3)
µi[j] ←− Fact(zi[j],µi−1, ei;θ1), (4) µL ←− Fout(zL,µL−1, eL;θ2). (5)
During the backward update, for i = L− 1, . . . , 1, we have
µi[j] ←− Bact(zi[j],µi+1, ei+1;θ3), (6) µ0[j] ←− Binp(z0[j],µ1, e1;θ4). (7)
Update functions F andB take the form of multi-layered fully-connected networks with ReLU activation functions or composites of these simple update networks. The terms θi denote the parameters of the networks. A detailed description of update functions is provided in the appendices.
We point out that our forward-backward update scheme does not depend on the underlying neural network structure and thus should be generalizable to network architectures that differ from the one we use for training. However, it does rely on the information used to compute convex relaxations, so underlying data distribution, features and bounding methods are assumed to be the same when the trained model is applied to different networks. Furthermore, our forward-backward update is memory efficient, as we are dealing with one layer at a time and only the updated embedding vectors of the layer are used to update the embedding vectors in the next (forward-pass) and the previous (backward-pass) layer. This makes it readily applicable to large networks.
Scores At the end of the forward-backward updates, embedding vectors for potential branching decision nodes (all input nodes and ambiguous activation nodes) are gathered and treated as inputs of a score function gs(·;θ5) : Rp → R, which takes the form of a fully-connected network with parameters θ5. It assigns a scalar score for each input embedding vector. The final branching decision is determined by picking the node with the largest score.

heading: 5 PARAMETER ESTIMATION
text: Training We train a GNN via supervised learning. To estimate Θ := (θ0,θ1,θ2,θ3,θ4,θ5), we propose a new hinge rank loss function that is specifically designed for our framework. Before we give details of the loss, we introduce a relative improvement measure m first. Given a domain D, for each branching decision node v, the two generated sub-problems have output lower bounds lLs1 and lLs2 . We measure the relative improvement of splitting at the node v over the output lower bound lLD as follows
mv := (min(l L s1 , 0) + min(l L s2 , 0)− 2 · l L D)/(−2 · lLD). (8)
Intuitively, m (0 ≤ m ≤ 1) measures the average relative sub-problem lower bound improvement to the maximum improvement possible, that is −lLD. Any potential branching decision node v can be compared and ranked via its relative improvement value mv . Since we are only interested in branching nodes with large improvement measures, ranking loss is a natural choice. A direct pairwise rank loss might be difficult to learn for NN verification problems, given the large number of branching decision nodes on each domain D. In addition, many branching decisions may give similar performance, so it is redundant and potentially harmful to the learning process if we learn a ranking among these similar nodes. To deal with these issues, we develop our loss by first dividing all potential branching nodes into M classes (M is much smaller than the total number of branching decision nodes) through the improvement value mv of a node. We denote the class label as Yv for a node v. Labels are assigned in an ascending order such that Yv >= Yv′ if mv > mv′ . We then compute the pairwise hinge-rank loss on these newly assigned labels as
lossD(Θ) = 1
K N∑ i=1 ( N∑ j=1 φ(gs(µj ;Θ)− gs(µi;Θ)) · 1Yj>Yi ) , (9)
where φ(z) = (1 − z)+ is the hinge function, N is the total number of branching decision nodes and K is the total number of pairs where Yj > Yi for any branching decision nodes vi, vj . The loss measures the average hinge loss on score difference (gs(µj ;Θ) − gs(µi;Θ)) for all pairs of branching decision nodes vi, vj such that Yj > Yi. Finally, we evaluate Θ by solving the following optimization problem:
Θ = argmin Θ
λ 2 ‖Θ‖2 + 1 n n∑ i lossDi(Θ), (10)
where the lossDi is the one introduced in Eq. (9) and n is the number of training samples.
Fail-safe Strategy We introduce a fail-safe strategy employed by our framework to ensure that consistent high-quality branching decisions are made throughout a BaB process. The proposed framework uses a GNN to imitate the behavior of the strong branching heuristic. Although computationally cheap, in some cases, the output decision by the learned graph neural network might be suboptimal. When this happens, it could lead to considerably deteriorated performance for two reasons. Firstly, we observed that for certain problems, which requires multiple splits to reach a conclusion on this problem, if a few low-quality branching decisions are made at the beginning or the middle stage of the branching process, the total number of splits required might increase substantially. The total BaB path is thus, to some extent, sensitive to the quality of each branching decision apart from those made near the end of the BaB process. Secondly, once a low-quality decision is made on a given problem, a decision of similar quality is likely to be made on the two newly generated sub-problems, leading to exponential decrease in performance. Features for newly generated sub-problems are normally similar to those of the parent problem, especially in the cases where the branching decision of the parent problem is made on the later layers and loose intermediate bounds are used. Thus, it is reasonable to expect the GNN fails again on the resulting sub-problems.
To deal with this issue, we keep track of the output lower bound improvement for each branching decision, as introduced in Eq. (8). We then set a pre-determined threshold parameter. If the improvement is below the threshold, a computationally cheap heuristic is called to make a branching decision. Generally, the back-up heuristic is able to give an above-threshold improvement and generate sub-problems sufficiently different from the parent problem to allow the learned GNN to recover from the next step onwards.
Online Learning Online learning is a strategy to fine-tune the network for a particular property after we have learnt Θ. It can be seen as an extension of the fail-safe strategy employed. Every time a heuristic branching decision node vh is used instead of the node vgnn chosen by the GNN, we can use vh and vgnn to update the GNN accordingly. Since a correct GNN model should output an embedding vector µh resulting in a higher score gs(µh;Θ) for the heuristic decision, a loss can be developed based on the two scores gs(µh;Θ) and gs(µgnn;Θ) to generate optimization signals for correcting the GNN behaviour. For example, the loss used in our experimental setting is:
lossonline(Θ) = gs(µgnn;Θ)− gs(µh;Θ) + γ · ((mh −mgnn) > t). (11)
The last term is used to amplify (γ > 0) the loss if the relative improvement made by the heuristic decision is more than t percent higher than that by the GNN. We update Θ of the GNN by taking one gradient step with a small learning rate of the following minimization problem.
Θ = argmin Θ
λ 2 ‖Θ‖2 + lossonline(Θ). (12)
Online learning is property specific: it uses the decisions made by heuristics to fine tune the GNN model so it can best accommodate the property at hand. As will be shown in our experiments, a small but significant improvement in performance is achieved when online learning is used.

heading: 6 EXPERIMENTS
text: We now validate the effectiveness of our proposed framework through comparative experiments against other available NN verification methods. A comprehensive study of NN verification methods has been done in Bunel et al. (2019). We thus design our experiments based on the results presented in Bunel et al. (2019).

heading: 6.1 SETUP
text: We are interested in verifying properties on large network architectures with convolutional layers. In Bunel et al. (2019), existing NN methods are compared on a robustly trained convolutional network on MNIST. We adopt a similar network structure but using a more challenging dataset, namely CIFAR-10, for an increased difficulty level. We compare against the following two methods: (i) MIPplanet, a mixed integer solver backed by the commercial solver Gurobi; and (ii) BaBSR, a BaB based method utilising a ReLU-split heuristic. Our choice is motivated by their superior performance over other methods for MNIST verification problems in the previous work (Bunel et al., 2019).
We provide the detailed experimental setup through four perspectives: bounding methods, branching strategies, network structures, and verification properties tested. (Bounding methods) We compute intermediate bounds using linear bounds relaxations (Figure 1(b)). For the output lower bound, we use Planet relaxation (Figure 1(c)) and solve the corresponding LP with Gurobi. For the output upper bound, we compute it by directly evaluating the network value at the input provided by the LP solution. (Branching strategy) We focus on ReLU split only in our experiments. As shown in Bunel et al. (2019), domain split only outperforms ReLU split on low input dimensional and small scale networks. Also, since one of the Baseline method BaBSR employs a ReLU-split heuristic, we consider ReLU split only for a fair comparison. However, we emphasize that our framework is readily applicable to work with a combined domain and ReLU split strategy. (Network structures) Three neural network structures will be studied. The base one is of the similar structure and size to the one used in Bunel et al. (2019). It has two convolutional layers, followed by two fully connected layers and is trained robustly using the method provided in Wong & Kolter (2018). This particular choice of network size is made because the time required for solving each LP increases substantially with the size of the network. To best evaluate the performance of the branching strategy, we have to work with a medium sized network so that within the given timeout, a sufficient amount of branching decisions can be made to allow effective comparisons. When testing the transferability of the framework, two larger networks will be tested but their sizes are still restricted by the LP bottleneck. A detailed description of the network architecture is provided in the appendices. (Verification properties) Finally, we consider the following verification properties. Given an image x for which the model correctly predicted the label yc, we randomly choose a label yc′ such that for a given , we want to prove (e(c) − e(c′))T f ′(x′) > 0, ∀x′ s.t ‖x−x′‖∞ ≤ . Here, f ′ is the original neural network, e(c) and e(c
′) are one-hot encoding vectors for labels yc and yc′ . We want to verify that for a given , the trained network will not make a mistake by labelling the image as yc′ . Since BaBSR is claimed to be the best performing method on convolutional networks, we use it to determine the
values, which govern the difficulty level of verification properties. Small values mean that most ReLU activation units are fixed so their associated verification properties are easy to prove while large values could lead to easy detection of counter-examples. The most challenging values are those at which a large number of activation units are ambiguous. We use binary search with BaBSR method to find suitable values. We only consider values that result in true properties and timed out properties. Binary search process is simplified by our choice of robustly trained models. Since these models are trained to be robust over a δ ball, the predetermined value δ can be used as a starting value for binary search.

heading: 6.2 TRAINING DATASET
text: In order to generate training data, we firstly pick 565 random images and for each image, we randomly select an incorrect class. For each property, the value is determined by running binary search with BaBSR and 800 seconds timeout, so the final set of properties consists of mainly easily solvable properties and a limited number of timed out properties.
We collect training data along a BaB process for solving a verification property. At each given domain, given the large number of potential branching decisions, we perform the strong branching heuristic on a selected subset of all potential branching decisions. The subset consists of branching decisions that are estimated to be of high quality by the BaBSR heuristic and randomly selected ones, which ensure a minimum 5% coverage on each layer.
To construct a training dataset that is representative enough of the whole problem space, we need to cover a large number of properties. In addition, within a BaB framework, it is important to include training data at different stages of a BaB process. However, running a complete BaB process with the strong branching heuristic for hundreds of properties is computationally expensive and considerably time consuming. We thus propose the following procedure for generating a training dataset to guarantee a wide coverage both in terms of the verification properties and BaB stages. For generated verification properties, we randomly select 25% of non-timeout property to conduct a complete BaB process with the strong branching heuristic. For the rest of the properties, we try to generate at least B = 20 training data for each verification property. Given the maximum number of branches q = 10 and an effective and computationally cheap heuristic, we first generate a random integer k from [0, q]. Then, we run a BaB process with the selected cheap heuristic for k steps. Finally, we call the strong branching heuristic to generate a training sample. We repeat the process until B training samples are generated or the BaB process terminated. A detailed algorithm is provided in the appendices.

heading: 6.3 BASE MODEL
text: We test our learned model on the same model structure but on properties of three different difficulty levels. Testing verification properties are generated by binary search with BaBSR and 3600s timeout. We categorise verification properties solved within 800s as easy, which is consistent with training data generated, between 800s and 2400s as medium and more than 2400s as hard. In total, we generated 467 easy properties, 773 medium properties and 426 hard properties.
Results are given in the Table 1. Methods are compared in three perspectives: the average time over all properties, average number of branches required over the properties that are solved by all methods (we exclude timed out properties) and also the ratio of timed out properties. Since the properties are generated based on BaBSR, the timed out ratios of BaBSR on easy and medium properties are not comparable with that of other methods. All other numbers should give a fair evaluation of the effectiveness of our branching strategy. BaBSR, GNN and GNN-online only differ in the branching strategy used.
On all three sets of properties, we see that our learned branching strategy has led to a more than 50% reduction in the total average number of branches required for a property. As a direct result, the average time required achieves at least a 50% reduction as well. Our framework is thus an effective scheme and enjoys horizontal transferability. A further performance improvement is obtained through instance-specific online learning. Among all 1666 tested verification properties, GNN with online-learning solves 61.52% of properties with fewer number of branches and 60.20% of properties in less time when compared to the standard GNN.
We also provide a time cactus plot (Figure 3a) for all properties on the Base model. Time cactus plots for each category of properties can be found in the appendices. All these time cactus plots look
similar. Although BaBSR performs better than the commercial solver encoded method MIPplanet overall, MIPplanet wins on a subset of properties. The learned model GNN, however, is capable of giving consistent high quality performance over all properties tested.

heading: 6.4 TRANSFERABILITY: LARGER MODELS
text: We also robustly trained two larger networks. One has the same layer structure as the Base model but has more hidden units on each layer, which we refer to as the Wide model. The other has a similar number of hidden units on each layer but more layers. We refer to it as the Deep model. The detailed network architecture is provided in the appendices. Apart from the network structure, everything else is kept the same as for the Base model experiments. We use BaBSR and a timeout of 7200s to generate 300 properties for the Wide model and 250 properties for the Deep model. For these two models, each LP called for solving a sub-problem output lower bound is much more time consuming, especially for the Deep model. This is reason that the average number of branches considered is much fewer that those of the Base model within the given time limit.
The model learned on the Base network is tested on verification properties of large networks. Experimental results are given in the Table 2 and time cactus plots (Figures 3b, 3c) are also provided. All results are similar to what we observed on the Base model, which show that our framework enjoys vertical transferability.

heading: 7 DISCUSSION
text: The key observation of our work is that the neural network we wish to verify can be used to design a GNN to improve branching strategies. This observation can be used in enhancing the performances of other aspects of BaB. Possible future works include employing GNNs to find fast-converging starting values for solving LPs on a neural network and utilising GNNs to develop a lazy verifier, that only solves the corresponding LP on a domain when it could lead to pruning.

heading: ACKNOWLEDGMENTS
text: This work is supported by Clarendon Fund Scholarship. We thank Rudy Bunel for useful discussions. We thank Florian Jaeckle for proofreading the article.

heading: APPENDIX A. BRANCH AND BOUND ALGORITHM
text: The following generic Branch and Bound Algorithm is provided in Bunel et al. (2019). Given a neural network net and a verification property problem we wish to verify, the BaB procedure examines the truthfulness of the property through an iterative procedure. During each step of BaB, we first use the pick out function (line 6) to choose a problem prob to branch on. The split function (line 7) determines the branching strategy and splits the chosen problem prob into sub-problems. We compute output upper and lower bounds on each sub-problem with functions compute UB and compute LB respectively. Newly computed output upper bounds are used to tighten the global upper bound, which allows more sub-problems to be pruned. We prune a sub-problem if its output lower bound is greater than or equal to the global upper bound, so the smaller the global upper bound the better it is. Newly calculated output lower bounds are used to tighten the global lower bound, which is defined as the minimum of the output lower bounds of all remained sub-problems after pruning. We consider the BaB procedure converges when the difference between the global upper bound and the global lower bound is smaller than .
In our case, our interested verification problem Eq. (1) is a satisfiability problem. We thus can simplify the BaB procedure by initialising the global upper bound global ub as 0. As a result, we prune all sub-problems whose output lower bounds are above 0. In addition, the BaB procedure is terminated early when a below 0 output upper bound of a sub-problem is obtained, which means a counterexample exits.
Algorithm 1 Branch and Bound 1: function BAB(net, problem, ) 2: global lb← compute LB(net, problem) 3: global ub← compute UB(net, problem) 4: probs← [(global lb, problem)] 5: while global ub− global lb > do 6: ( , prob)← pick out(probs) 7: [subprob 1, . . . , subprob s]← split(prob) 8: for i = 1 . . . s do 9: sub lb← compute LB(net, subprob i) 10: sub ub← compute UB(net, subprob i) 11: if sub ub < global ub then 12: global ub← sub ub 13: prune probs(probs, global ub) 14: end if 15: if sub lb < global ub then 16: probs.append((sub lb, subprob i)) 17: end if 18: end for 19: global lb← min{lb | (lb, prob) ∈ probs} 20: end while 21: return global ub 22: end function

heading: APPENDIX B. IMPLEMENTATION OF FORWARD AND BACKWARD PASSES
text: We give implementation details of forward and backward updates for embedding vectors for the model used in the experiments section. Choices of forward and backward update functions are based on the bounding methods used. In our experiments, we used linear bound relaxations for computing intermediate bounds and Planet relaxation for computing the final output lower bound. We start with a graph neural network mimicking the structure of the network we want to verify. We denote domain lower and upper bounds as l0 and u0 respectively. Similarly, we denote the intermediate bounds (pre-activation) for layers i = 1, . . . , L − 1 as li and ui. Since an LP solver is called for the final output lower bound, we have primal values for all nodes of V and dual values for all ambiguous nodes of V . Finally, let W 1, . . . ,WL be the layer weights and b1, . . . , bL be the layer biases of the network f , which we wish to verify.

heading: B.1 FORWARD PASS
text: Unless otherwise stated, all functions F∗ are 2-layer fully connected network with ReLU activation units.
B.1.1 INPUT NODES
We update the embedding vectors of input nodes only during the first round of forward pass. That is we update µ0[j] when it is zero for all j. After that, input nodes embedding vectors are updated only in backward pass. For each input node, we form the feature vector z0[j] as a vector of l0[j], u0[j] and its associated primal solution. The input node embedding vectors are computed as
µ0[j] = Finp(z0[j];θ0). (13)

heading: B.1.2 ACTIVATION NODES
text: The update function Fact can be broken down into three parts: 1) compute information from local features 2) compute information from neighbourhood embedding vectors and 3) combine information from 1) and 2) to update current layer’s embedding vectors.
Information from local features Since we compute the final lower bound with the Planet relaxation (Figure 1(c)), we introduce a new feature related to the relaxation: the intercept of the relaxation triangle, shown in Figure 4. We denote an intercept as β and compute it as
βi[j] = −li[j] · ui[j] ui[j] − li[j] . (14)
The intercept of a relaxation triangle can be used as a measure of the amount of relaxation introduced at the current ambiguous node.
Therefore, the local feature vector zi[j] of an ambiguous node x′i[j] consists of li[j], ui[j], βi[j], its associated layer bias value, primal values (one for pre-activation variable and one for post-activation variable) and dual values. We obtain information from local features via
Ri[j] =
{ Fact−lf (zi[j];θ 0 1) if x ′ i[j] is ambiguous,
0 otherwise. (15)
where Ri[j] ∈ Rp.
Information from neighbourhood embedding vectors During the forward pass, we focus on embedding vectors of the previous layer only. To update an embedding vector on layer i, we first combine embedding vectors of the previous layer with edge weights via
Ei[j] = ∑ k W ikj · µi−1[k]. (16)
To compute the information from neighbourhood embedding vectors to an arbitrary activation node x′i[j], we consider each activation unit as a gate. We observe that the amount of the information
from neighbourhood embedding vectors that remains after passing through a gate is dependent on the its lower bound li[j] and upper bound ui[j]. When li[j] and ui[j] are of different signs, x′i[j] is an ambiguous node. With relaxation, for any input value between li[j] and ui[j], the maximum output achievable after passing an activation unit is shown by the red slope in Figure 5(a). The red slope si[j] is computed as
si[j](x̂i[j]) = ui[j]
ui[j] − li[j] · x̂i[j] + βi[j]. (17)
Thus, the amount of information from neighbourhood embedding vectors that remains after passing through an ambiguous gate is related to the ratio α := ui[j]ui[j]−li[j] . When ui[j] is no greater than zero, the activation node x′i[j] completely blocks all information. For any input value, the output value is zero after passing the activation unit, as shown by the red line in Figure 5(b). We have α = 0 in this case. Finally, when li[j] is no less than 0, the activation node x′i[j] allows a complete passing of information and α = 1. It is shown by the red line in Figure 5(c). We incorporate these observations into our evaluations and compute the information from neighbourhood embedding vectors as
Ni[j] = fact−nb([α · Ei[j], α′ · Ei[j]];θ11), (18)
where α′ = 1−α when 0 < α < 1 and α′ = α otherwise. Here, we use [a, b] to denote the concatenation of two vectors a, b ∈ Rp into a vector of R2p. We introduce α′ to be more informative. We do not consider the information that relate to the intercept βi[j] in the ambiguous case for the sake of simplicity. Improved performance could be expected if the βi[j] related information is incorporated as well.
Combine previous information Finally, we combine the information from local features and the information from neighbourhood embedding vectors to update the embedding vectors of activation nodes. Specifically,
µi[j] = Fact−com([Ri[j], Ni[j]];θ 2 1). (19)

heading: B.1.3 OUTPUT NODE
text: Embedding vectors of output nodes are updated in a similar fashion to that of activation nodes. We first compute information from local features.
RLj = Fout−lf (zLj ;θ 0 2) (20)
For output nodes, the vector of local features zL consists of output lower bound, output upper bound, primal solution and layer bias. Fout−lf is a one-layer fully-connected network with ReLU activation units. We then compute information from neighbourhood embedding vectors. Since the output node does not have an activation unit associated with it, we directly compute the information of neighbourhood embedding vectors as
EL[j] = ∑ k WLkj · µL−1[k]. (21)
Finally, we update the embedding vector of the output node as µLj = Fout−com([RL[j], EL[j]];θ 1 2). (22)

heading: B.2 BACKWARD PASS
text: During backward message passing, for i = L−1, . . . , 1, we update embedding vectors for activation nodes and input node. Again, all functions B∗ are 2-layer fully-connected networks unless specified otherwise.

heading: B.2.1 ACTIVATION NODES
text: Similar to updates of embedding vectors carried out for activation nodes in a forward pass, we update embedding vectors of activation nodes using the same three steps in the backward pass, but with minor modifications.
Information from local features We use the same feature zi[j] as the one used in the forward pass and compute the information from local features as
Rbi[j] =
{ Bact−lf1(zi[j];θ 0 3) if x ′ i[j] is ambiguous,
0 otherwise. (23)
We recall that a dual value indicates how the final objective function is affected if its associated constraint is relaxed by a unit. To better measure the importance of each relaxation to the final objective function, we further update the information from local features by
Rb ′
i[j] =
{ Bact−lf2([di[j] Rbi[j], R b i[j]];θ 1 3) if R b i[j] 6= 0
0 otherwise. (24)
Here, di[j] is the vector of dual values corresponding to the activation node x′i[j]. We use to mean that we multiply Rbi[j] by each element value of di[j] and concatenate them as a singe vector.
Information from neighbourhood embedding vectors During the backward pass, we focus on embedding vectors of the next layer only. In order to update an embedding vector on layer i, we compute the neighbourhood embedding vectors as
Ebi[j] = ∑ k W i+1jk · µi+1[k]. (25)
We point out that there might be an issue with computing Ei[j] if the layer i + 1 is a convolutional layer in the backward pass. For a convolutional layer, depending on the padding number, stride number and dilation number, each node x′i[j] may connect to a different number of nodes on the layer i + 1. Thus, to obtain a consistent measure of Ei[j], we divide Ei[j] by the number of connecting node on the layer i+ 1, denoted as Eb ′
i[j] and use the averaged E b′ i[j] instead. Let
Eb∗i[j] =
{ Eb ′
i[j] if layer i+1 convolutional, Ebi[j] otherwise.
(26)
The following steps are the same as the forward pass. We first evaluate
N bi[j] = Bact−nb([α · E b∗ i[j], α ′ · Eb∗i[j]];θ 2 3), (27)
and the update embedding vectors as
µi[j] = Bact−com([R b′ i[j], N b i[j]];θ 3 3). (28)
B.2.2 INPUT NODES
Finally, we update the input nodes. We use the feature vector zb0, which consists of domain upper bound and domain lower bound. Information from local features is evaluated as
R0j = Binp−lf (z b 0[j];θ 0 4). (29)
We compute the information from neighbourhood embedding vectors in the same manner as we do for activation nodes in the backward pass, shown in Eq (26). Denote the computed information as Eb∗0[j]. The embedding vectors of input nodes are updated by
µ0[j] = Binp−com([R b′ 0[j], E b∗ 0[j]];θ 1 4). (30)

heading: APPENDIX C. ALGORITHM FOR GENERATING TRAINING DATASET
text: Algorithm 2 outlines the procedure for generating the training dataset. The algorithm ensures the generated training date have a wide coverage both in terms of the verification properties and BaB stages while at the same time is computationally efficient. Specifically, we randomly pick 25% of all properties that do not time out and run a complete BaB procedure on each of them with the strong branching heuristic to generate training samples (line 3-5). For the remaining properties, we attempt to generate B training samples for each of them. To cover different stages of a BaB process of a property, we use a computationally cheap heuristic together with the strong branching heuristic. Given a property, we first use the cheap heuristic for k steps (line 10-15) to reach a new stage of the BaB procedure and then call the strong branching heuristic to generate a training sample (line 16). We repeat the process until B training samples are generated or the BaB processs terminates.
Algorithm 2 Generating Training Dataset 1: Provided: total P properties; minimum B training data for each property; a maximum q
branches between strong branching decisions 2: for p = 1, . . . , P do: 3: α←− random number from [0, 1] 4: if p is not a timed out property and α ≤ 0.25 then 5: Running a complete BaB process with the Strong Branching Heuristic 6: else 7: b = 0 8: while b ≤ B do 9: k ←− random integer from[0, q]
10: while k > 0 do 11: Call a computationally cheap heuristic 12: if BaB process terminates then return 13: end if 14: k = k − 1 15: end while 16: Call the strong branching heuristic and generate a training sample 17: if BaB process terminates then return 18: end if 19: b = b+ 1 20: end while 21: end if 22: end for

heading: APPENDIX D. EXPERIMENT DETAILS
text: All the hyper-parameters used in the experiments are determined by testing a small set of numbers over the validation set. Due to the limited number of tests, we believe better sets of hyper-parameters could be found.

heading: D.1 TRAINING DETAILS
text: Training dataset To generate a training dataset, 565 random images are selected. Binary serach with BaBSR and 800 seconds timeout are used to determine on the Base model. Among 565 verification properties determined, we use 430 properties to generate 17958 training samples and the rest of properties to generate 5923 validation samples. Training samples and validation samples are generated using Algorithm 2 with B = 20 and q = 10.
For a typical epsilon value, each sub-domain generally contains 1300 ambiguous ReLU nodes. Among them, approximately 140 ReLU nodes are chosen for strong branching heuristics, which leads to roughly 200 seconds for generating a training sample. We point out that the total amount of time required for generating a training sample equals the 2*(per LP solve time)*(number of ambiguous ReLU nodes chosen). Although both the second and the third terms increase with the size of the model used for generating training dataset, the vertical transferability of our GNN enables us to efficiently generate training dataset by working with a small substitute of the model we are interested in. In our case, we trained on the Base model and generalised to Wide and Deep model.
Training We initialise a GNN by assigning each node a 64-dimensional zero embedding vector. GNN updates embedding vectors through two rounds of forward and backward updates. To train the GNN, we use hinge rank loss (Eq. (9)) with M = 10. Parameters Θ are computed and updated through Adam optimizer with weight decay rate λ = 1e−4 and learning rate 1e−4. If the validation loss does not decrease for 10 consecutive epochs, we decrease the learning rate by a factor of 5. If the validation loss does not decrease for 20 consecutive epochs, we terminate the learning procedure. The batch size is set to 2. In our experiments, each training epoch took less than 400 seconds and the GNN converges within 60 epochs.
In terms of the training accuracy, we first evaluate each branching decision using the metric defined by Eq. (8) 1. Since there are several branching choices that give similar performance at each subdomain, we considered all branching choices that have mv above 0.9 as correct decisions. Under this assumption, our trained GNN achieves 85.8% accuracy on the training dataset and 83.1% accuracy on the validation dataset.
D.2 VERIFICATION EXPERIMENT DETAILS
We ran all verification experiments in parallel on 16 CPU cores, with one property being verified on one CPU core. We observed that although we specifically set the thread number to be one for MIPplanet (backed by the commercial solver Gurobi), the time required for solving a property depends on the total number of CPUs used. For a machine with 20 cpu cores, MIPplanet requires much less time on average for proving the same set of properties on fewer (say 4) CPU cores in parallel than on many (say 16) CPU cores in parallel (the rest of CPU cores remain idle). Since BaBSR, GNN and GNN-online all use Gurobi for the bounding problems, similar time variations, depending on the number of CPU cores used, are observed. We ran each method in the same setting and on 16 CPUs in parallel, so our reported results and time are comparable. However, we remind readers to take the time variation into consideration when replicating our experiments or using our results for comparison.
Fail-safe strategy Since, to the best of our knowledge, the branching heurisitc of BaBSR is the best performing one on convolutional neural networks so far, we choose it for our fail-safe strategy. The threshold is set to be 0.2. Every time when the relative improvementmgnn of a GNN branching decision vgnn is less than 0.2, we call the heuristic to make a new branching decision vh. We solve
1we have tried various other metrics, including picking the minimum of the two subdomain lower bounds and the maximum of the two lower bounds. Among these metrics, metric defined by Eq. (8) performs the best.
the corresponding LPs for the new branching decision and compute its relative improvement mh. The node with higher relative improvement is chosen to be the final branching decision.
Online learning We take a conservative approach in terms of online learning. We refer to a GNN decision as a failed decision if the relative improvement offered by heuristic branching is better than the one offered by the GNN. We record all GNN failed decisions and only update the GNN model online when the same failed decision is made at least twice. To update the GNN model, we use Adam optimizer with weight decay rate λ = 1e−4 and learning rate 1e−4. The GNN model is updated with one gradient step only with respect to the optimization problem Eq. (12), where γ = 1 and t = 0.1 in the loss function lossonline, defined in Eq. (11).

heading: D.3 BASELINES
text: We decided our baselines based on the experiment results of Bunel et al. (2019). In Bunel et al. (2019), methods including MIPplanet, BaBSR, planet (Ehlers, 2017), reluBaB and reluplex (Katz et al., 2017) are compared on a small convolutional MNIST network. Among them, BaBSR and MIPplanet significantly outperform other methods. We thus evaluate our methods against these two methods only in the experiments section. In order to strengthen our baseline, we compare against two additional methods here.
Neurify (Wang et al., 2018a) Similar to BaBSR, Neurify splits on ReLU activation nodes. It makes a branching decision by computing gradient scores to prioritise ReLU nodes. Since the updated version of Neurify’s released code supports verification, we conducted a comparison experiment between between Neurify and BaBSR for inclusiveness.
Neurify does not support CIFAR dataset. To evaluate the performance of Neurify, we obtained the trained ROBUST MNIST model and corresponding verification properties from Bunel et al. (2019). We ranked all verification properties in terms of the BaBSR solving time and selected the first 200 properties, which are solved by BaBSR within one minute, as our test properties. For a fair comparison, we have restricted Neurify to use one CPU core only and set the timeout limit to be two minutes. Among all test properties, Neurify timed out on 183 out of 200 properties. BaBSR thus outperforms Neurify significantly. Combining with the results of Bunel et al. (2019), BaBSR is indeed a fairly strong baseline to be compared against.
MIP based algorithm (Tjeng et al., 2019) We also compared our MIPplanet baseline against a new MIP based algorithm (Tjeng et al., 2019), published in ICLR 2019. To test these two methods, we randomly selected 100 verification properties from the CIFAR Base experiment with timeout 3600s. In terms of solving time, MIPplanet requires 1732.18 seconds on average while the new MIP algorithm requires 2736.60 seconds. Specifically, MIPplanet outperforms the new MIP algorithm on 78 out of 100 properties. MIPplanet is therefore a strong baseline for comparison.
As a caveat, we mention that the main difference between MIPplanet and the algorithm of (Tjeng et al., 2019) is the intermediate bound computation, which is complementary to our focus. If better intermediate bounds are shown to help verification, we can still use our approach to get better branching decisions corresponding to those bounds.

heading: D.4 MODEL ARCHITECTURE
text: We provide the architecture detail of the neural networks verified in the experiments in the following table.

heading: APPENDIX E. ADDITIONAL EXPERIMENT RESULTS

heading: E.1 FAIL-SAFE HEURISTIC DEPENDENCE
text: In all our experiments, we have compared against BaBSR, which employs only the fail-safe heuristic for branching. In other words, removing the GNN and using only the fail-safe heuristic is equivalent to BaBSR. The fact that GNN significantly outperforms BaBSR demonstrates that GNN is doing most of the job. To better evaluate the GNN’s reliance on a fail-safe heuristic, we study the ratio of times that a GNN branching decision is used for each verification property of a given model. Results are listed in Table 4. On all three models, GNN accounts for more than 90% of branching decisions employed on average, ensuring the effectiveness of our GNN framework.

heading: E.2 GNN FEATURE ANALYSIS
text: We evaluate the importance of different features used in GNN. We note that two types of features are used in GNN. The first type (including intermediates bounds, network weights and biases) can be collected at negligible costs. The other type is LP features (primal and dual values) that are acquired by solving a strong LP relaxation, which are expensive to compute but potentially highly informative. To evaluate their effect, we trained a new GNN with LP features removed and tested the new GNN on 260 randomly selected verification properties on the Base model. Among the selected properties, 140 are categorised as easy, 70 as medium and 50 as hard. We denote the model trained on all features as GNN and the newly trained model as GNN-R (we use R to indicate reduced features).
From Table 5, we observe that removing primal and dual information deteriorates the GNN performance, but GNN-R still outperforms the baseline heuristic BaBSR. We believe cheap features are the most important. Depending on the cost of LP, potential users can either remove expensive LP features or train a GNN with a smaller architecture.

heading: E.3 MIPPLANET BRANCHING NUMBER
text: MIPplanet is implemented with the commercial solver Gurobi. Since Gurobi outputs internal branch number, we recorded MIPplanet branch number for a subset of verification properties for each model. In detail, we randomly selected 120 properties of various difficulty levels for the Base model and 27 properties each for the Wide and Deep model. Results are summarised in Table 6.
One key observation we made is that Gurobi branch number is not positively related to the solving time. For instance, on timed out properties of the Wide model, MIPplanet branch number varies between 1 and 7479. We suspect Gurobi performs cutting before branching, so time spent on branching varies between properties, leading to inconsistent branch number and solving time. As the result, the MIPplanet branch number is not comparable with that of BaBSR, GNN and GNN-online. This is also the reason that we did not include MIPplant branch number in Table 1 and Table 2.

heading: E.4 LP SOLVING TIME AND GNN COMPUTING TIME
text: We mention that LP solving time is the main bottleneck for branch-and-bound based verification methods. Although both GNN evaluation time and LP solving time increase with the size of network, LP solving time grows at a significantly faster speed. For instance, in CIFAR experiments, GNN requires on average 0.02, 0.03, 0.08 seconds to make a branching decision on Base, Wide and Deep model respectively but the corresponding one LP solving time on average are roughly 1.1, 4.9, 9.6 seconds. GNN evaluation is almost negligible for large neural networks when compared to LP solving time.

heading: E.5 GEOMETRIC MEAN
text: For all our experiments, we based our analyses on the statistics of average solving time and branching number. To ensure the reported numbers are not biased by potential outliers, we measure methods’ performance with the geometric mean as well and summarize results in Table 7 and Table 8. Statistics of geometric mean are consistent with that of arithmetic mean, validating the analyses of the main paper.

heading: E.6 ADDITIONAL PLOTS
text: We provide cactus plots for the Base model on easy, medium and hard difficulty level properties respectively.

heading: APPENDIX F. MNIST DATASET
text: We replicate the CIFAR experiments on the MNIST dataset to test the generalization ability of our GNN framework.

heading: F.1 MODEL ARCHITECTURE AND VERIFICATION PROPERTIES
text: We trained three different networks on MNIST with the method provided in Wong & Kolter (2018). The base model is mainly used for generating the training dataset and testing the horizontal generalization ability of the trained GNN. The Wide and Deep models are used for evaluating the vertical generalization. Verification properties are found via binary search with BaBSR. We set the binary search time limit to be 1800 seconds for the Base model and 3600 seconds for the other two models.

heading: F.2 TRAINING DETAILS
text: Training dataset Training dataset is generated on the Base model. We point out that we explicitly choose a Base model of small network size for efficient and fast training data generation.
To generate a training dataset, 538 random MNIST images are selected. Binary serach with BaBSR and 600 seconds timeout are used to determine on the Base model. Among 538 verification properties determined, we use 403 properties to generate 18231 training samples and the rest of properties to generate 5921 validation samples. Training samples and validation samples are generated using Algorithm 2 with B = 20 and q = 10. For a typical epsilon value, each sub-domain generally contains 480 ambiguous ReLU nodes. Among them, approximately 80 ReLU nodes are chosen for strong branching heuristics, which leads to roughly 45 seconds for generating a training sample.
Training The same set of parameters and training procedure are used for training a GNN for MNIST dataset. The GNN converges in 70 epochs with each epoch took less than 400 seconds. The trained GNN reached 86.5% accuracy on the training dataset and 83.1% accuracy on the validation dataset.

heading: F.3 EXPERIMENT RESULTS
text: We first note that we use verification properties with timeout 1800 seconds on the Base model to allow for an integrated evaluation of GNN on both its performance and its horizontal transferability. Vertical transferability is tested on the Wide and Deep model.
We observe that MIPplanet outperforms all BaB based methods on verification properties of the Base model. Given that the network size of the Base model is particularly small (1226 hidden units only), we believe that MIP algorithms backed by commercial solvers could be the most effective tool on verification problems of small size. Our conjecture is further confirmed by the fact that MIPplanet timed out on almost all properties of both the Wide and Deep models. On all three models, GNN consistently outperforms BaBSR, demonstrating the transferability of our framework. Finally, when online learning is considered, we found it is effective in fine-tuning the trained GNN and enabling further performance improvements, especially on the Wide model.

heading: F.4 FAIL-SAFE HEURISTIC DEPENDENCE
text: Results of Table 11 ensure that the trained GNN is indeed account for the most branching decisions.

heading: F.5 GEOMETRIC MEAN
text: The consistency between the results of Table 12 and Table 10 confirm that our analyses based on arithmetic mean are not biased by outliers.

heading: F.6 TRANSFERABILITY BETWEEN DATASETS
text: To evaluate whether our framework can generalize further to support transferring between datasets, we tested CIFAR trained GNN on MNIST verification properties. In detail, we have tested the GNN on 20 randomly picked verification properties of MNIST Base model. We found that BABSR outperforms CIFAR trained GNN on all properties, so the CIFAR trained GNN model does not transfer to MNIST dataset. This is expected as MNIST and CIFAR images differ significantly from each other.

Review summaries:
This paper deals with complete formal verification of Neural Network, based on the Branch and Bound framework. The authors focus on branching strategies, which have been shown to be a critical design decision in order to obtain good performance. The tactic employed here is to learn a Graph Neural Network (which allows to transfer the heuristic from small networks to large networks), using supervised training to imitate strong branching. The authors also discuss fallback mechanism to prevent bad failures case, as well as an online fine-tuning strategy that provide better performance.

The paper proposes learning a branching heuristic to be used inside a branch-and-bound algorithm used for solving integer programming problems corresponding to neural network verification. The heuristic is parameterized as a neural network and trained to imitate an existing heuristic called Strong Branching which is computationally expensive but produces smaller branch-and-bound trees than other heuristics. A graph neural network architecture is used to take the neural network being verified as input, and a message passing schedule that follows a forward pass and a backward pass along the network being verified is used. An online learning variant is also considered that fine tunes the learned heuristic at test time as a problem instance is being solved. Results for verifying large convolutional neural networks on CIFAR-10 show approximately 2x improvement in average running time of the branch-and-bound algorithm.